Changes from branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90 at r8627 to trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90 at r4990
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90
r8627 r4990 18 18 USE oce ! ocean dynamics and tracers 19 19 USE dom_oce ! ocean space and time domain 20 USE phycst ! physical constants21 20 USE ldfdyn_oce ! ocean dynamics: lateral physics 22 21 ! 23 22 USE in_out_manager ! I/O manager 24 USE iom ! I/O library25 23 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 26 24 USE wrk_nemo ! Memory Allocation … … 77 75 INTEGER, INTENT(in) :: kt ! ocean time-step index 78 76 ! 79 INTEGER :: ji, jj, jk ! dummy loop indices 80 REAL(wp) :: zua, zva, zbt, ze2u, ze2v, zzz ! local scalar 81 REAL(wp), POINTER, DIMENSION(:,: ) :: zcu, zcv 82 REAL(wp), POINTER, DIMENSION(:,:,:) :: zuf, zut, zlu, zlv 83 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d ! 2D workspace 77 INTEGER :: ji, jj, jk ! dummy loop indices 78 REAL(wp) :: zua, zva, zbt, ze2u, ze2v ! temporary scalar 79 REAL(wp), POINTER, DIMENSION(:,: ) :: zcu, zcv 80 REAL(wp), POINTER, DIMENSION(:,:,:) :: zuf, zut, zlu, zlv 84 81 !!---------------------------------------------------------------------- 85 82 ! … … 115 112 DO jj = 2, jpjm1 116 113 DO ji = fs_2, fs_jpim1 ! vector opt. 117 zlu(ji,jj,jk) = - ( zuf(ji ,jj,jk) -zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) &118 & + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj ,jk) ) /e1u(ji,jj)119 120 zlv(ji,jj,jk) = + ( zuf(ji,jj ,jk) -zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) &121 & + ( hdivb(ji,jj+1,jk) - hdivb(ji ,jj,jk) ) /e2v(ji,jj)114 zlu(ji,jj,jk) = - ( zuf(ji,jj,jk) - zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) & 115 & + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj,jk) ) / e1u(ji,jj) 116 117 zlv(ji,jj,jk) = + ( zuf(ji,jj,jk) - zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) & 118 & + ( hdivb(ji,jj+1,jk) - hdivb(ji,jj,jk) ) / e2v(ji,jj) 122 119 END DO 123 120 END DO … … 126 123 DO ji = fs_2, fs_jpim1 ! vector opt. 127 124 zlu(ji,jj,jk) = - ( rotb (ji ,jj,jk) - rotb (ji,jj-1,jk) ) / e2u(ji,jj) & 128 & 125 & + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj ,jk) ) / e1u(ji,jj) 129 126 130 127 zlv(ji,jj,jk) = + ( rotb (ji,jj ,jk) - rotb (ji-1,jj,jk) ) / e1v(ji,jj) & 131 & 128 & + ( hdivb(ji,jj+1,jk) - hdivb(ji ,jj,jk) ) / e2v(ji,jj) 132 129 END DO 133 130 END DO … … 136 133 CALL lbc_lnk( zlu, 'U', -1. ) ; CALL lbc_lnk( zlv, 'V', -1. ) ! Boundary conditions 137 134 138 IF( iom_use('dispkexyfo') ) THEN ! ocean kinetic energy dissipation per unit area 139 ! ! due to xy friction (xy=lateral) 140 ! see NEMO_book appendix C, §C.7.2 (N.B. here averaged at t-points) 141 ! local dissipation of KE at t-point due to bilaplacian operator is given by : 142 ! + ahmu mi( zlu**2 ) + mj( ahmv zlv**2 ) 143 ! Note that a sign + is used as in v3.6 ahm is negative for bilaplacian viscosity 144 ! 145 ! NB: ahm is negative when bilaplacian is used 146 ALLOCATE( z2d(jpi,jpj) ) 147 z2d(:,:) = 0._wp 148 DO jk = 1, jpkm1 149 DO jj = 2, jpjm1 150 DO ji = 2, jpim1 151 z2d(ji,jj) = z2d(ji,jj) & 152 & + ( fsahmu(ji,jj,jk)*zlu(ji,jj,jk)**2 + fsahmu(ji-1,jj,jk)*zlu(ji-1,jj,jk)**2 & 153 & + fsahmv(ji,jj,jk)*zlv(ji,jj,jk)**2 + fsahmv(ji,jj-1,jk)*zlv(ji,jj-1,jk)**2 ) * tmask(ji,jj,jk) 154 END DO 155 END DO 156 END DO 157 zzz = 0.5_wp * rau0 158 z2d(:,:) = zzz * z2d(:,:) 159 CALL lbc_lnk( z2d,'T', 1. ) 160 CALL iom_put( 'dispkexyfo', z2d ) 161 DEALLOCATE( z2d ) 162 ENDIF 163 164 165 ! Third derivative 166 ! ---------------- 167 ! 135 168 136 DO jk = 1, jpkm1 169 ! 137 138 ! Third derivative 139 ! ---------------- 140 170 141 ! Multiply by the eddy viscosity coef. (at u- and v-points) 171 zlu(:,:,jk) = zlu(:,:,jk) * fsahmu(:,:,jk) 172 zlv(:,:,jk) = zlv(:,:,jk) * fsahmv(:,:,jk) 173 ! 142 zlu(:,:,jk) = zlu(:,:,jk) * ( fsahmu(:,:,jk) * (1-nkahm_smag) + nkahm_smag) 143 144 zlv(:,:,jk) = zlv(:,:,jk) * ( fsahmv(:,:,jk) * (1-nkahm_smag) + nkahm_smag) 145 174 146 ! Contravariant "laplacian" 175 147 zcu(:,:) = e1u(:,:) * zlu(:,:,jk) … … 180 152 DO ji = 1, fs_jpim1 ! vector opt. 181 153 zuf(ji,jj,jk) = fmask(ji,jj,jk) * ( zcv(ji+1,jj ) - zcv(ji,jj) & 182 & 183 & * fse3f(ji,jj,jk) * r1_e12f(ji,jj)154 & - zcu(ji ,jj+1) + zcu(ji,jj) ) & 155 & * fse3f(ji,jj,jk) / ( e1f(ji,jj)*e2f(ji,jj) ) 184 156 END DO 185 157 END DO … … 203 175 END DO 204 176 177 205 178 ! boundary conditions on the laplacian curl and div (zuf,zut) 206 179 !!bug gm no need to do this 2 following lbc... … … 208 181 CALL lbc_lnk( zut, 'T', 1. ) 209 182 210 DO jk = 1, jpkm1 ! Bilaplacian 183 DO jk = 1, jpkm1 184 185 ! Bilaplacian 186 ! ----------- 187 211 188 DO jj = 2, jpjm1 212 189 DO ji = fs_2, fs_jpim1 ! vector opt. … … 220 197 & + ( zut(ji,jj+1,jk) - zut(ji ,jj,jk) ) / e2v(ji,jj) 221 198 ! add it to the general momentum trends 222 ua(ji,jj,jk) = ua(ji,jj,jk) + zua 223 va(ji,jj,jk) = va(ji,jj,jk) + zva 199 ua(ji,jj,jk) = ua(ji,jj,jk) + zua * ( fsahmu(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag )) 200 va(ji,jj,jk) = va(ji,jj,jk) + zva * ( fsahmv(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag )) 224 201 END DO 225 202 END DO 203 226 204 ! ! =============== 227 205 END DO ! End of slab
Note: See TracChangeset
for help on using the changeset viewer.