Changeset 13899 for NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA
- Timestamp:
- 2020-11-27T17:26:33+01:00 (4 years ago)
- Location:
- NEMO/branches/2020/tickets_icb_1900
- Files:
-
- 22 edited
-
. (modified) (1 prop)
-
src/OCE/TRA/eosbn2.F90 (modified) (16 diffs)
-
src/OCE/TRA/traadv_cen.F90 (modified) (9 diffs)
-
src/OCE/TRA/traadv_fct.F90 (modified) (27 diffs)
-
src/OCE/TRA/traadv_mus.F90 (modified) (8 diffs)
-
src/OCE/TRA/traadv_qck.F90 (modified) (13 diffs)
-
src/OCE/TRA/traadv_ubs.F90 (modified) (13 diffs)
-
src/OCE/TRA/traatf.F90 (modified) (2 diffs)
-
src/OCE/TRA/traatf_qco.F90 (modified) (2 diffs)
-
src/OCE/TRA/trabbc.F90 (modified) (1 diff)
-
src/OCE/TRA/trabbl.F90 (modified) (8 diffs)
-
src/OCE/TRA/tradmp.F90 (modified) (4 diffs)
-
src/OCE/TRA/traisf.F90 (modified) (1 diff)
-
src/OCE/TRA/traldf_iso.F90 (modified) (12 diffs)
-
src/OCE/TRA/traldf_lap_blp.F90 (modified) (3 diffs)
-
src/OCE/TRA/traldf_triad.F90 (modified) (13 diffs)
-
src/OCE/TRA/tramle.F90 (modified) (10 diffs)
-
src/OCE/TRA/tranpc.F90 (modified) (1 diff)
-
src/OCE/TRA/traqsr.F90 (modified) (12 diffs)
-
src/OCE/TRA/trasbc.F90 (modified) (6 diffs)
-
src/OCE/TRA/trazdf.F90 (modified) (6 diffs)
-
src/OCE/TRA/zpshde.F90 (modified) (9 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/tickets_icb_1900
- Property svn:externals
-
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/eosbn2.F90
r13237 r13899 238 238 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 239 239 ! 240 DO_3D _11_11(1, jpkm1 )240 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 241 241 ! 242 242 zh = pdep(ji,jj,jk) * r1_Z0 ! depth … … 274 274 CASE( np_seos ) !== simplified EOS ==! 275 275 ! 276 DO_3D _11_11(1, jpkm1 )276 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 277 277 zt = pts (ji,jj,jk,jp_tem) - 10._wp 278 278 zs = pts (ji,jj,jk,jp_sal) - 35._wp … … 338 338 END DO 339 339 ! 340 DO_3D _11_11(1, jpkm1 )340 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 341 341 ! 342 342 ! compute density (2*nn_sto_eos) times: … … 388 388 ! Non-stochastic equation of state 389 389 ELSE 390 DO_3D _11_11(1, jpkm1 )390 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 391 391 ! 392 392 zh = pdep(ji,jj,jk) * r1_Z0 ! depth … … 426 426 CASE( np_seos ) !== simplified EOS ==! 427 427 ! 428 DO_3D _11_11(1, jpkm1 )428 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 429 429 zt = pts (ji,jj,jk,jp_tem) - 10._wp 430 430 zs = pts (ji,jj,jk,jp_sal) - 35._wp … … 480 480 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 481 481 ! 482 DO_2D _11_11482 DO_2D( 1, 1, 1, 1 ) 483 483 ! 484 484 zh = pdep(ji,jj) * r1_Z0 ! depth … … 515 515 CASE( np_seos ) !== simplified EOS ==! 516 516 ! 517 DO_2D _11_11517 DO_2D( 1, 1, 1, 1 ) 518 518 ! 519 519 zt = pts (ji,jj,jp_tem) - 10._wp … … 563 563 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 564 564 ! 565 DO_3D _11_11(1, jpkm1 )565 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 566 566 ! 567 567 zh = gdept(ji,jj,jk,Kmm) * r1_Z0 ! depth … … 616 616 CASE( np_seos ) !== simplified EOS ==! 617 617 ! 618 DO_3D _11_11(1, jpkm1 )618 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 619 619 zt = pts (ji,jj,jk,jp_tem) - 10._wp ! pot. temperature anomaly (t-T0) 620 620 zs = pts (ji,jj,jk,jp_sal) - 35._wp ! abs. salinity anomaly (s-S0) … … 670 670 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 671 671 ! 672 DO_2D _11_11672 DO_2D( 1, 1, 1, 1 ) 673 673 ! 674 674 zh = pdep(ji,jj) * r1_Z0 ! depth … … 723 723 CASE( np_seos ) !== simplified EOS ==! 724 724 ! 725 DO_2D _11_11725 DO_2D( 1, 1, 1, 1 ) 726 726 ! 727 727 zt = pts (ji,jj,jp_tem) - 10._wp ! pot. temperature anomaly (t-T0) … … 873 873 IF( ln_timing ) CALL timing_start('bn2') 874 874 ! 875 DO_3D _11_11( 2, jpkm1 )875 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) ! interior points only (2=< jk =< jpkm1 ); surface and bottom value set to zero one for all in istate.F90 876 876 zrw = ( gdepw(ji,jj,jk ,Kmm) - gdept(ji,jj,jk,Kmm) ) & 877 877 & / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) ) … … 921 921 z1_T0 = 1._wp/40._wp 922 922 ! 923 DO_2D _11_11923 DO_2D( 1, 1, 1, 1 ) 924 924 ! 925 925 zt = ctmp (ji,jj) * z1_T0 … … 974 974 ! 975 975 z1_S0 = 1._wp / 35.16504_wp 976 DO_2D _11_11976 DO_2D( 1, 1, 1, 1 ) 977 977 zs= SQRT( ABS( psal(ji,jj) ) * z1_S0 ) ! square root salinity 978 978 ptf(ji,jj) = ((((1.46873e-03_wp*zs-9.64972e-03_wp)*zs+2.28348e-02_wp)*zs & … … 1081 1081 CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! 1082 1082 ! 1083 DO_3D _11_11(1, jpkm1 )1083 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 1084 1084 ! 1085 1085 zh = gdept(ji,jj,jk,Kmm) * r1_Z0 ! depth … … 1140 1140 CASE( np_seos ) !== Vallis (2006) simplified EOS ==! 1141 1141 ! 1142 DO_3D _11_11(1, jpkm1 )1142 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 1143 1143 zt = pts(ji,jj,jk,jp_tem) - 10._wp ! temperature anomaly (t-T0) 1144 1144 zs = pts (ji,jj,jk,jp_sal) - 35._wp ! abs. salinity anomaly (s-S0) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traadv_cen.F90
r13237 r13899 104 104 ! 105 105 CASE( 2 ) !* 2nd order centered 106 DO_3D _10_10(1, jpkm1 )106 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 107 107 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ) 108 108 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) ) … … 112 112 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 113 113 ztv(:,:,jpk) = 0._wp 114 DO_3D _00_00( 1, jpkm1 )114 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! masked gradient 115 115 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 116 116 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) … … 118 118 CALL lbc_lnk_multi( 'traadv_cen', ztu, 'U', -1.0_wp , ztv, 'V', -1.0_wp ) ! Lateral boundary cond. 119 119 ! 120 DO_3D _00_10( 1, jpkm1 )120 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! Horizontal advective fluxes 121 121 zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ! C2 interpolation of T at u- & v-points (x2) 122 122 zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) … … 128 128 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * zC4t_v 129 129 END_3D 130 CALL lbc_lnk_multi( 'traadv_cen', zwx, 'U', -1. , zwy, 'V', -1. ) 130 131 ! 131 132 CASE DEFAULT 132 CALL ctl_stop( 'traadv_ fct: wrong value for nn_fct' )133 CALL ctl_stop( 'traadv_cen: wrong value for nn_cen' ) 133 134 END SELECT 134 135 ! … … 136 137 ! 137 138 CASE( 2 ) !* 2nd order centered 138 DO_3D _00_00(2, jpk )139 DO_3D( 0, 0, 0, 0, 2, jpk ) 139 140 zwz(ji,jj,jk) = 0.5 * pW(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) * wmask(ji,jj,jk) 140 141 END_3D … … 142 143 CASE( 4 ) !* 4th order compact 143 144 CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw ) ! ztw = interpolated value of T at w-point 144 DO_3D _00_00(2, jpkm1 )145 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 145 146 zwz(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 146 147 END_3D … … 150 151 IF( ln_linssh ) THEN !* top value (linear free surf. only as zwz is multiplied by wmask) 151 152 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 152 DO_2D _11_11153 DO_2D( 1, 1, 1, 1 ) 153 154 zwz(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm) 154 155 END_2D … … 158 159 ENDIF 159 160 ! 160 DO_3D _00_00( 1, jpkm1 )161 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !-- Divergence of advective fluxes --! 161 162 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) & 162 163 & - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & … … 165 166 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 166 167 END_3D 167 ! ! trend diagnostics168 ! ! trend diagnostics 168 169 IF( l_trd ) THEN 169 170 CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) ) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traadv_fct.F90
r13237 r13899 139 139 IF( ll_zAimp ) THEN 140 140 ALLOCATE(zwdia(jpi,jpj,jpk), zwinf(jpi,jpj,jpk),zwsup(jpi,jpj,jpk)) 141 DO_3D _00_00(1, jpkm1 )141 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 142 142 zwdia(ji,jj,jk) = 1._wp + p2dt * ( MAX( wi(ji,jj,jk) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) & 143 143 & / e3t(ji,jj,jk,Krhs) … … 151 151 ! !== upstream advection with initial mass fluxes & intermediate update ==! 152 152 ! !* upstream tracer flux in the i and j direction 153 DO_3D _10_10(1, jpkm1 )153 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 154 154 ! upstream scheme 155 155 zfp_ui = pU(ji,jj,jk) + ABS( pU(ji,jj,jk) ) … … 160 160 zwy(ji,jj,jk) = 0.5 * ( zfp_vj * pt(ji,jj,jk,jn,Kbb) + zfm_vj * pt(ji ,jj+1,jk,jn,Kbb) ) 161 161 END_3D 162 ! !* upstream tracer flux in the k direction *!163 DO_3D _11_11( 2, jpkm1)162 ! !* upstream tracer flux in the k direction *! 163 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) ! Interior value ( multiplied by wmask) 164 164 zfp_wk = pW(ji,jj,jk) + ABS( pW(ji,jj,jk) ) 165 165 zfm_wk = pW(ji,jj,jk) - ABS( pW(ji,jj,jk) ) 166 166 zwz(ji,jj,jk) = 0.5 * ( zfp_wk * pt(ji,jj,jk,jn,Kbb) + zfm_wk * pt(ji,jj,jk-1,jn,Kbb) ) * wmask(ji,jj,jk) 167 167 END_3D 168 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked)169 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface170 DO_2D _11_11168 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as zwz has been w-masked) 169 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 170 DO_2D( 1, 1, 1, 1 ) 171 171 zwz(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) ! linear free surface 172 172 END_2D 173 ELSE ! no cavities: only at the ocean surface 174 zwz(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 173 ELSE ! no cavities: only at the ocean surface 174 DO_2D( 1, 1, 1, 1 ) 175 zwz(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kbb) 176 END_2D 175 177 ENDIF 176 178 ENDIF 177 179 ! 178 DO_3D _00_00( 1, jpkm1 )179 ! ! total intermediate advective trends180 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !* trend and after field with monotonic scheme 181 ! ! total intermediate advective trends 180 182 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 181 183 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 182 184 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) 183 ! ! update and guess with monotonic sheme185 ! ! update and guess with monotonic sheme 184 186 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra & 185 187 & / e3t(ji,jj,jk,Kmm ) * tmask(ji,jj,jk) … … 192 194 ! 193 195 ztw(:,:,1) = 0._wp ; ztw(:,:,jpk) = 0._wp ; 194 DO_3D _00_00( 2, jpkm1)196 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! Interior value ( multiplied by wmask) 195 197 zfp_wk = wi(ji,jj,jk) + ABS( wi(ji,jj,jk) ) 196 198 zfm_wk = wi(ji,jj,jk) - ABS( wi(ji,jj,jk) ) … … 198 200 zwz(ji,jj,jk) = zwz(ji,jj,jk) + ztw(ji,jj,jk) ! update vertical fluxes 199 201 END_3D 200 DO_3D _00_00(1, jpkm1 )202 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 201 203 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) & 202 204 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) … … 216 218 ! 217 219 CASE( 2 ) !- 2nd order centered 218 DO_3D _10_10(1, jpkm1 )220 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 219 221 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj,jk,jn,Kmm) ) - zwx(ji,jj,jk) 220 222 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj+1,jk,jn,Kmm) ) - zwy(ji,jj,jk) … … 225 227 zltv(:,:,jpk) = 0._wp 226 228 DO jk = 1, jpkm1 ! Laplacian 227 DO_2D _10_10229 DO_2D( 1, 0, 1, 0 ) ! 1st derivative (gradient) 228 230 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 229 231 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 230 232 END_2D 231 DO_2D _00_00233 DO_2D( 0, 0, 0, 0 ) ! 2nd derivative * 1/ 6 232 234 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) + ztu(ji-1,jj,jk) ) * r1_6 233 235 zltv(ji,jj,jk) = ( ztv(ji,jj,jk) + ztv(ji,jj-1,jk) ) * r1_6 … … 236 238 CALL lbc_lnk_multi( 'traadv_fct', zltu, 'T', 1.0_wp , zltv, 'T', 1.0_wp ) ! Lateral boundary cond. (unchanged sgn) 237 239 ! 238 DO_3D _10_10( 1, jpkm1 )240 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) ! Horizontal advective fluxes 239 241 zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ! 2 x C2 interpolation of T at u- & v-points 240 242 zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) 241 ! ! C4 minus upstream advective fluxes243 ! ! C4 minus upstream advective fluxes 242 244 zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 243 245 zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) … … 247 249 ztu(:,:,jpk) = 0._wp ! Bottom value : flux set to zero 248 250 ztv(:,:,jpk) = 0._wp 249 DO_3D _10_10( 1, jpkm1)251 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) ! 1st derivative (gradient) 250 252 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 251 253 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) … … 253 255 CALL lbc_lnk_multi( 'traadv_fct', ztu, 'U', -1.0_wp , ztv, 'V', -1.0_wp ) ! Lateral boundary cond. (unchanged sgn) 254 256 ! 255 DO_3D _00_00( 1, jpkm1 )257 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! Horizontal advective fluxes 256 258 zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj ,jk,jn,Kmm) ! 2 x C2 interpolation of T at u- & v-points (x2) 257 259 zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji ,jj+1,jk,jn,Kmm) … … 269 271 ! 270 272 CASE( 2 ) !- 2nd order centered 271 DO_3D _00_00(2, jpkm1 )273 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 272 274 zwz(ji,jj,jk) = ( pW(ji,jj,jk) * 0.5_wp * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) & 273 275 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) … … 276 278 CASE( 4 ) !- 4th order COMPACT 277 279 CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw ) ! zwt = COMPACT interpolation of T at w-point 278 DO_3D _00_00(2, jpkm1 )280 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 279 281 zwz(ji,jj,jk) = ( pW(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 280 282 END_3D … … 286 288 ! 287 289 IF ( ll_zAimp ) THEN 288 DO_3D _00_00( 1, jpkm1 )289 ! ! total intermediate advective trends290 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !* trend and after field with monotonic scheme 291 ! ! total intermediate advective trends 290 292 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 291 293 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & … … 296 298 CALL tridia_solver( zwdia, zwsup, zwinf, ztw, ztw , 0 ) 297 299 ! 298 DO_3D _00_00( 2, jpkm1)300 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! Interior value ( multiplied by wmask) 299 301 zfp_wk = wi(ji,jj,jk) + ABS( wi(ji,jj,jk) ) 300 302 zfm_wk = wi(ji,jj,jk) - ABS( wi(ji,jj,jk) ) … … 311 313 ! !== final trend with corrected fluxes ==! 312 314 ! 313 DO_3D _00_00(1, jpkm1 )315 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 314 316 ztra = - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 315 317 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & … … 322 324 ! 323 325 ztw(:,:,1) = 0._wp ; ztw(:,:,jpk) = 0._wp 324 DO_3D _00_00( 2, jpkm1)326 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! Interior value ( multiplied by wmask) 325 327 zfp_wk = wi(ji,jj,jk) + ABS( wi(ji,jj,jk) ) 326 328 zfm_wk = wi(ji,jj,jk) - ABS( wi(ji,jj,jk) ) … … 328 330 zwz(ji,jj,jk) = zwz(ji,jj,jk) + ztw(ji,jj,jk) ! Update vertical fluxes for trend diagnostic 329 331 END_3D 330 DO_3D _00_00(1, jpkm1 )332 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 331 333 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji ,jj ,jk+1) ) & 332 334 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) … … 407 409 DO jk = 1, jpkm1 408 410 ikm1 = MAX(jk-1,1) 409 DO_2D _00_00411 DO_2D( 0, 0, 0, 0 ) 410 412 411 413 ! search maximum in neighbourhood … … 441 443 ! 3. monotonic flux in the i & j direction (paa & pbb) 442 444 ! ---------------------------------------- 443 DO_3D _00_00(1, jpkm1 )445 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 444 446 zau = MIN( 1._wp, zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) ) 445 447 zbu = MIN( 1._wp, zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) ) … … 452 454 pbb(ji,jj,jk) = pbb(ji,jj,jk) * ( zcv * zav + ( 1._wp - zcv) * zbv ) 453 455 454 ! monotonic flux in the k direction, i.e. pcc455 ! -------------------------------------------456 ! monotonic flux in the k direction, i.e. pcc 457 ! ------------------------------------------- 456 458 za = MIN( 1., zbetdo(ji,jj,jk+1), zbetup(ji,jj,jk) ) 457 459 zb = MIN( 1., zbetup(ji,jj,jk+1), zbetdo(ji,jj,jk) ) … … 479 481 !!---------------------------------------------------------------------- 480 482 481 DO_3D _11_11( 3, jpkm1 )483 DO_3D( 1, 1, 1, 1, 3, jpkm1 ) !== build the three diagonal matrix ==! 482 484 zwd (ji,jj,jk) = 4._wp 483 485 zwi (ji,jj,jk) = 1._wp … … 493 495 END_3D 494 496 ! 495 jk = 2 496 DO_2D _11_11497 jk = 2 ! Switch to second order centered at top 498 DO_2D( 1, 1, 1, 1 ) 497 499 zwd (ji,jj,jk) = 1._wp 498 500 zwi (ji,jj,jk) = 0._wp … … 502 504 ! 503 505 ! !== tridiagonal solve ==! 504 DO_2D _11_11506 DO_2D( 1, 1, 1, 1 ) ! first recurrence 505 507 zwt(ji,jj,2) = zwd(ji,jj,2) 506 508 END_2D 507 DO_3D _11_11(3, jpkm1 )509 DO_3D( 1, 1, 1, 1, 3, jpkm1 ) 508 510 zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) /zwt(ji,jj,jk-1) 509 511 END_3D 510 512 ! 511 DO_2D _11_11513 DO_2D( 1, 1, 1, 1 ) ! second recurrence: Zk = Yk - Ik / Tk-1 Zk-1 512 514 pt_out(ji,jj,2) = zwrm(ji,jj,2) 513 515 END_2D 514 DO_3D _11_11(3, jpkm1 )516 DO_3D( 1, 1, 1, 1, 3, jpkm1 ) 515 517 pt_out(ji,jj,jk) = zwrm(ji,jj,jk) - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) *pt_out(ji,jj,jk-1) 516 518 END_3D 517 519 518 DO_2D _11_11520 DO_2D( 1, 1, 1, 1 ) ! third recurrence: Xk = (Zk - Sk Xk+1 ) / Tk 519 521 pt_out(ji,jj,jpkm1) = pt_out(ji,jj,jpkm1) / zwt(ji,jj,jpkm1) 520 522 END_2D 521 DO_3DS _11_11(jpk-2, 2, -1 )523 DO_3DS( 1, 1, 1, 1, jpk-2, 2, -1 ) 522 524 pt_out(ji,jj,jk) = ( pt_out(ji,jj,jk) - zws(ji,jj,jk) * pt_out(ji,jj,jk+1) ) / zwt(ji,jj,jk) 523 525 END_3D … … 544 546 ! !== build the three diagonal matrix & the RHS ==! 545 547 ! 546 DO_3D _00_00( 3, jpkm1)548 DO_3D( 0, 0, 0, 0, 3, jpkm1 ) ! interior (from jk=3 to jpk-1) 547 549 zwd (ji,jj,jk) = 3._wp * wmask(ji,jj,jk) + 1._wp ! diagonal 548 550 zwi (ji,jj,jk) = wmask(ji,jj,jk) ! lower diagonal … … 563 565 END IF 564 566 ! 565 DO_2D _00_00567 DO_2D( 0, 0, 0, 0 ) ! 2nd order centered at top & bottom 566 568 ikt = mikt(ji,jj) + 1 ! w-point below the 1st wet point 567 569 ikb = MAX(mbkt(ji,jj), 2) ! - above the last wet point … … 580 582 ! !== tridiagonal solver ==! 581 583 ! 582 DO_2D _00_00584 DO_2D( 0, 0, 0, 0 ) !* 1st recurrence: Tk = Dk - Ik Sk-1 / Tk-1 583 585 zwt(ji,jj,2) = zwd(ji,jj,2) 584 586 END_2D 585 DO_3D _00_00(3, jpkm1 )587 DO_3D( 0, 0, 0, 0, 3, jpkm1 ) 586 588 zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) /zwt(ji,jj,jk-1) 587 589 END_3D 588 590 ! 589 DO_2D _00_00591 DO_2D( 0, 0, 0, 0 ) !* 2nd recurrence: Zk = Yk - Ik / Tk-1 Zk-1 590 592 pt_out(ji,jj,2) = zwrm(ji,jj,2) 591 593 END_2D 592 DO_3D _00_00(3, jpkm1 )594 DO_3D( 0, 0, 0, 0, 3, jpkm1 ) 593 595 pt_out(ji,jj,jk) = zwrm(ji,jj,jk) - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) *pt_out(ji,jj,jk-1) 594 596 END_3D 595 597 596 DO_2D _00_00598 DO_2D( 0, 0, 0, 0 ) !* 3d recurrence: Xk = (Zk - Sk Xk+1 ) / Tk 597 599 pt_out(ji,jj,jpkm1) = pt_out(ji,jj,jpkm1) / zwt(ji,jj,jpkm1) 598 600 END_2D 599 DO_3DS _00_00(jpk-2, 2, -1 )601 DO_3DS( 0, 0, 0, 0, jpk-2, 2, -1 ) 600 602 pt_out(ji,jj,jk) = ( pt_out(ji,jj,jk) - zws(ji,jj,jk) * pt_out(ji,jj,jk+1) ) / zwt(ji,jj,jk) 601 603 END_3D … … 636 638 kstart = 1 + klev 637 639 ! 638 DO_2D _00_00640 DO_2D( 0, 0, 0, 0 ) !* 1st recurrence: Tk = Dk - Ik Sk-1 / Tk-1 639 641 zwt(ji,jj,kstart) = pD(ji,jj,kstart) 640 642 END_2D 641 DO_3D _00_00(kstart+1, jpkm1 )643 DO_3D( 0, 0, 0, 0, kstart+1, jpkm1 ) 642 644 zwt(ji,jj,jk) = pD(ji,jj,jk) - pL(ji,jj,jk) * pU(ji,jj,jk-1) /zwt(ji,jj,jk-1) 643 645 END_3D 644 646 ! 645 DO_2D _00_00647 DO_2D( 0, 0, 0, 0 ) !* 2nd recurrence: Zk = Yk - Ik / Tk-1 Zk-1 646 648 pt_out(ji,jj,kstart) = pRHS(ji,jj,kstart) 647 649 END_2D 648 DO_3D _00_00(kstart+1, jpkm1 )650 DO_3D( 0, 0, 0, 0, kstart+1, jpkm1 ) 649 651 pt_out(ji,jj,jk) = pRHS(ji,jj,jk) - pL(ji,jj,jk) / zwt(ji,jj,jk-1) *pt_out(ji,jj,jk-1) 650 652 END_3D 651 653 652 DO_2D _00_00654 DO_2D( 0, 0, 0, 0 ) !* 3d recurrence: Xk = (Zk - Sk Xk+1 ) / Tk 653 655 pt_out(ji,jj,jpkm1) = pt_out(ji,jj,jpkm1) / zwt(ji,jj,jpkm1) 654 656 END_2D 655 DO_3DS _00_00(jpk-2, kstart, -1 )657 DO_3DS( 0, 0, 0, 0, jpk-2, kstart, -1 ) 656 658 pt_out(ji,jj,jk) = ( pt_out(ji,jj,jk) - pU(ji,jj,jk) * pt_out(ji,jj,jk+1) ) / zwt(ji,jj,jk) 657 659 END_3D -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traadv_mus.F90
r13237 r13899 132 132 zwx(:,:,jpk) = 0._wp ! bottom values 133 133 zwy(:,:,jpk) = 0._wp 134 DO_3D _10_10(1, jpkm1 )134 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 135 135 zwx(ji,jj,jk) = umask(ji,jj,jk) * ( pt(ji+1,jj,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 136 136 zwy(ji,jj,jk) = vmask(ji,jj,jk) * ( pt(ji,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) … … 141 141 zslpx(:,:,jpk) = 0._wp ! bottom values 142 142 zslpy(:,:,jpk) = 0._wp 143 DO_3D _01_01(1, jpkm1 )143 DO_3D( 0, 1, 0, 1, 1, jpkm1 ) 144 144 zslpx(ji,jj,jk) = ( zwx(ji,jj,jk) + zwx(ji-1,jj ,jk) ) & 145 145 & * ( 0.25 + SIGN( 0.25_wp, zwx(ji,jj,jk) * zwx(ji-1,jj ,jk) ) ) … … 148 148 END_3D 149 149 ! 150 DO_3D _01_01( 1, jpkm1 )150 DO_3D( 0, 1, 0, 1, 1, jpkm1 ) !-- Slopes limitation 151 151 zslpx(ji,jj,jk) = SIGN( 1.0_wp, zslpx(ji,jj,jk) ) * MIN( ABS( zslpx(ji ,jj,jk) ), & 152 152 & 2.*ABS( zwx (ji-1,jj,jk) ), & … … 157 157 END_3D 158 158 ! 159 DO_3D _00_00( 1, jpkm1 )159 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !-- MUSCL horizontal advective fluxes 160 160 ! MUSCL fluxes 161 161 z0u = SIGN( 0.5_wp, pU(ji,jj,jk) ) … … 175 175 CALL lbc_lnk_multi( 'traadv_mus', zwx, 'U', -1.0_wp , zwy, 'V', -1.0_wp ) ! lateral boundary conditions (changed sign) 176 176 ! 177 DO_3D _00_00( 1, jpkm1 )177 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !-- Tracer advective trend 178 178 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 179 179 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) ) & … … 200 200 ! !-- Slopes of tracer 201 201 zslpx(:,:,1) = 0._wp ! surface values 202 DO_3D _11_11(2, jpkm1 )202 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) 203 203 zslpx(ji,jj,jk) = ( zwx(ji,jj,jk) + zwx(ji,jj,jk+1) ) & 204 204 & * ( 0.25 + SIGN( 0.25_wp, zwx(ji,jj,jk) * zwx(ji,jj,jk+1) ) ) 205 205 END_3D 206 DO_3D _11_11( 2, jpkm1 )206 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) !-- Slopes limitation 207 207 zslpx(ji,jj,jk) = SIGN( 1.0_wp, zslpx(ji,jj,jk) ) * MIN( ABS( zslpx(ji,jj,jk ) ), & 208 208 & 2.*ABS( zwx (ji,jj,jk+1) ), & 209 209 & 2.*ABS( zwx (ji,jj,jk ) ) ) 210 210 END_3D 211 DO_3D _00_00( 1, jpk-2 )211 DO_3D( 0, 0, 0, 0, 1, jpk-2 ) !-- vertical advective flux 212 212 z0w = SIGN( 0.5_wp, pW(ji,jj,jk+1) ) 213 213 zalpha = 0.5 + z0w … … 219 219 IF( ln_linssh ) THEN ! top values, linear free surface only 220 220 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 221 DO_2D _11_11221 DO_2D( 1, 1, 1, 1 ) 222 222 zwx(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) 223 223 END_2D … … 227 227 ENDIF 228 228 ! 229 DO_3D _00_00( 1, jpkm1 )229 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !-- vertical advective trend 230 230 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) & 231 231 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traadv_qck.F90
r13237 r13899 142 142 ! 143 143 !!gm why not using a SHIFT instruction... 144 DO_3D _00_00( 1, jpkm1 )144 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !--- Computation of the ustream and downstream value of the tracer and the mask 145 145 zfc(ji,jj,jk) = pt(ji-1,jj,jk,jn,Kbb) ! Upstream in the x-direction for the tracer 146 146 zfd(ji,jj,jk) = pt(ji+1,jj,jk,jn,Kbb) ! Downstream in the x-direction for the tracer … … 151 151 ! Horizontal advective fluxes 152 152 ! --------------------------- 153 DO_3D _00_00(1, jpkm1 )153 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 154 154 zdir = 0.5 + SIGN( 0.5_wp, pU(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 155 155 zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji+1,jj,jk) ! FU in the x-direction for T 156 156 END_3D 157 157 ! 158 DO_3D _00_00(1, jpkm1 )158 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 159 159 zdir = 0.5 + SIGN( 0.5_wp, pU(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 160 160 zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) … … 170 170 ! 171 171 ! Mask at the T-points in the x-direction (mask=0 or mask=1) 172 DO_3D _00_00(1, jpkm1 )172 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 173 173 zfu(ji,jj,jk) = tmask(ji-1,jj,jk) + tmask(ji,jj,jk) + tmask(ji+1,jj,jk) - 2. 174 174 END_3D … … 179 179 DO jk = 1, jpkm1 180 180 ! 181 DO_2D _00_00181 DO_2D( 0, 0, 0, 0 ) 182 182 zdir = 0.5 + SIGN( 0.5_wp, pU(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 183 183 !--- If the second ustream point is a land point … … 192 192 ! 193 193 ! Computation of the trend 194 DO_3D _00_00(1, jpkm1 )194 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 195 195 zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 196 196 ! horizontal advective trends … … 233 233 ! 234 234 !--- Computation of the ustream and downstream value of the tracer and the mask 235 DO_2D _00_00235 DO_2D( 0, 0, 0, 0 ) 236 236 ! Upstream in the x-direction for the tracer 237 237 zfc(ji,jj,jk) = pt(ji,jj-1,jk,jn,Kbb) … … 247 247 ! --------------------------- 248 248 ! 249 DO_3D _00_00(1, jpkm1 )249 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 250 250 zdir = 0.5 + SIGN( 0.5_wp, pV(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 251 251 zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji,jj+1,jk) ! FU in the x-direction for T 252 252 END_3D 253 253 ! 254 DO_3D _00_00(1, jpkm1 )254 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 255 255 zdir = 0.5 + SIGN( 0.5_wp, pV(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 256 256 zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) … … 267 267 ! 268 268 ! Mask at the T-points in the x-direction (mask=0 or mask=1) 269 DO_3D _00_00(1, jpkm1 )269 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 270 270 zfu(ji,jj,jk) = tmask(ji,jj-1,jk) + tmask(ji,jj,jk) + tmask(ji,jj+1,jk) - 2. 271 271 END_3D … … 275 275 DO jk = 1, jpkm1 276 276 ! 277 DO_2D _00_00277 DO_2D( 0, 0, 0, 0 ) 278 278 zdir = 0.5 + SIGN( 0.5_wp, pV(ji,jj,jk) ) ! if pU > 0 : zdir = 1 otherwise zdir = 0 279 279 !--- If the second ustream point is a land point … … 288 288 ! 289 289 ! Computation of the trend 290 DO_3D _00_00(1, jpkm1 )290 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 291 291 zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 292 292 ! horizontal advective trends … … 327 327 ! ! =========== 328 328 ! 329 DO_3D _00_00( 2, jpkm1)329 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !* Interior point (w-masked 2nd order centered flux) 330 330 zwz(ji,jj,jk) = 0.5 * pW(ji,jj,jk) * ( pt(ji,jj,jk-1,jn,Kmm) + pt(ji,jj,jk,jn,Kmm) ) * wmask(ji,jj,jk) 331 331 END_3D 332 332 IF( ln_linssh ) THEN !* top value (only in linear free surf. as zwz is multiplied by wmask) 333 333 IF( ln_isfcav ) THEN ! ice-shelf cavities (top of the ocean) 334 DO_2D _11_11334 DO_2D( 1, 1, 1, 1 ) 335 335 zwz(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm) ! linear free surface 336 336 END_2D … … 340 340 ENDIF 341 341 ! 342 DO_3D _00_00( 1, jpkm1 )342 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Tracer flux divergence added to the general trend ==! 343 343 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) & 344 344 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) … … 369 369 !---------------------------------------------------------------------- 370 370 ! 371 DO_3D _11_11(1, jpkm1 )371 DO_3D( 1, 1, 1, 1, 1, jpkm1 ) 372 372 zc = puc(ji,jj,jk) ! Courant number 373 373 zcurv = pfd(ji,jj,jk) + pfu(ji,jj,jk) - 2. * pfc(ji,jj,jk) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traadv_ubs.F90
r13237 r13899 124 124 ! ! =========== 125 125 ! 126 DO jk = 1, jpkm1 !== horizontal laplacian of before tracer ==!127 DO_2D _10_10126 DO jk = 1, jpkm1 !== horizontal laplacian of before tracer ==! 127 DO_2D( 1, 0, 1, 0 ) ! First derivative (masked gradient) 128 128 zeeu = e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) * umask(ji,jj,jk) 129 129 zeev = e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * vmask(ji,jj,jk) … … 131 131 ztv(ji,jj,jk) = zeev * ( pt(ji ,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 132 132 END_2D 133 DO_2D _00_00133 DO_2D( 0, 0, 0, 0 ) ! Second derivative (divergence) 134 134 zcoef = 1._wp / ( 6._wp * e3t(ji,jj,jk,Kmm) ) 135 135 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zcoef … … 140 140 CALL lbc_lnk( 'traadv_ubs', zltu, 'T', 1.0_wp ) ; CALL lbc_lnk( 'traadv_ubs', zltv, 'T', 1.0_wp ) ! Lateral boundary cond. (unchanged sgn) 141 141 ! 142 DO_3D _10_10( 1, jpkm1)143 zfp_ui = pU(ji,jj,jk) + ABS( pU(ji,jj,jk) ) ! upstream transport (x2)142 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) !== Horizontal advective fluxes ==! (UBS) 143 zfp_ui = pU(ji,jj,jk) + ABS( pU(ji,jj,jk) ) ! upstream transport (x2) 144 144 zfm_ui = pU(ji,jj,jk) - ABS( pU(ji,jj,jk) ) 145 145 zfp_vj = pV(ji,jj,jk) + ABS( pV(ji,jj,jk) ) … … 156 156 ! 157 157 DO jk = 1, jpkm1 !== add the horizontal advective trend ==! 158 DO_2D _00_00158 DO_2D( 0, 0, 0, 0 ) 159 159 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) & 160 160 & - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk) & … … 166 166 ! 167 167 zltu(:,:,:) = pt(:,:,:,jn,Krhs) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case 168 ! ! and/or in trend diagnostic (l_trd=T)168 ! ! and/or in trend diagnostic (l_trd=T) 169 169 ! 170 170 IF( l_trd ) THEN ! trend diagnostics … … 187 187 IF( l_trd ) zltv(:,:,:) = pt(:,:,:,jn,Krhs) ! store pt(:,:,:,:,Krhs) if trend diag. 188 188 ! 189 ! !* upstream advection with initial mass fluxes & intermediate update ==!190 DO_3D _11_11(2, jpkm1 )189 ! !* upstream advection with initial mass fluxes & intermediate update ==! 190 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) 191 191 zfp_wk = pW(ji,jj,jk) + ABS( pW(ji,jj,jk) ) 192 192 zfm_wk = pW(ji,jj,jk) - ABS( pW(ji,jj,jk) ) 193 193 ztw(ji,jj,jk) = 0.5_wp * ( zfp_wk * pt(ji,jj,jk,jn,Kbb) + zfm_wk * pt(ji,jj,jk-1,jn,Kbb) ) * wmask(ji,jj,jk) 194 194 END_3D 195 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as ztw has been w-masked)196 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface197 DO_2D _11_11195 IF( ln_linssh ) THEN ! top ocean value (only in linear free surface as ztw has been w-masked) 196 IF( ln_isfcav ) THEN ! top of the ice-shelf cavities and at the ocean surface 197 DO_2D( 1, 1, 1, 1 ) 198 198 ztw(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) ! linear free surface 199 199 END_2D 200 ELSE ! no cavities: only at the ocean surface200 ELSE ! no cavities: only at the ocean surface 201 201 ztw(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kbb) 202 202 ENDIF 203 203 ENDIF 204 204 ! 205 DO_3D _00_00( 1, jpkm1 )205 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !* trend and after field with monotonic scheme 206 206 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) & 207 207 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) … … 212 212 ! 213 213 ! !* anti-diffusive flux : high order minus low order 214 DO_3D _11_11(2, jpkm1 )214 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) 215 215 ztw(ji,jj,jk) = ( 0.5_wp * pW(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) & 216 216 & - ztw(ji,jj,jk) ) * wmask(ji,jj,jk) … … 223 223 CASE( 4 ) ! 4th order COMPACT 224 224 CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw ) ! 4th order compact interpolation of T at w-point 225 DO_3D _00_00(2, jpkm1 )225 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 226 226 ztw(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 227 227 END_3D … … 230 230 END SELECT 231 231 ! 232 DO_3D _00_00( 1, jpkm1 )232 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! final trend with corrected fluxes 233 233 pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) & 234 234 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 235 235 END_3D 236 236 ! 237 IF( l_trd ) THEN ! vertical advective trend diagnostics238 DO_3D _00_00( 1, jpkm1)237 IF( l_trd ) THEN ! vertical advective trend diagnostics 238 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! (compute -w.dk[ptn]= -dk[w.ptn] + ptn.dk[w]) 239 239 zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk) & 240 240 & + pt(ji,jj,jk,jn,Kmm) * ( pW(ji,jj,jk) - pW(ji,jj,jk+1) ) & … … 286 286 DO jk = 1, jpkm1 ! search maximum in neighbourhood 287 287 ikm1 = MAX(jk-1,1) 288 DO_2D _00_00288 DO_2D( 0, 0, 0, 0 ) 289 289 zbetup(ji,jj,jk) = MAX( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), & 290 290 & pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), & … … 298 298 DO jk = 1, jpkm1 ! search minimum in neighbourhood 299 299 ikm1 = MAX(jk-1,1) 300 DO_2D _00_00300 DO_2D( 0, 0, 0, 0 ) 301 301 zbetdo(ji,jj,jk) = MIN( pbef(ji ,jj ,jk ), paft(ji ,jj ,jk ), & 302 302 & pbef(ji ,jj ,ikm1), pbef(ji ,jj ,jk+1), & … … 310 310 ! Positive and negative part of fluxes and beta terms 311 311 ! --------------------------------------------------- 312 DO_3D _00_00(1, jpkm1 )312 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 313 313 ! positive & negative part of the flux 314 314 zpos = MAX( 0., pcc(ji ,jj ,jk+1) ) - MIN( 0., pcc(ji ,jj ,jk ) ) … … 322 322 ! monotonic flux in the k direction, i.e. pcc 323 323 ! ------------------------------------------- 324 DO_3D _00_00(2, jpkm1 )324 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 325 325 za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) ) 326 326 zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) ) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traatf.F90
r13237 r13899 210 210 DO jn = 1, kjpt 211 211 ! 212 DO_3D _00_00(1, jpkm1 )212 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 213 213 ztn = pt(ji,jj,jk,jn,Kmm) 214 214 ztd = pt(ji,jj,jk,jn,Kaa) - 2._wp * ztn + pt(ji,jj,jk,jn,Kbb) ! time laplacian on tracers … … 275 275 zfact2 = zfact1 * r1_rho0 276 276 DO jn = 1, kjpt 277 DO_3D _00_00(1, jpkm1 )277 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 278 278 ze3t_b = e3t(ji,jj,jk,Kbb) 279 279 ze3t_n = e3t(ji,jj,jk,Kmm) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traatf_qco.F90
r13237 r13899 203 203 DO jn = 1, kjpt 204 204 ! 205 DO_3D _00_00(1, jpkm1 )205 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 206 206 ztn = pt(ji,jj,jk,jn,Kmm) 207 207 ztd = pt(ji,jj,jk,jn,Kaa) - 2._wp * ztn + pt(ji,jj,jk,jn,Kbb) ! time laplacian on tracers … … 268 268 zfact2 = zfact1 * r1_rho0 269 269 DO jn = 1, kjpt 270 DO_3D _00_00(1, jpkm1 )270 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 271 271 ze3t_b = e3t(ji,jj,jk,Kbb) 272 272 ze3t_n = e3t(ji,jj,jk,Kmm) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/trabbc.F90
r13237 r13899 91 91 ENDIF 92 92 ! ! Add the geothermal trend on temperature 93 DO_2D _00_0093 DO_2D( 0, 0, 0, 0 ) 94 94 pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs) = pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs) & 95 95 & + qgh_trd0(ji,jj) / e3t(ji,jj,mbkt(ji,jj),Kmm) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/trabbl.F90
r13237 r13899 192 192 DO jn = 1, kjpt ! tracer loop 193 193 ! ! =========== 194 DO_2D _11_11194 DO_2D( 1, 1, 1, 1 ) 195 195 ik = mbkt(ji,jj) ! bottom T-level index 196 196 zptb(ji,jj) = pt(ji,jj,ik,jn) ! bottom before T and S 197 197 END_2D 198 198 ! 199 DO_2D _00_00199 DO_2D( 0, 0, 0, 0 ) ! Compute the trend 200 200 ik = mbkt(ji,jj) ! bottom T-level index 201 201 pt_rhs(ji,jj,ik,jn) = pt_rhs(ji,jj,ik,jn) & … … 343 343 ENDIF 344 344 ! !* bottom variables (T, S, alpha, beta, depth, velocity) 345 DO_2D _11_11345 DO_2D( 1, 1, 1, 1 ) 346 346 ik = mbkt(ji,jj) ! bottom T-level index 347 347 zts (ji,jj,jp_tem) = ts(ji,jj,ik,jp_tem,Kbb) ! bottom before T and S … … 358 358 IF( nn_bbl_ldf == 1 ) THEN ! diffusive bbl ! 359 359 ! !-------------------! 360 DO_2D _10_10360 DO_2D( 1, 0, 1, 0 ) ! (criteria for non zero flux: grad(rho).grad(h) < 0 ) 361 361 ! ! i-direction 362 362 za = zab(ji+1,jj,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at u-point … … 388 388 ! 389 389 CASE( 1 ) != use of upper velocity 390 DO_2D _10_10390 DO_2D( 1, 0, 1, 0 ) ! criteria: grad(rho).grad(h)<0 and grad(rho).grad(h)<0 391 391 ! ! i-direction 392 392 za = zab(ji+1,jj,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at u-point … … 417 417 CASE( 2 ) != bbl velocity = F( delta rho ) 418 418 zgbbl = grav * rn_gambbl 419 DO_2D _10_10419 DO_2D( 1, 0, 1, 0 ) ! criteria: rho_up > rho_down 420 420 ! ! i-direction 421 421 ! down-slope T-point i/k-index (deep) & up-slope T-point i/k-index (shelf) … … 505 505 IF( tra_bbl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'tra_bbl_init : unable to allocate arrays' ) 506 506 ! 507 IF( nn_bbl_adv == 1 ) WRITE(numout,*) ' * Advective BBL using upper velocity' 508 IF( nn_bbl_adv == 2 ) WRITE(numout,*) ' * Advective BBL using velocity = F( delta rho)' 507 IF(lwp) THEN 508 IF( nn_bbl_adv == 1 ) WRITE(numout,*) ' * Advective BBL using upper velocity' 509 IF( nn_bbl_adv == 2 ) WRITE(numout,*) ' * Advective BBL using velocity = F( delta rho)' 510 ENDIF 509 511 ! 510 512 ! !* vertical index of "deep" bottom u- and v-points 511 DO_2D _10_10513 DO_2D( 1, 0, 1, 0 ) ! (the "shelf" bottom k-indices are mbku and mbkv) 512 514 mbku_d(ji,jj) = MAX( mbkt(ji+1,jj ) , mbkt(ji,jj) ) ! >= 1 as mbkt=1 over land 513 515 mbkv_d(ji,jj) = MAX( mbkt(ji ,jj+1) , mbkt(ji,jj) ) … … 520 522 ! !* sign of grad(H) at u- and v-points; zero if grad(H) = 0 521 523 mgrhu(:,:) = 0 ; mgrhv(:,:) = 0 522 DO_2D _10_10524 DO_2D( 1, 0, 1, 0 ) 523 525 IF( gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 524 526 mgrhu(ji,jj) = INT( SIGN( 1.0_wp, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) ) ) … … 530 532 END_2D 531 533 ! 532 DO_2D _10_10534 DO_2D( 1, 0, 1, 0 ) !* bbl thickness at u- (v-) point; minimum of top & bottom e3u_0 (e3v_0) 533 535 e3u_bbl_0(ji,jj) = MIN( e3u_0(ji,jj,mbkt(ji+1,jj )), e3u_0(ji,jj,mbkt(ji,jj)) ) 534 536 e3v_bbl_0(ji,jj) = MIN( e3v_0(ji,jj,mbkt(ji ,jj+1)), e3v_0(ji,jj,mbkt(ji,jj)) ) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/tradmp.F90
r12377 r13899 112 112 CASE( 0 ) !* newtonian damping throughout the water column *! 113 113 DO jn = 1, jpts 114 DO_3D _00_00(1, jpkm1 )114 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 115 115 pts(ji,jj,jk,jn,Krhs) = pts(ji,jj,jk,jn,Krhs) & 116 116 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jn) - pts(ji,jj,jk,jn,Kbb) ) … … 119 119 ! 120 120 CASE ( 1 ) !* no damping in the turbocline (avt > 5 cm2/s) *! 121 DO_3D _00_00(1, jpkm1 )121 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 122 122 IF( avt(ji,jj,jk) <= avt_c ) THEN 123 123 pts(ji,jj,jk,jp_tem,Krhs) = pts(ji,jj,jk,jp_tem,Krhs) & … … 129 129 ! 130 130 CASE ( 2 ) !* no damping in the mixed layer *! 131 DO_3D _00_00(1, jpkm1 )131 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 132 132 IF( gdept(ji,jj,jk,Kmm) >= hmlp (ji,jj) ) THEN 133 133 pts(ji,jj,jk,jp_tem,Krhs) = pts(ji,jj,jk,jp_tem,Krhs) & … … 208 208 ! ! Read in mask from file 209 209 CALL iom_open ( cn_resto, imask) 210 CALL iom_get ( imask, jpdom_auto glo, 'resto', resto )210 CALL iom_get ( imask, jpdom_auto, 'resto', resto ) 211 211 CALL iom_close( imask ) 212 212 ENDIF -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traisf.F90
r13237 r13899 108 108 ! 109 109 ! update pts(:,:,:,:,Krhs) 110 DO_2D _11_11110 DO_2D( 1, 1, 1, 1 ) 111 111 ! 112 112 ikt = ktop(ji,jj) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_iso.F90
r13237 r13899 141 141 IF( kpass == 1 ) THEN !== first pass only ==! 142 142 ! 143 DO_3D _00_00(2, jpkm1 )143 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 144 144 ! 145 145 zmsku = wmask(ji,jj,jk) / MAX( umask(ji ,jj,jk-1) + umask(ji-1,jj,jk) & … … 158 158 ! 159 159 IF( ln_traldf_msc ) THEN ! stabilizing vertical diffusivity coefficient 160 DO_3D _00_00(2, jpkm1 )160 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 161 161 akz(ji,jj,jk) = 0.25_wp * ( & 162 162 & ( pahu(ji ,jj,jk) + pahu(ji ,jj,jk-1) ) / ( e1u(ji ,jj) * e1u(ji ,jj) ) & … … 167 167 ! 168 168 IF( ln_traldf_blp ) THEN ! bilaplacian operator 169 DO_3D _10_10(2, jpkm1 )169 DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 170 170 akz(ji,jj,jk) = 16._wp & 171 171 & * ah_wslp2 (ji,jj,jk) & … … 175 175 END_3D 176 176 ELSEIF( ln_traldf_lap ) THEN ! laplacian operator 177 DO_3D _10_10(2, jpkm1 )177 DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 178 178 ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 179 179 zcoef0 = rDt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) … … 200 200 201 201 ! Horizontal tracer gradient 202 DO_3D _10_10(1, jpkm1 )202 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 203 203 zdit(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 204 204 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 205 205 END_3D 206 206 IF( ln_zps ) THEN ! botton and surface ocean correction of the horizontal gradient 207 DO_2D _10_10207 DO_2D( 1, 0, 1, 0 ) ! bottom correction (partial bottom cell) 208 208 zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 209 209 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 210 210 END_2D 211 211 IF( ln_isfcav ) THEN ! first wet level beneath a cavity 212 DO_2D _10_10212 DO_2D( 1, 0, 1, 0 ) 213 213 IF( miku(ji,jj) > 1 ) zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn) 214 214 IF( mikv(ji,jj) > 1 ) zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn) … … 229 229 ELSE ; zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk) 230 230 ENDIF 231 DO_2D _10_10231 DO_2D( 1, 0, 1, 0 ) !== Horizontal fluxes 232 232 zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) 233 233 zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) … … 250 250 END_2D 251 251 ! 252 DO_2D _00_00252 DO_2D( 0, 0, 0, 0 ) !== horizontal divergence and add to pta 253 253 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) & 254 254 & + zsign * ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk) + zftv(ji,jj,jk) - zftv(ji,jj-1,jk) ) & … … 266 266 ztfw(:,:, 1 ) = 0._wp ; ztfw(:,:,jpk) = 0._wp 267 267 268 DO_3D _00_00( 2, jpkm1)268 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! interior (2=<jk=<jpk-1) 269 269 ! 270 270 zmsku = wmask(ji,jj,jk) / MAX( umask(ji ,jj,jk-1) + umask(ji-1,jj,jk) & … … 288 288 ! !== add the vertical 33 flux ==! 289 289 IF( ln_traldf_lap ) THEN ! laplacian case: eddy coef = ah_wslp2 - akz 290 DO_3D _00_00(2, jpkm1 )290 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 291 291 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk) & 292 292 & * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) ) & … … 297 297 SELECT CASE( kpass ) 298 298 CASE( 1 ) ! 1st pass : eddy coef = ah_wslp2 299 DO_3D _00_00(2, jpkm1 )299 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 300 300 ztfw(ji,jj,jk) = & 301 301 & ztfw(ji,jj,jk) + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj) & … … 303 303 END_3D 304 304 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt and pt2 gradients, resp. 305 DO_3D _00_00(2, jpkm1 )305 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 306 306 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk) & 307 307 & * ( ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) ) & … … 311 311 ENDIF 312 312 ! 313 DO_3D _00_00( 1, jpkm1 )313 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Divergence of vertical fluxes added to pta ==! 314 314 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) & 315 315 & / e3t(ji,jj,jk,Kmm) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_lap_blp.F90
r13237 r13899 99 99 ELSE ; zsign = -1._wp 100 100 ENDIF 101 DO_3D _10_10(1, jpkm1 )101 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 102 102 zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) !!gm * umask(ji,jj,jk) pah masked! 103 103 zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) !!gm * vmask(ji,jj,jk) … … 108 108 ! ! =========== ! 109 109 ! 110 DO_3D _10_10( 1, jpkm1 )110 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) !== First derivative (gradient) ==! 111 111 ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) 112 112 ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) 113 113 END_3D 114 IF( ln_zps ) THEN ! set gradient at bottom/top ocean level115 DO_2D _10_10114 IF( ln_zps ) THEN ! set gradient at bottom/top ocean level 115 DO_2D( 1, 0, 1, 0 ) ! bottom 116 116 ztu(ji,jj,mbku(ji,jj)) = zaheeu(ji,jj,mbku(ji,jj)) * pgu(ji,jj,jn) 117 117 ztv(ji,jj,mbkv(ji,jj)) = zaheev(ji,jj,mbkv(ji,jj)) * pgv(ji,jj,jn) 118 118 END_2D 119 IF( ln_isfcav ) THEN ! top in ocean cavities only120 DO_2D _10_10119 IF( ln_isfcav ) THEN ! top in ocean cavities only 120 DO_2D( 1, 0, 1, 0 ) 121 121 IF( miku(ji,jj) > 1 ) ztu(ji,jj,miku(ji,jj)) = zaheeu(ji,jj,miku(ji,jj)) * pgui(ji,jj,jn) 122 122 IF( mikv(ji,jj) > 1 ) ztv(ji,jj,mikv(ji,jj)) = zaheev(ji,jj,mikv(ji,jj)) * pgvi(ji,jj,jn) … … 125 125 ENDIF 126 126 ! 127 DO_3D _00_00( 1, jpkm1 )127 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Second derivative (divergence) added to the general tracer trends ==! 128 128 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) & 129 129 & + ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) & -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_triad.F90
r13237 r13899 137 137 DO ip = 0, 1 ! i-k triads 138 138 DO kp = 0, 1 139 DO_3D _10_10(1, jpkm1 )139 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 140 140 ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,Kmm) 141 141 zbu = e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) … … 157 157 DO jp = 0, 1 ! j-k triads 158 158 DO kp = 0, 1 159 DO_3D _10_10(1, jpkm1 )159 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 160 160 ze3wr = 1.0_wp / e3w(ji,jj+jp,jk+kp,Kmm) 161 161 zbv = e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) … … 179 179 ! 180 180 IF( ln_traldf_blp ) THEN ! bilaplacian operator 181 DO_3D _10_10(2, jpkm1 )181 DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 182 182 akz(ji,jj,jk) = 16._wp & 183 183 & * ah_wslp2 (ji,jj,jk) & … … 187 187 END_3D 188 188 ELSEIF( ln_traldf_lap ) THEN ! laplacian operator 189 DO_3D _10_10(2, jpkm1 )189 DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 190 190 ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 191 191 zcoef0 = rDt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) … … 211 211 zftv(:,:,:) = 0._wp 212 212 ! 213 DO_3D _10_10( 1, jpkm1 )213 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) !== before lateral T & S gradients at T-level jk ==! 214 214 zdit(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 215 215 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 216 216 END_3D 217 217 IF( ln_zps .AND. l_grad_zps ) THEN ! partial steps: correction at top/bottom ocean level 218 DO_2D _10_10218 DO_2D( 1, 0, 1, 0 ) ! bottom level 219 219 zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 220 220 zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 221 221 END_2D 222 222 IF( ln_isfcav ) THEN ! top level (ocean cavities only) 223 DO_2D _10_10223 DO_2D( 1, 0, 1, 0 ) 224 224 IF( miku(ji,jj) > 1 ) zdit(ji,jj,miku(ji,jj) ) = pgui(ji,jj,jn) 225 225 IF( mikv(ji,jj) > 1 ) zdjt(ji,jj,mikv(ji,jj) ) = pgvi(ji,jj,jn) … … 246 246 DO ip = 0, 1 !== Horizontal & vertical fluxes 247 247 DO kp = 0, 1 248 DO_2D _10_10248 DO_2D( 1, 0, 1, 0 ) 249 249 ze1ur = r1_e1u(ji,jj) 250 250 zdxt = zdit(ji,jj,jk) * ze1ur … … 267 267 DO jp = 0, 1 268 268 DO kp = 0, 1 269 DO_2D _10_10269 DO_2D( 1, 0, 1, 0 ) 270 270 ze2vr = r1_e2v(ji,jj) 271 271 zdyt = zdjt(ji,jj,jk) * ze2vr … … 289 289 DO ip = 0, 1 !== Horizontal & vertical fluxes 290 290 DO kp = 0, 1 291 DO_2D _10_10291 DO_2D( 1, 0, 1, 0 ) 292 292 ze1ur = r1_e1u(ji,jj) 293 293 zdxt = zdit(ji,jj,jk) * ze1ur … … 310 310 DO jp = 0, 1 311 311 DO kp = 0, 1 312 DO_2D _10_10312 DO_2D( 1, 0, 1, 0 ) 313 313 ze2vr = r1_e2v(ji,jj) 314 314 zdyt = zdjt(ji,jj,jk) * ze2vr … … 329 329 ENDIF 330 330 ! !== horizontal divergence and add to the general trend ==! 331 DO_2D _00_00331 DO_2D( 0, 0, 0, 0 ) 332 332 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) & 333 333 & + zsign * ( zftu(ji-1,jj ,jk) - zftu(ji,jj,jk) & … … 340 340 ! !== add the vertical 33 flux ==! 341 341 IF( ln_traldf_lap ) THEN ! laplacian case: eddy coef = ah_wslp2 - akz 342 DO_3D _10_00(2, jpkm1 )342 DO_3D( 1, 0, 0, 0, 2, jpkm1 ) 343 343 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) & 344 344 & * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) ) & … … 348 348 SELECT CASE( kpass ) 349 349 CASE( 1 ) ! 1st pass : eddy coef = ah_wslp2 350 DO_3D _10_00(2, jpkm1 )350 DO_3D( 1, 0, 0, 0, 2, jpkm1 ) 351 351 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) & 352 352 & * ah_wslp2(ji,jj,jk) * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) 353 353 END_3D 354 354 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt and pt2 gradients, resp. 355 DO_3D _10_00(2, jpkm1 )355 DO_3D( 1, 0, 0, 0, 2, jpkm1 ) 356 356 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) & 357 357 & * ( ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) ) & … … 361 361 ENDIF 362 362 ! 363 DO_3D _00_00( 1, jpkm1 )363 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Divergence of vertical fluxes added to pta ==! 364 364 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) & 365 365 & + zsign * ( ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk) ) & -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/tramle.F90
r13237 r13899 100 100 inml_mle(:,:) = mbkt(:,:) + 1 ! init. to number of ocean w-level (T-level + 1) 101 101 IF ( nla10 > 0 ) THEN ! avoid case where first level is thicker than 10m 102 DO_3DS _11_11( jpkm1, nlb10, -1)102 DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! from the bottom to nlb10 (10m) 103 103 IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle ) inml_mle(ji,jj) = jk ! Mixed layer 104 104 END_3D … … 110 110 zbm (:,:) = 0._wp 111 111 zn2 (:,:) = 0._wp 112 DO_3D _11_11( 1, ikmax )112 DO_3D( 1, 1, 1, 1, 1, ikmax ) ! MLD and mean buoyancy and N2 over the mixed layer 113 113 zc = e3t(ji,jj,jk,Kmm) * REAL( MIN( MAX( 0, inml_mle(ji,jj)-jk ) , 1 ) ) ! zc being 0 outside the ML t-points 114 114 zmld(ji,jj) = zmld(ji,jj) + zc … … 119 119 SELECT CASE( nn_mld_uv ) ! MLD at u- & v-pts 120 120 CASE ( 0 ) != min of the 2 neighbour MLDs 121 DO_2D _10_10121 DO_2D( 1, 0, 1, 0 ) 122 122 zhu(ji,jj) = MIN( zmld(ji+1,jj), zmld(ji,jj) ) 123 123 zhv(ji,jj) = MIN( zmld(ji,jj+1), zmld(ji,jj) ) 124 124 END_2D 125 125 CASE ( 1 ) != average of the 2 neighbour MLDs 126 DO_2D _10_10126 DO_2D( 1, 0, 1, 0 ) 127 127 zhu(ji,jj) = ( zmld(ji+1,jj) + zmld(ji,jj) ) * 0.5_wp 128 128 zhv(ji,jj) = ( zmld(ji,jj+1) + zmld(ji,jj) ) * 0.5_wp 129 129 END_2D 130 130 CASE ( 2 ) != max of the 2 neighbour MLDs 131 DO_2D _10_10131 DO_2D( 1, 0, 1, 0 ) 132 132 zhu(ji,jj) = MAX( zmld(ji+1,jj), zmld(ji,jj) ) 133 133 zhv(ji,jj) = MAX( zmld(ji,jj+1), zmld(ji,jj) ) … … 146 146 ! 147 147 IF( nn_mle == 0 ) THEN ! Fox-Kemper et al. 2010 formulation 148 DO_2D _10_10148 DO_2D( 1, 0, 1, 0 ) 149 149 zpsim_u(ji,jj) = rn_ce * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 150 150 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) & … … 157 157 ! 158 158 ELSEIF( nn_mle == 1 ) THEN ! New formulation (Lf = 5km fo/ff with fo=Coriolis parameter at latitude rn_lat) 159 DO_2D _10_10159 DO_2D( 1, 0, 1, 0 ) 160 160 zpsim_u(ji,jj) = rc_f * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 161 161 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) … … 167 167 ! 168 168 IF( nn_conv == 1 ) THEN ! No MLE in case of convection 169 DO_2D _10_10169 DO_2D( 1, 0, 1, 0 ) 170 170 IF( MIN( zn2(ji,jj) , zn2(ji+1,jj) ) < 0._wp ) zpsim_u(ji,jj) = 0._wp 171 171 IF( MIN( zn2(ji,jj) , zn2(ji,jj+1) ) < 0._wp ) zpsim_v(ji,jj) = 0._wp … … 174 174 ! 175 175 ! !== structure function value at uw- and vw-points ==! 176 DO_2D _10_10176 DO_2D( 1, 0, 1, 0 ) 177 177 zhu(ji,jj) = 1._wp / zhu(ji,jj) ! hu --> 1/hu 178 178 zhv(ji,jj) = 1._wp / zhv(ji,jj) … … 182 182 zpsi_vw(:,:,:) = 0._wp 183 183 ! 184 DO_3D _10_10( 2, ikmax )184 DO_3D( 1, 0, 1, 0, 2, ikmax ) ! start from 2 : surface value = 0 185 185 zcuw = 1._wp - ( gdepw(ji+1,jj,jk,Kmm) + gdepw(ji,jj,jk,Kmm) ) * zhu(ji,jj) 186 186 zcvw = 1._wp - ( gdepw(ji,jj+1,jk,Kmm) + gdepw(ji,jj,jk,Kmm) ) * zhv(ji,jj) … … 196 196 ! !== transport increased by the MLE induced transport ==! 197 197 DO jk = 1, ikmax 198 DO_2D _10_10198 DO_2D( 1, 0, 1, 0 ) ! CAUTION pu,pv must be defined at row/column i=1 / j=1 199 199 pu(ji,jj,jk) = pu(ji,jj,jk) + ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji,jj,jk+1) ) 200 200 pv(ji,jj,jk) = pv(ji,jj,jk) + ( zpsi_vw(ji,jj,jk) - zpsi_vw(ji,jj,jk+1) ) 201 201 END_2D 202 DO_2D _00_00202 DO_2D( 0, 0, 0, 0 ) 203 203 pw(ji,jj,jk) = pw(ji,jj,jk) - ( zpsi_uw(ji,jj,jk) - zpsi_uw(ji-1,jj,jk) & 204 204 & + zpsi_vw(ji,jj,jk) - zpsi_vw(ji,jj-1,jk) ) … … 283 283 IF( ierr /= 0 ) CALL ctl_stop( 'tra_adv_mle_init: failed to allocate arrays' ) 284 284 z1_t2 = 1._wp / ( rn_time * rn_time ) 285 DO_2D _01_01285 DO_2D( 0, 1, 0, 1 ) ! "coriolis+ time^-1" at u- & v-points 286 286 zfu = ( ff_f(ji,jj) + ff_f(ji,jj-1) ) * 0.5_wp 287 287 zfv = ( ff_f(ji,jj) + ff_f(ji-1,jj) ) * 0.5_wp -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/tranpc.F90
r13237 r13899 103 103 inpcc = 0 104 104 ! 105 DO_2D _00_00105 DO_2D( 0, 0, 0, 0 ) ! interior column only 106 106 ! 107 107 IF( tmask(ji,jj,2) == 1 ) THEN ! At least 2 ocean points -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traqsr.F90
r13237 r13899 63 63 REAL(wp) :: xsi1r ! inverse of rn_si1 64 64 ! 65 REAL(wp) , DIMENSION(3,61):: rkrgb ! tabulated attenuation coefficients for RGB absorption65 REAL(wp) , PUBLIC, DIMENSION(3,61) :: rkrgb ! tabulated attenuation coefficients for RGB absorption 66 66 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_chl ! structure of input Chl (file informations, fields read) 67 67 … … 138 138 IF(lwp) WRITE(numout,*) ' nit000-1 qsr tracer content forcing field read in the restart file' 139 139 z1_2 = 0.5_wp 140 CALL iom_get( numror, jpdom_auto glo, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios ) ! before heat content trend due to Qsr flux140 CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b, ldxios = lrxios ) ! before heat content trend due to Qsr flux 141 141 ELSE ! No restart or restart not found: Euler forward time stepping 142 142 z1_2 = 1._wp … … 172 172 ! most expensive calculations) 173 173 ! 174 DO_2D _00_00174 DO_2D( 0, 0, 0, 0 ) 175 175 ! zlogc = log(zchl) 176 176 zlogc = LOG ( MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) ) … … 191 191 192 192 ! 193 DO_3D _00_00 (1, nksr + 1 )193 DO_3D( 0, 0, 0, 0, 1, nksr + 1 ) 194 194 ! zchl = ALOG( ze0(ji,jj) ) 195 195 zlogc = ze0(ji,jj) … … 221 221 ! 222 222 zcoef = ( 1. - rn_abs ) / 3._wp !* surface equi-partition in R-G-B 223 DO_2D _00_00223 DO_2D( 0, 0, 0, 0 ) 224 224 ze0(ji,jj) = rn_abs * qsr(ji,jj) 225 225 ze1(ji,jj) = zcoef * qsr(ji,jj) … … 231 231 END_2D 232 232 ! 233 ! * interior equi-partition in R-G-B depending on vertical profile of Chl234 DO_3D _00_00 (2, nksr + 1 )233 ! !* interior equi-partition in R-G-B depending on vertical profile of Chl 234 DO_3D( 0, 0, 0, 0, 2, nksr + 1 ) 235 235 ze3t = e3t(ji,jj,jk-1,Kmm) 236 236 irgb = NINT( ztmp3d(ji,jj,jk) ) … … 246 246 END_3D 247 247 ! 248 DO_3D _00_00( 1, nksr )248 DO_3D( 0, 0, 0, 0, 1, nksr ) !* now qsr induced heat content 249 249 qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( ztmp3d(ji,jj,jk) - ztmp3d(ji,jj,jk+1) ) 250 250 END_3D … … 256 256 zz0 = rn_abs * r1_rho0_rcp ! surface equi-partition in 2-bands 257 257 zz1 = ( 1. - rn_abs ) * r1_rho0_rcp 258 DO_3D _00_00( 1, nksr )258 DO_3D( 0, 0, 0, 0, 1, nksr ) ! solar heat absorbed at T-point in the top 400m 259 259 zc0 = zz0 * EXP( -gdepw(ji,jj,jk ,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk ,Kmm)*xsi1r ) 260 260 zc1 = zz0 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi1r ) … … 265 265 ! 266 266 ! !-----------------------------! 267 DO_3D_00_00( 1, nksr ) 267 ! ! update to the temp. trend ! 268 ! !-----------------------------! 269 DO_3D( 0, 0, 0, 0, 1, nksr ) 268 270 pts(ji,jj,jk,jp_tem,Krhs) = pts(ji,jj,jk,jp_tem,Krhs) & 269 271 & + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) & … … 272 274 ! 273 275 ! sea-ice: store the 1st ocean level attenuation coefficient 274 DO_2D _00_00276 DO_2D( 0, 0, 0, 0 ) 275 277 IF( qsr(ji,jj) /= 0._wp ) THEN ; fraqsr_1lev(ji,jj) = qsr_hc(ji,jj,1) / ( r1_rho0_rcp * qsr(ji,jj) ) 276 278 ELSE ; fraqsr_1lev(ji,jj) = 1._wp … … 417 419 IF( .NOT.lk_top ) CALL ctl_stop( 'No bio model : ln_qsr_bio = true impossible ' ) 418 420 ! 421 CALL trc_oce_rgb( rkrgb ) ! tabulated attenuation coef. 422 ! 423 nksr = trc_oce_ext_lev( r_si2, 33._wp ) ! level of light extinction 424 ! 425 IF(lwp) WRITE(numout,*) ' level of light extinction = ', nksr, ' ref depth = ', gdepw_1d(nksr+1), ' m' 426 ! 419 427 END SELECT 420 428 ! … … 423 431 ! 1st ocean level attenuation coefficient (used in sbcssm) 424 432 IF( iom_varid( numror, 'fraqsr_1lev', ldstop = .FALSE. ) > 0 ) THEN 425 CALL iom_get( numror, jpdom_auto glo, 'fraqsr_1lev' , fraqsr_1lev, ldxios = lrxios )433 CALL iom_get( numror, jpdom_auto, 'fraqsr_1lev' , fraqsr_1lev, ldxios = lrxios ) 426 434 ELSE 427 435 fraqsr_1lev(:,:) = 1._wp ! default : no penetration -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/trasbc.F90
r13237 r13899 112 112 zfact = 0.5_wp 113 113 sbc_tsc(:,:,:) = 0._wp 114 CALL iom_get( numror, jpdom_auto glo, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem), ldxios = lrxios ) ! before heat content sbc trend115 CALL iom_get( numror, jpdom_auto glo, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal), ldxios = lrxios ) ! before salt content sbc trend114 CALL iom_get( numror, jpdom_auto, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem), ldxios = lrxios ) ! before heat content sbc trend 115 CALL iom_get( numror, jpdom_auto, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal), ldxios = lrxios ) ! before salt content sbc trend 116 116 ELSE ! No restart or restart not found: Euler forward time stepping 117 117 zfact = 1._wp … … 124 124 ENDIF 125 125 ! !== Now sbc tracer content fields ==! 126 DO_2D _01_00126 DO_2D( 0, 1, 0, 0 ) 127 127 sbc_tsc(ji,jj,jp_tem) = r1_rho0_rcp * qns(ji,jj) ! non solar heat flux 128 128 sbc_tsc(ji,jj,jp_sal) = r1_rho0 * sfx(ji,jj) ! salt flux due to freezing/melting 129 129 END_2D 130 130 IF( ln_linssh ) THEN !* linear free surface 131 DO_2D _01_00131 DO_2D( 0, 1, 0, 0 ) !==>> add concentration/dilution effect due to constant volume cell 132 132 sbc_tsc(ji,jj,jp_tem) = sbc_tsc(ji,jj,jp_tem) + r1_rho0 * emp(ji,jj) * pts(ji,jj,1,jp_tem,Kmm) 133 133 sbc_tsc(ji,jj,jp_sal) = sbc_tsc(ji,jj,jp_sal) + r1_rho0 * emp(ji,jj) * pts(ji,jj,1,jp_sal,Kmm) 134 END_2D 134 END_2D !==>> output c./d. term 135 135 IF( iom_use('emp_x_sst') ) CALL iom_put( "emp_x_sst", emp (:,:) * pts(:,:,1,jp_tem,Kmm) ) 136 136 IF( iom_use('emp_x_sss') ) CALL iom_put( "emp_x_sss", emp (:,:) * pts(:,:,1,jp_sal,Kmm) ) … … 138 138 ! 139 139 DO jn = 1, jpts !== update tracer trend ==! 140 DO_2D _01_00140 DO_2D( 0, 1, 0, 0 ) 141 141 pts(ji,jj,1,jn,Krhs) = pts(ji,jj,1,jn,Krhs) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) & 142 142 & / e3t(ji,jj,1,Kmm) … … 157 157 IF( ln_rnf ) THEN ! input of heat and salt due to river runoff 158 158 zfact = 0.5_wp 159 DO_2D _01_00159 DO_2D( 0, 1, 0, 0 ) 160 160 IF( rnf(ji,jj) /= 0._wp ) THEN 161 161 zdep = zfact / h_rnf(ji,jj) … … 182 182 ! 183 183 IF( ln_linssh ) THEN 184 DO_2D _01_00184 DO_2D( 0, 1, 0, 0 ) 185 185 ztim = ssh_iau(ji,jj) / e3t(ji,jj,1,Kmm) 186 186 pts(ji,jj,1,jp_tem,Krhs) = pts(ji,jj,1,jp_tem,Krhs) + pts(ji,jj,1,jp_tem,Kmm) * ztim … … 188 188 END_2D 189 189 ELSE 190 DO_2D _01_00190 DO_2D( 0, 1, 0, 0 ) 191 191 ztim = ssh_iau(ji,jj) / ( ht(ji,jj) + 1. - ssmask(ji, jj) ) 192 192 pts(ji,jj,:,jp_tem,Krhs) = pts(ji,jj,:,jp_tem,Krhs) + pts(ji,jj,:,jp_tem,Kmm) * ztim -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/trazdf.F90
r13237 r13899 161 161 IF( l_ldfslp ) THEN ! isoneutral diffusion: add the contribution 162 162 IF( ln_traldf_msc ) THEN ! MSC iso-neutral operator 163 DO_3D _00_00(2, jpkm1 )163 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 164 164 zwt(ji,jj,jk) = zwt(ji,jj,jk) + akz(ji,jj,jk) 165 165 END_3D 166 166 ELSE ! standard or triad iso-neutral operator 167 DO_3D _00_00(2, jpkm1 )167 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 168 168 zwt(ji,jj,jk) = zwt(ji,jj,jk) + ah_wslp2(ji,jj,jk) 169 169 END_3D … … 173 173 ! Diagonal, lower (i), upper (s) (including the bottom boundary condition since avt is masked) 174 174 IF( ln_zad_Aimp ) THEN ! Adaptive implicit vertical advection 175 DO_3D _00_00(1, jpkm1 )175 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 176 176 zzwi = - p2dt * zwt(ji,jj,jk ) / e3w(ji,jj,jk ,Kmm) 177 177 zzws = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) … … 182 182 END_3D 183 183 ELSE 184 DO_3D _00_00(1, jpkm1 )184 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 185 185 zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk ) / e3w(ji,jj,jk,Kmm) 186 186 zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) … … 208 208 ! used as a work space array: its value is modified. 209 209 ! 210 DO_2D _00_00210 DO_2D( 0, 0, 0, 0 ) !* 1st recurrence: Tk = Dk - Ik Sk-1 / Tk-1 (increasing k) ! done one for all passive tracers (so included in the IF instruction) 211 211 zwt(ji,jj,1) = zwd(ji,jj,1) 212 212 END_2D 213 DO_3D _00_00(2, jpkm1 )213 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 214 214 zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwt(ji,jj,jk-1) 215 215 END_3D … … 217 217 ENDIF 218 218 ! 219 DO_2D _00_00219 DO_2D( 0, 0, 0, 0 ) !* 2nd recurrence: Zk = Yk - Ik / Tk-1 Zk-1 220 220 pt(ji,jj,1,jn,Kaa) = e3t(ji,jj,1,Kbb) * pt(ji,jj,1,jn,Kbb) & 221 221 & + p2dt * e3t(ji,jj,1,Kmm) * pt(ji,jj,1,jn,Krhs) 222 222 END_2D 223 DO_3D _00_00(2, jpkm1 )223 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 224 224 zrhs = e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) & 225 225 & + p2dt * e3t(ji,jj,jk,Kmm) * pt(ji,jj,jk,jn,Krhs) ! zrhs=right hand side … … 227 227 END_3D 228 228 ! 229 DO_2D _00_00229 DO_2D( 0, 0, 0, 0 ) !* 3d recurrence: Xk = (Zk - Sk Xk+1 ) / Tk (result is the after tracer) 230 230 pt(ji,jj,jpkm1,jn,Kaa) = pt(ji,jj,jpkm1,jn,Kaa) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1) 231 231 END_2D 232 DO_3DS _00_00(jpk-2, 1, -1 )232 DO_3DS( 0, 0, 0, 0, jpk-2, 1, -1 ) 233 233 pt(ji,jj,jk,jn,Kaa) = ( pt(ji,jj,jk,jn,Kaa) - zws(ji,jj,jk) * pt(ji,jj,jk+1,jn,Kaa) ) & 234 234 & / zwt(ji,jj,jk) * tmask(ji,jj,jk) -
NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/zpshde.F90
r13237 r13899 107 107 DO jn = 1, kjpt !== Interpolation of tracers at the last ocean level ==! 108 108 ! 109 DO_2D _10_10109 DO_2D( 1, 0, 1, 0 ) 110 110 iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points 111 111 ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 ) ! if level first is a p-step, ik.m1=1 … … 151 151 pgru(:,:) = 0._wp 152 152 pgrv(:,:) = 0._wp ! depth of the partial step level 153 DO_2D _10_10153 DO_2D( 1, 0, 1, 0 ) 154 154 iku = mbku(ji,jj) 155 155 ikv = mbkv(ji,jj) … … 167 167 CALL eos( ztj, zhj, zrj ) ! at the partial step depth output in zri, zrj 168 168 ! 169 DO_2D _10_10169 DO_2D( 1, 0, 1, 0 ) ! Gradient of density at the last level 170 170 iku = mbku(ji,jj) 171 171 ikv = mbkv(ji,jj) … … 262 262 DO jn = 1, kjpt !== Interpolation of tracers at the last ocean level ==! 263 263 ! 264 DO_2D _10_10264 DO_2D( 1, 0, 1, 0 ) 265 265 266 266 iku = mbku(ji,jj); ikum1 = MAX( iku - 1 , 1 ) ! last and before last ocean level at u- & v-points … … 308 308 pgru(:,:)=0.0_wp ; pgrv(:,:)=0.0_wp ; 309 309 ! 310 DO_2D _10_10310 DO_2D( 1, 0, 1, 0 ) 311 311 312 312 iku = mbku(ji,jj) … … 329 329 CALL eos( ztj, zhj, zrj ) 330 330 331 DO_2D _10_10331 DO_2D( 1, 0, 1, 0 ) ! Gradient of density at the last level 332 332 iku = mbku(ji,jj) 333 333 ikv = mbkv(ji,jj) … … 351 351 ! 352 352 DO jn = 1, kjpt !== Interpolation of tracers at the last ocean level ==! ! 353 DO_2D _10_10353 DO_2D( 1, 0, 1, 0 ) 354 354 iku = miku(ji,jj); ikup1 = miku(ji,jj) + 1 355 355 ikv = mikv(ji,jj); ikvp1 = mikv(ji,jj) + 1 … … 400 400 ! 401 401 pgrui(:,:) =0.0_wp; pgrvi(:,:) =0.0_wp; 402 DO_2D _10_10402 DO_2D( 1, 0, 1, 0 ) 403 403 404 404 iku = miku(ji,jj) … … 420 420 CALL eos( ztj, zhj, zrj ) ! at the partial step depth output in zri, zrj 421 421 ! 422 DO_2D _10_10422 DO_2D( 1, 0, 1, 0 ) ! Gradient of density at the last level 423 423 iku = miku(ji,jj) 424 424 ikv = mikv(ji,jj)
Note: See TracChangeset
for help on using the changeset viewer.
