Changeset 10928 for NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf_iso.F90
- Timestamp:
- 2019-05-03T17:44:56+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf_iso.F90
r10425 r10928 60 60 61 61 62 SUBROUTINE dyn_ldf_iso( kt )62 SUBROUTINE dyn_ldf_iso( kt, Kbb, Kmm, puu, pvv, Krhs ) 63 63 !!---------------------------------------------------------------------- 64 64 !! *** ROUTINE dyn_ldf_iso *** … … 81 81 !! horizontal fluxes associated with the rotated lateral mixing: 82 82 !! u-component: 83 !! ziut = ( ahmt + rn_ahm_b ) e2t * e3t / e1t di[ u b]84 !! - ahmt e2t * mi-1(uslp) dk[ mi(mk(u b)) ]85 !! zjuf = ( ahmf + rn_ahm_b ) e1f * e3f / e2f dj[ u b]86 !! - ahmf e1f * mi(vslp) dk[ mj(mk(u b)) ]83 !! ziut = ( ahmt + rn_ahm_b ) e2t * e3t / e1t di[ uu ] 84 !! - ahmt e2t * mi-1(uslp) dk[ mi(mk(uu)) ] 85 !! zjuf = ( ahmf + rn_ahm_b ) e1f * e3f / e2f dj[ uu ] 86 !! - ahmf e1f * mi(vslp) dk[ mj(mk(uu)) ] 87 87 !! v-component: 88 !! zivf = ( ahmf + rn_ahm_b ) e2t * e3t / e1t di[ v b]89 !! - ahmf e2t * mj(uslp) dk[ mi(mk(v b)) ]90 !! zjvt = ( ahmt + rn_ahm_b ) e1f * e3f / e2f dj[ ub]91 !! - ahmt e1f * mj-1(vslp) dk[ mj(mk(v b)) ]88 !! zivf = ( ahmf + rn_ahm_b ) e2t * e3t / e1t di[ vv ] 89 !! - ahmf e2t * mj(uslp) dk[ mi(mk(vv)) ] 90 !! zjvt = ( ahmt + rn_ahm_b ) e1f * e3f / e2f dj[ vv ] 91 !! - ahmt e1f * mj-1(vslp) dk[ mj(mk(vv)) ] 92 92 !! take the horizontal divergence of the fluxes: 93 93 !! diffu = 1/(e1u*e2u*e3u) { di [ ziut ] + dj-1[ zjuf ] } 94 94 !! diffv = 1/(e1v*e2v*e3v) { di-1[ zivf ] + dj [ zjvt ] } 95 !! Add this trend to the general trend (u a,va):96 !! u a = ua+ diffu95 !! Add this trend to the general trend (uu(rhs),vv(rhs)): 96 !! uu(rhs) = uu(rhs) + diffu 97 97 !! CAUTION: here the isopycnal part is with a coeff. of aht. This 98 98 !! should be modified for applications others than orca_r2 (!!bug) 99 99 !! 100 100 !! ** Action : 101 !! -( ua,va) updated with the before geopotential harmonic mixing trend101 !! -(puu(:,:,:,Krhs),pvv(:,:,:,Krhs)) updated with the before geopotential harmonic mixing trend 102 102 !! -(akzu,akzv) to accompt for the diagonal vertical component 103 103 !! of the rotated operator in dynzdf module 104 104 !!---------------------------------------------------------------------- 105 INTEGER, INTENT( in ) :: kt ! ocean time-step index 105 INTEGER , INTENT( in ) :: kt ! ocean time-step index 106 INTEGER , INTENT( in ) :: Kbb, Kmm, Krhs ! ocean time level indices 107 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 106 108 ! 107 109 INTEGER :: ji, jj, jk ! dummy loop indices … … 128 130 DO jj = 2, jpjm1 129 131 DO ji = 2, jpim1 130 uslp (ji,jj,jk) = - ( gdept _b(ji+1,jj,jk) - gdept_b(ji ,jj ,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk)131 vslp (ji,jj,jk) = - ( gdept _b(ji,jj+1,jk) - gdept_b(ji ,jj ,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk)132 wslpi(ji,jj,jk) = - ( gdepw _b(ji+1,jj,jk) - gdepw_b(ji-1,jj,jk) ) * r1_e1t(ji,jj) * tmask(ji,jj,jk) * 0.5133 wslpj(ji,jj,jk) = - ( gdepw _b(ji,jj+1,jk) - gdepw_b(ji,jj-1,jk) ) * r1_e2t(ji,jj) * tmask(ji,jj,jk) * 0.5132 uslp (ji,jj,jk) = - ( gdept(ji+1,jj,jk,Kbb) - gdept(ji ,jj ,jk,Kbb) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 133 vslp (ji,jj,jk) = - ( gdept(ji,jj+1,jk,Kbb) - gdept(ji ,jj ,jk,Kbb) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk) 134 wslpi(ji,jj,jk) = - ( gdepw(ji+1,jj,jk,Kbb) - gdepw(ji-1,jj,jk,Kbb) ) * r1_e1t(ji,jj) * tmask(ji,jj,jk) * 0.5 135 wslpj(ji,jj,jk) = - ( gdepw(ji,jj+1,jk,Kbb) - gdepw(ji,jj-1,jk,Kbb) ) * r1_e2t(ji,jj) * tmask(ji,jj,jk) * 0.5 134 136 END DO 135 137 END DO … … 151 153 ! zdkv(jk=1)=zdkv(jk=2) 152 154 153 zdk1u(:,:) = ( ub(:,:,jk) -ub(:,:,jk+1) ) * umask(:,:,jk+1)154 zdk1v(:,:) = ( vb(:,:,jk) -vb(:,:,jk+1) ) * vmask(:,:,jk+1)155 zdk1u(:,:) = ( puu(:,:,jk,Kbb) -puu(:,:,jk+1,Kbb) ) * umask(:,:,jk+1) 156 zdk1v(:,:) = ( pvv(:,:,jk,Kbb) -pvv(:,:,jk+1,Kbb) ) * vmask(:,:,jk+1) 155 157 156 158 IF( jk == 1 ) THEN … … 158 160 zdkv(:,:) = zdk1v(:,:) 159 161 ELSE 160 zdku(:,:) = ( ub(:,:,jk-1) - ub(:,:,jk) ) * umask(:,:,jk)161 zdkv(:,:) = ( vb(:,:,jk-1) - vb(:,:,jk) ) * vmask(:,:,jk)162 zdku(:,:) = ( puu(:,:,jk-1,Kbb) - puu(:,:,jk,Kbb) ) * umask(:,:,jk) 163 zdkv(:,:) = ( pvv(:,:,jk-1,Kbb) - pvv(:,:,jk,Kbb) ) * vmask(:,:,jk) 162 164 ENDIF 163 165 … … 171 173 DO jj = 2, jpjm1 172 174 DO ji = fs_2, jpi ! vector opt. 173 zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * MIN( e3u _n(ji,jj,jk), e3u_n(ji-1,jj,jk) ) * r1_e1t(ji,jj)175 zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * MIN( e3u(ji,jj,jk,Kmm), e3u(ji-1,jj,jk,Kmm) ) * r1_e1t(ji,jj) 174 176 175 177 zmskt = 1._wp / MAX( umask(ji-1,jj,jk )+umask(ji,jj,jk+1) & … … 178 180 zcof1 = - zaht_0 * e2t(ji,jj) * zmskt * 0.5 * ( uslp(ji-1,jj,jk) + uslp(ji,jj,jk) ) 179 181 180 ziut(ji,jj) = ( zabe1 * ( ub(ji,jj,jk) - ub(ji-1,jj,jk) ) &182 ziut(ji,jj) = ( zabe1 * ( puu(ji,jj,jk,Kbb) - puu(ji-1,jj,jk,Kbb) ) & 181 183 & + zcof1 * ( zdku (ji,jj) + zdk1u(ji-1,jj) & 182 184 & +zdk1u(ji,jj) + zdku (ji-1,jj) ) ) * tmask(ji,jj,jk) … … 186 188 DO jj = 2, jpjm1 187 189 DO ji = fs_2, jpi ! vector opt. 188 zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * e3t _n(ji,jj,jk) * r1_e1t(ji,jj)190 zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e1t(ji,jj) 189 191 190 192 zmskt = 1._wp / MAX( umask(ji-1,jj,jk ) + umask(ji,jj,jk+1) & … … 193 195 zcof1 = - zaht_0 * e2t(ji,jj) * zmskt * 0.5 * ( uslp(ji-1,jj,jk) + uslp(ji,jj,jk) ) 194 196 195 ziut(ji,jj) = ( zabe1 * ( ub(ji,jj,jk) - ub(ji-1,jj,jk) ) &197 ziut(ji,jj) = ( zabe1 * ( puu(ji,jj,jk,Kbb) - puu(ji-1,jj,jk,Kbb) ) & 196 198 & + zcof1 * ( zdku (ji,jj) + zdk1u(ji-1,jj) & 197 199 & +zdk1u(ji,jj) + zdku (ji-1,jj) ) ) * tmask(ji,jj,jk) … … 203 205 DO jj = 1, jpjm1 204 206 DO ji = 1, fs_jpim1 ! vector opt. 205 zabe2 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e1f(ji,jj) * e3f _n(ji,jj,jk) * r1_e2f(ji,jj)207 zabe2 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e1f(ji,jj) * e3f(ji,jj,jk) * r1_e2f(ji,jj) 206 208 207 209 zmskf = 1._wp / MAX( umask(ji,jj+1,jk )+umask(ji,jj,jk+1) & … … 210 212 zcof2 = - zaht_0 * e1f(ji,jj) * zmskf * 0.5 * ( vslp(ji+1,jj,jk) + vslp(ji,jj,jk) ) 211 213 212 zjuf(ji,jj) = ( zabe2 * ( ub(ji,jj+1,jk) - ub(ji,jj,jk) ) &214 zjuf(ji,jj) = ( zabe2 * ( puu(ji,jj+1,jk,Kbb) - puu(ji,jj,jk,Kbb) ) & 213 215 & + zcof2 * ( zdku (ji,jj+1) + zdk1u(ji,jj) & 214 216 & +zdk1u(ji,jj+1) + zdku (ji,jj) ) ) * fmask(ji,jj,jk) … … 224 226 DO jj = 2, jpjm1 225 227 DO ji = 1, fs_jpim1 ! vector opt. 226 zabe1 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e2f(ji,jj) * e3f _n(ji,jj,jk) * r1_e1f(ji,jj)228 zabe1 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e2f(ji,jj) * e3f(ji,jj,jk) * r1_e1f(ji,jj) 227 229 228 230 zmskf = 1._wp / MAX( vmask(ji+1,jj,jk )+vmask(ji,jj,jk+1) & … … 231 233 zcof1 = - zaht_0 * e2f(ji,jj) * zmskf * 0.5 * ( uslp(ji,jj+1,jk) + uslp(ji,jj,jk) ) 232 234 233 zivf(ji,jj) = ( zabe1 * ( vb(ji+1,jj,jk) - vb(ji,jj,jk) ) &235 zivf(ji,jj) = ( zabe1 * ( pvv(ji+1,jj,jk,Kbb) - pvv(ji,jj,jk,Kbb) ) & 234 236 & + zcof1 * ( zdkv (ji,jj) + zdk1v(ji+1,jj) & 235 237 & + zdk1v(ji,jj) + zdkv (ji+1,jj) ) ) * fmask(ji,jj,jk) … … 241 243 DO jj = 2, jpj 242 244 DO ji = 1, fs_jpim1 ! vector opt. 243 zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * MIN( e3v _n(ji,jj,jk), e3v_n(ji,jj-1,jk) ) * r1_e2t(ji,jj)245 zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * MIN( e3v(ji,jj,jk,Kmm), e3v(ji,jj-1,jk,Kmm) ) * r1_e2t(ji,jj) 244 246 245 247 zmskt = 1._wp / MAX( vmask(ji,jj-1,jk )+vmask(ji,jj,jk+1) & … … 248 250 zcof2 = - zaht_0 * e1t(ji,jj) * zmskt * 0.5 * ( vslp(ji,jj-1,jk) + vslp(ji,jj,jk) ) 249 251 250 zjvt(ji,jj) = ( zabe2 * ( vb(ji,jj,jk) - vb(ji,jj-1,jk) ) &252 zjvt(ji,jj) = ( zabe2 * ( pvv(ji,jj,jk,Kbb) - pvv(ji,jj-1,jk,Kbb) ) & 251 253 & + zcof2 * ( zdkv (ji,jj-1) + zdk1v(ji,jj) & 252 254 & +zdk1v(ji,jj-1) + zdkv (ji,jj) ) ) * tmask(ji,jj,jk) … … 256 258 DO jj = 2, jpj 257 259 DO ji = 1, fs_jpim1 ! vector opt. 258 zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * e3t _n(ji,jj,jk) * r1_e2t(ji,jj)260 zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e2t(ji,jj) 259 261 260 262 zmskt = 1./MAX( vmask(ji,jj-1,jk )+vmask(ji,jj,jk+1) & … … 263 265 zcof2 = - zaht_0 * e1t(ji,jj) * zmskt * 0.5 * ( vslp(ji,jj-1,jk) + vslp(ji,jj,jk) ) 264 266 265 zjvt(ji,jj) = ( zabe2 * ( vb(ji,jj,jk) - vb(ji,jj-1,jk) ) &267 zjvt(ji,jj) = ( zabe2 * ( pvv(ji,jj,jk,Kbb) - pvv(ji,jj-1,jk,Kbb) ) & 266 268 & + zcof2 * ( zdkv (ji,jj-1) + zdk1v(ji,jj) & 267 269 & +zdk1v(ji,jj-1) + zdkv (ji,jj) ) ) * tmask(ji,jj,jk) … … 275 277 DO jj = 2, jpjm1 276 278 DO ji = 2, jpim1 !!gm Question vectop possible??? !!bug 277 ua(ji,jj,jk) = ua(ji,jj,jk) + ( ziut(ji+1,jj) - ziut(ji,jj ) &278 & + zjuf(ji ,jj) - zjuf(ji,jj-1) ) * r1_e1e2u(ji,jj) / e3u _n(ji,jj,jk)279 va(ji,jj,jk) = va(ji,jj,jk) + ( zivf(ji,jj ) - zivf(ji-1,jj) &280 & + zjvt(ji,jj+1) - zjvt(ji,jj ) ) * r1_e1e2v(ji,jj) / e3v _n(ji,jj,jk)279 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + ( ziut(ji+1,jj) - ziut(ji,jj ) & 280 & + zjuf(ji ,jj) - zjuf(ji,jj-1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 281 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + ( zivf(ji,jj ) - zivf(ji-1,jj) & 282 & + zjvt(ji,jj+1) - zjvt(ji,jj ) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 281 283 END DO 282 284 END DO … … 286 288 287 289 ! print sum trends (used for debugging) 288 IF(ln_ctl) CALL prt_ctl( tab3d_1= ua, clinfo1=' ldfh - Ua: ', mask1=umask, &289 & tab3d_2= va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )290 IF(ln_ctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' ldfh - Ua: ', mask1=umask, & 291 & tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 290 292 291 293 … … 306 308 DO ji = 2, jpi 307 309 ! i-gradient of u at jj 308 zdiu (ji,jk) = tmask(ji,jj ,jk) * ( ub(ji,jj ,jk) - ub(ji-1,jj ,jk) )310 zdiu (ji,jk) = tmask(ji,jj ,jk) * ( puu(ji,jj ,jk,Kbb) - puu(ji-1,jj ,jk,Kbb) ) 309 311 ! j-gradient of u and v at jj 310 zdju (ji,jk) = fmask(ji,jj ,jk) * ( ub(ji,jj+1,jk) - ub(ji ,jj ,jk) )311 zdjv (ji,jk) = tmask(ji,jj ,jk) * ( vb(ji,jj ,jk) - vb(ji ,jj-1,jk) )312 zdju (ji,jk) = fmask(ji,jj ,jk) * ( puu(ji,jj+1,jk,Kbb) - puu(ji ,jj ,jk,Kbb) ) 313 zdjv (ji,jk) = tmask(ji,jj ,jk) * ( pvv(ji,jj ,jk,Kbb) - pvv(ji ,jj-1,jk,Kbb) ) 312 314 ! j-gradient of u and v at jj+1 313 zdj1u(ji,jk) = fmask(ji,jj-1,jk) * ( ub(ji,jj ,jk) - ub(ji ,jj-1,jk) )314 zdj1v(ji,jk) = tmask(ji,jj+1,jk) * ( vb(ji,jj+1,jk) - vb(ji ,jj ,jk) )315 zdj1u(ji,jk) = fmask(ji,jj-1,jk) * ( puu(ji,jj ,jk,Kbb) - puu(ji ,jj-1,jk,Kbb) ) 316 zdj1v(ji,jk) = tmask(ji,jj+1,jk) * ( pvv(ji,jj+1,jk,Kbb) - pvv(ji ,jj ,jk,Kbb) ) 315 317 END DO 316 318 END DO … … 318 320 DO ji = 1, jpim1 319 321 ! i-gradient of v at jj 320 zdiv (ji,jk) = fmask(ji,jj ,jk) * ( vb(ji+1,jj,jk) - vb(ji ,jj ,jk) )322 zdiv (ji,jk) = fmask(ji,jj ,jk) * ( pvv(ji+1,jj,jk,Kbb) - pvv(ji ,jj ,jk,Kbb) ) 321 323 END DO 322 324 END DO … … 391 393 DO jk = 1, jpkm1 392 394 DO ji = 2, jpim1 393 ua(ji,jj,jk) = ua(ji,jj,jk) + ( zfuw(ji,jk) - zfuw(ji,jk+1) ) * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk)394 va(ji,jj,jk) = va(ji,jj,jk) + ( zfvw(ji,jk) - zfvw(ji,jk+1) ) * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk)395 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + ( zfuw(ji,jk) - zfuw(ji,jk+1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 396 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + ( zfvw(ji,jk) - zfvw(ji,jk+1) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 395 397 END DO 396 398 END DO
Note: See TracChangeset
for help on using the changeset viewer.