Changeset 14045 for NEMO/trunk/src/OCE/TRA
- Timestamp:
- 2020-12-03T13:04:25+01:00 (3 years ago)
- Location:
- NEMO/trunk/src/OCE/TRA
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/TRA/traatf.F90
r13982 r14045 117 117 IF( l_trdtra ) THEN 118 118 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 119 ztrdt(:,:, jpk) = 0._wp120 ztrds(:,:, jpk) = 0._wp119 ztrdt(:,:,:) = 0._wp 120 ztrds(:,:,:) = 0._wp 121 121 IF( ln_traldf_iso ) THEN ! diagnose the "pure" Kz diffusive trend 122 122 CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_zdfp, ztrdt ) -
NEMO/trunk/src/OCE/TRA/tramle.F90
r13982 r14045 20 20 USE lib_mpp ! MPP library 21 21 USE lbclnk ! lateral boundary condition / mpp link 22 23 ! where OSMOSIS_OBL is used with integrated FK 24 USE zdf_oce, ONLY : ln_zdfosm 25 USE zdfosm, ONLY : ln_osm_mle, hmle, dbdx_mle, dbdy_mle, mld_prof 22 26 23 27 IMPLICIT NONE … … 99 103 !!---------------------------------------------------------------------- 100 104 ! 101 ! !== MLD used for MLE ==! 102 ! ! compute from the 10m density to deal with the diurnal cycle 103 DO_2D( 1, 1, 1, 1 ) 104 inml_mle(ji,jj) = mbkt(ji,jj) + 1 ! init. to number of ocean w-level (T-level + 1) 105 END_2D 106 IF ( nla10 > 0 ) THEN ! avoid case where first level is thicker than 10m 107 DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! from the bottom to nlb10 (10m) 108 IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle ) inml_mle(ji,jj) = jk ! Mixed layer 105 ! 106 IF(ln_osm_mle.and.ln_zdfosm) THEN 107 ikmax = MIN( MAXVAL( mld_prof(:,:) ), jpkm1 ) ! max level of the computation 108 ! 109 ! 110 SELECT CASE( nn_mld_uv ) ! MLD at u- & v-pts 111 CASE ( 0 ) != min of the 2 neighbour MLDs 112 DO_2D( 1, 0, 1, 0 ) 113 zhu(ji,jj) = MIN( hmle(ji+1,jj), hmle(ji,jj) ) 114 zhv(ji,jj) = MIN( hmle(ji,jj+1), hmle(ji,jj) ) 115 END_2D 116 CASE ( 1 ) != average of the 2 neighbour MLDs 117 DO_2D( 1, 0, 1, 0 ) 118 zhu(ji,jj) = MAX( hmle(ji+1,jj), hmle(ji,jj) ) 119 zhv(ji,jj) = MAX( hmle(ji,jj+1), hmle(ji,jj) ) 120 END_2D 121 CASE ( 2 ) != max of the 2 neighbour MLDs 122 DO_2D( 1, 0, 1, 0 ) 123 zhu(ji,jj) = MAX( hmle(ji+1,jj), hmle(ji,jj) ) 124 zhv(ji,jj) = MAX( hmle(ji,jj+1), hmle(ji,jj) ) 125 END_2D 126 END SELECT 127 IF( nn_mle == 0 ) THEN ! Fox-Kemper et al. 2010 formulation 128 DO_2D( 1, 0, 1, 0 ) 129 zpsim_u(ji,jj) = rn_ce * zhu(ji,jj) * zhu(ji,jj) * e2u(ji,jj) & 130 & * dbdx_mle(ji,jj) * MIN( 111.e3_wp , e1u(ji,jj) ) & 131 & / ( MAX( rn_lf * rfu(ji,jj) , SQRT( rb_c * zhu(ji,jj) ) ) ) 132 ! 133 zpsim_v(ji,jj) = rn_ce * zhv(ji,jj) * zhv(ji,jj) * e1v(ji,jj) & 134 & * dbdy_mle(ji,jj) * MIN( 111.e3_wp , e2v(ji,jj) ) & 135 & / ( MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) ) ) 136 END_2D 137 ! 138 ELSEIF( nn_mle == 1 ) THEN ! New formulation (Lf = 5km fo/ff with fo=Coriolis parameter at latitude rn_lat) 139 DO_2D( 1, 0, 1, 0 ) 140 zpsim_u(ji,jj) = rc_f * zhu(ji,jj) * zhu(ji,jj) * e2u(ji,jj) & 141 & * dbdx_mle(ji,jj) * MIN( 111.e3_wp , e1u(ji,jj) ) 142 ! 143 zpsim_v(ji,jj) = rc_f * zhv(ji,jj) * zhv(ji,jj) * e1v(ji,jj) & 144 & * dbdy_mle(ji,jj) * MIN( 111.e3_wp , e2v(ji,jj) ) 145 END_2D 146 ENDIF 147 148 ELSE !do not use osn_mle 149 ! !== MLD used for MLE ==! 150 ! ! compute from the 10m density to deal with the diurnal cycle 151 DO_2D( 1, 1, 1, 1 ) 152 inml_mle(ji,jj) = mbkt(ji,jj) + 1 ! init. to number of ocean w-level (T-level + 1) 153 END_2D 154 IF ( nla10 > 0 ) THEN ! avoid case where first level is thicker than 10m 155 DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! from the bottom to nlb10 (10m) 156 IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle ) inml_mle(ji,jj) = jk ! Mixed layer 157 END_3D 158 ENDIF 159 ikmax = MIN( MAXVAL( inml_mle(:,:) ), jpkm1 ) ! max level of the computation 160 ! 161 ! 162 zmld(:,:) = 0._wp !== Horizontal shape of the MLE ==! 163 zbm (:,:) = 0._wp 164 zn2 (:,:) = 0._wp 165 DO_3D( 1, 1, 1, 1, 1, ikmax ) ! MLD and mean buoyancy and N2 over the mixed layer 166 zc = e3t(ji,jj,jk,Kmm) * REAL( MIN( MAX( 0, inml_mle(ji,jj)-jk ) , 1 ) ) ! zc being 0 outside the ML t-points 167 zmld(ji,jj) = zmld(ji,jj) + zc 168 zbm (ji,jj) = zbm (ji,jj) + zc * (rho0 - rhop(ji,jj,jk) ) * r1_rho0 169 zn2 (ji,jj) = zn2 (ji,jj) + zc * (rn2(ji,jj,jk)+rn2(ji,jj,jk+1))*0.5_wp 109 170 END_3D 110 ENDIF 111 ikmax = MIN( MAXVAL( inml_mle(:,:) ), jpkm1 ) ! max level of the computation 112 ! 113 ! 114 zmld(:,:) = 0._wp !== Horizontal shape of the MLE ==! 115 zbm (:,:) = 0._wp 116 zn2 (:,:) = 0._wp 117 DO_3D( 1, 1, 1, 1, 1, ikmax ) ! MLD and mean buoyancy and N2 over the mixed layer 118 zc = e3t(ji,jj,jk,Kmm) * REAL( MIN( MAX( 0, inml_mle(ji,jj)-jk ) , 1 ) ) ! zc being 0 outside the ML t-points 119 zmld(ji,jj) = zmld(ji,jj) + zc 120 zbm (ji,jj) = zbm (ji,jj) + zc * (rho0 - rhop(ji,jj,jk) ) * r1_rho0 121 zn2 (ji,jj) = zn2 (ji,jj) + zc * (rn2(ji,jj,jk)+rn2(ji,jj,jk+1))*0.5_wp 122 END_3D 123 124 SELECT CASE( nn_mld_uv ) ! MLD at u- & v-pts 125 CASE ( 0 ) != min of the 2 neighbour MLDs 126 DO_2D( 1, 0, 1, 0 ) 127 zhu(ji,jj) = MIN( zmld(ji+1,jj), zmld(ji,jj) ) 128 zhv(ji,jj) = MIN( zmld(ji,jj+1), zmld(ji,jj) ) 171 172 SELECT CASE( nn_mld_uv ) ! MLD at u- & v-pts 173 CASE ( 0 ) != min of the 2 neighbour MLDs 174 DO_2D( 1, 0, 1, 0 ) 175 zhu(ji,jj) = MIN( zmld(ji+1,jj), zmld(ji,jj) ) 176 zhv(ji,jj) = MIN( zmld(ji,jj+1), zmld(ji,jj) ) 177 END_2D 178 CASE ( 1 ) != average of the 2 neighbour MLDs 179 DO_2D( 1, 0, 1, 0 ) 180 zhu(ji,jj) = ( zmld(ji+1,jj) + zmld(ji,jj) ) * 0.5_wp 181 zhv(ji,jj) = ( zmld(ji,jj+1) + zmld(ji,jj) ) * 0.5_wp 182 END_2D 183 CASE ( 2 ) != max of the 2 neighbour MLDs 184 DO_2D( 1, 0, 1, 0 ) 185 zhu(ji,jj) = MAX( zmld(ji+1,jj), zmld(ji,jj) ) 186 zhv(ji,jj) = MAX( zmld(ji,jj+1), zmld(ji,jj) ) 187 END_2D 188 END SELECT 189 ! ! convert density into buoyancy 190 DO_2D( 1, 1, 1, 1 ) 191 zbm(ji,jj) = + grav * zbm(ji,jj) / MAX( e3t(ji,jj,1,Kmm), zmld(ji,jj) ) 129 192 END_2D 130 CASE ( 1 ) != average of the 2 neighbour MLDs 131 DO_2D( 1, 0, 1, 0 ) 132 zhu(ji,jj) = ( zmld(ji+1,jj) + zmld(ji,jj) ) * 0.5_wp 133 zhv(ji,jj) = ( zmld(ji,jj+1) + zmld(ji,jj) ) * 0.5_wp 134 END_2D 135 CASE ( 2 ) != max of the 2 neighbour MLDs 136 DO_2D( 1, 0, 1, 0 ) 137 zhu(ji,jj) = MAX( zmld(ji+1,jj), zmld(ji,jj) ) 138 zhv(ji,jj) = MAX( zmld(ji,jj+1), zmld(ji,jj) ) 139 END_2D 140 END SELECT 141 ! ! convert density into buoyancy 142 DO_2D( 1, 1, 1, 1 ) 143 zbm(ji,jj) = + grav * zbm(ji,jj) / MAX( e3t(ji,jj,1,Kmm), zmld(ji,jj) ) 144 END_2D 145 ! 146 ! 147 ! !== Magnitude of the MLE stream function ==! 148 ! 149 ! di[bm] Ds 150 ! Psi = Ce H^2 ---------------- e2u mu(z) where fu Lf = MAX( fu*rn_fl , (Db H)^1/2 ) 151 ! e1u Lf fu and the e2u for the "transport" 152 ! (not *e3u as divided by e3u at the end) 153 ! 154 IF( nn_mle == 0 ) THEN ! Fox-Kemper et al. 2010 formulation 155 DO_2D( 1, 0, 1, 0 ) 156 zpsim_u(ji,jj) = rn_ce * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 157 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) & 158 & / ( MAX( rn_lf * rfu(ji,jj) , SQRT( rb_c * zhu(ji,jj) ) ) ) 193 ! 194 ! 195 ! !== Magnitude of the MLE stream function ==! 196 ! 197 ! di[bm] Ds 198 ! Psi = Ce H^2 ---------------- e2u mu(z) where fu Lf = MAX( fu*rn_fl , (Db H)^1/2 ) 199 ! e1u Lf fu and the e2u for the "transport" 200 ! (not *e3u as divided by e3u at the end) 201 ! 202 IF( nn_mle == 0 ) THEN ! Fox-Kemper et al. 2010 formulation 203 DO_2D( 1, 0, 1, 0 ) 204 zpsim_u(ji,jj) = rn_ce * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 205 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) & 206 & / ( MAX( rn_lf * rfu(ji,jj) , SQRT( rb_c * zhu(ji,jj) ) ) ) 159 207 ! 160 zpsim_v(ji,jj) = rn_ce * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) &161 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) &162 & / ( MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) ) )163 END_2D164 !165 ELSEIF( nn_mle == 1 ) THEN ! New formulation (Lf = 5km fo/ff with fo=Coriolis parameter at latitude rn_lat)166 DO_2D( 1, 0, 1, 0 )167 zpsim_u(ji,jj) = rc_f * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) &168 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) )208 zpsim_v(ji,jj) = rn_ce * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) & 209 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) & 210 & / ( MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) ) ) 211 END_2D 212 ! 213 ELSEIF( nn_mle == 1 ) THEN ! New formulation (Lf = 5km fo/ff with fo=Coriolis parameter at latitude rn_lat) 214 DO_2D( 1, 0, 1, 0 ) 215 zpsim_u(ji,jj) = rc_f * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 216 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) 169 217 ! 170 zpsim_v(ji,jj) = rc_f * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) & 171 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) 172 END_2D 173 ENDIF 174 ! 175 IF( nn_conv == 1 ) THEN ! No MLE in case of convection 176 DO_2D( 1, 0, 1, 0 ) 177 IF( MIN( zn2(ji,jj) , zn2(ji+1,jj) ) < 0._wp ) zpsim_u(ji,jj) = 0._wp 178 IF( MIN( zn2(ji,jj) , zn2(ji,jj+1) ) < 0._wp ) zpsim_v(ji,jj) = 0._wp 179 END_2D 180 ENDIF 181 ! 182 ! !== structure function value at uw- and vw-points ==! 183 DO_2D( 1, 0, 1, 0 ) 184 zhu(ji,jj) = 1._wp / zhu(ji,jj) ! hu --> 1/hu 185 zhv(ji,jj) = 1._wp / zhv(ji,jj) 186 END_2D 187 ! 188 zpsi_uw(:,:,:) = 0._wp 189 zpsi_vw(:,:,:) = 0._wp 190 ! 218 zpsim_v(ji,jj) = rc_f * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) & 219 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) 220 END_2D 221 ENDIF 222 ! 223 IF( nn_conv == 1 ) THEN ! No MLE in case of convection 224 DO_2D( 1, 0, 1, 0 ) 225 IF( MIN( zn2(ji,jj) , zn2(ji+1,jj) ) < 0._wp ) zpsim_u(ji,jj) = 0._wp 226 IF( MIN( zn2(ji,jj) , zn2(ji,jj+1) ) < 0._wp ) zpsim_v(ji,jj) = 0._wp 227 END_2D 228 ENDIF 229 ! 230 ENDIF ! end of ln_osm_mle conditional 231 ! !== structure function value at uw- and vw-points ==! 232 DO_2D( 1, 0, 1, 0 ) 233 zhu(ji,jj) = 1._wp / MAX(zhu(ji,jj), rsmall) ! hu --> 1/hu 234 zhv(ji,jj) = 1._wp / MAX(zhv(ji,jj), rsmall) 235 END_2D 236 ! 237 zpsi_uw(:,:,:) = 0._wp 238 zpsi_vw(:,:,:) = 0._wp 239 ! 191 240 DO_3D( 1, 0, 1, 0, 2, ikmax ) ! start from 2 : surface value = 0 192 241 zcuw = 1._wp - ( gdepw(ji+1,jj,jk,Kmm) + gdepw(ji,jj,jk,Kmm) ) * zhu(ji,jj) … … 220 269 ENDIF 221 270 ! 222 DO_2D( 0, 0, 0, 0 ) 223 zLf_NH(ji,jj) = SQRT( rb_c * zmld(ji,jj) ) * r1_ft(ji,jj) ! Lf = N H / f 224 END_2D 271 IF (ln_osm_mle.and.ln_zdfosm) THEN 272 DO_2D( 0, 0, 0, 0 ) 273 zLf_NH(ji,jj) = SQRT( rb_c * hmle(ji,jj) ) * r1_ft(ji,jj) ! Lf = N H / f 274 END_2D 275 ELSE 276 DO_2D( 0, 0, 0, 0 ) 277 zLf_NH(ji,jj) = SQRT( rb_c * zmld(ji,jj) ) * r1_ft(ji,jj) ! Lf = N H / f 278 END_2D 279 ENDIF 225 280 ! 226 281 ! divide by cross distance to give streamfunction with dimensions m^2/s … … 239 294 ! 240 295 END SUBROUTINE tra_mle_trp 241 242 296 243 297 SUBROUTINE tra_mle_init … … 301 355 IF( ierr /= 0 ) CALL ctl_stop( 'tra_adv_mle_init: failed to allocate arrays' ) 302 356 z1_t2 = 1._wp / ( rn_time * rn_time ) 303 DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls) ! "coriolis+ time^-1" at u- & v-points357 DO_2D( 0, 1, 0, 1 ) ! "coriolis+ time^-1" at u- & v-points 304 358 zfu = ( ff_f(ji,jj) + ff_f(ji,jj-1) ) * 0.5_wp 305 359 zfv = ( ff_f(ji,jj) + ff_f(ji-1,jj) ) * 0.5_wp … … 307 361 rfv(ji,jj) = SQRT( zfv * zfv + z1_t2 ) 308 362 END_2D 309 IF (nn_hls.EQ.1)CALL lbc_lnk_multi( 'tramle', rfu, 'U', 1.0_wp , rfv, 'V', 1.0_wp )363 CALL lbc_lnk_multi( 'tramle', rfu, 'U', 1.0_wp , rfv, 'V', 1.0_wp ) 310 364 ! 311 365 ELSEIF( nn_mle == 1 ) THEN ! MLE array allocation & initialisation
Note: See TracChangeset
for help on using the changeset viewer.