- Timestamp:
- 2020-12-03T20:14:08+01:00 (4 years ago)
- Location:
- NEMO/branches/2020/dev_r13787_doc_latex_recovery
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/TRA/tramle.F90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13787_doc_latex_recovery
- Property svn:externals
-
old new 8 8 9 9 # SETTE 10 ^/utils/CI/sette @13559sette10 ^/utils/CI/sette_wave@13990 sette
-
- Property svn:externals
-
NEMO/branches/2020/dev_r13787_doc_latex_recovery/src/OCE/TRA/tramle.F90
r13497 r14066 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 … … 79 83 !! Fox-Kemper and Ferrari, JPO, 38, 1166-1179, 2008 80 84 !!---------------------------------------------------------------------- 81 INTEGER 82 INTEGER 83 INTEGER 84 CHARACTER(len=3) 85 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pu ! in : 3 ocean transport components86 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pv ! out: same 3 transport components87 REAL(wp), DIMENSION( jpi,jpj,jpk), INTENT(inout) :: pw ! increased by the MLE induced transport85 INTEGER , INTENT(in ) :: kt ! ocean time-step index 86 INTEGER , INTENT(in ) :: kit000 ! first time step index 87 INTEGER , INTENT(in ) :: Kmm ! ocean time level index 88 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 89 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pu ! in : 3 ocean transport components 90 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pv ! out: same 3 transport components 91 REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(inout) :: pw ! increased by the MLE induced transport 88 92 ! 89 93 INTEGER :: ji, jj, jk ! dummy loop indices … … 91 95 REAL(wp) :: zcuw, zmuw, zc ! local scalar 92 96 REAL(wp) :: zcvw, zmvw ! - - 93 INTEGER , DIMENSION(jpi,jpj) :: inml_mle 94 REAL(wp), DIMENSION(jpi,jpj) :: zpsim_u, zpsim_v, zmld, zbm, zhu, zhv, zn2, zLf_NH, zLf_MH 95 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpsi_uw, zpsi_vw 96 !!---------------------------------------------------------------------- 97 ! 98 ! !== MLD used for MLE ==! 99 ! ! compute from the 10m density to deal with the diurnal cycle 100 inml_mle(:,:) = mbkt(:,:) + 1 ! init. to number of ocean w-level (T-level + 1) 101 IF ( nla10 > 0 ) THEN ! avoid case where first level is thicker than 10m 102 DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! from the bottom to nlb10 (10m) 103 IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle ) inml_mle(ji,jj) = jk ! Mixed layer 97 INTEGER , DIMENSION(A2D(nn_hls)) :: inml_mle 98 REAL(wp), DIMENSION(A2D(nn_hls)) :: zpsim_u, zpsim_v, zmld, zbm, zhu, zhv, zn2, zLf_MH 99 REAL(wp), DIMENSION(A2D(nn_hls),jpk) :: zpsi_uw, zpsi_vw 100 ! TEMP: [tiling] These changes not necessary if using XIOS (subdomain support) 101 REAL(wp), DIMENSION(:,:), ALLOCATABLE, SAVE :: zLf_NH 102 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, SAVE :: zpsiu_mle, zpsiv_mle 103 !!---------------------------------------------------------------------- 104 ! 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 104 170 END_3D 105 ENDIF 106 ikmax = MIN( MAXVAL( inml_mle(:,:) ), jpkm1 ) ! max level of the computation 107 ! 108 ! 109 zmld(:,:) = 0._wp !== Horizontal shape of the MLE ==! 110 zbm (:,:) = 0._wp 111 zn2 (:,:) = 0._wp 112 DO_3D( 1, 1, 1, 1, 1, ikmax ) ! MLD and mean buoyancy and N2 over the mixed layer 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 zmld(ji,jj) = zmld(ji,jj) + zc 115 zbm (ji,jj) = zbm (ji,jj) + zc * (rho0 - rhop(ji,jj,jk) ) * r1_rho0 116 zn2 (ji,jj) = zn2 (ji,jj) + zc * (rn2(ji,jj,jk)+rn2(ji,jj,jk+1))*0.5_wp 117 END_3D 118 119 SELECT CASE( nn_mld_uv ) ! MLD at u- & v-pts 120 CASE ( 0 ) != min of the 2 neighbour MLDs 121 DO_2D( 1, 0, 1, 0 ) 122 zhu(ji,jj) = MIN( zmld(ji+1,jj), zmld(ji,jj) ) 123 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) ) 124 192 END_2D 125 CASE ( 1 ) != average of the 2 neighbour MLDs 126 DO_2D( 1, 0, 1, 0 ) 127 zhu(ji,jj) = ( zmld(ji+1,jj) + zmld(ji,jj) ) * 0.5_wp 128 zhv(ji,jj) = ( zmld(ji,jj+1) + zmld(ji,jj) ) * 0.5_wp 129 END_2D 130 CASE ( 2 ) != max of the 2 neighbour MLDs 131 DO_2D( 1, 0, 1, 0 ) 132 zhu(ji,jj) = MAX( zmld(ji+1,jj), zmld(ji,jj) ) 133 zhv(ji,jj) = MAX( zmld(ji,jj+1), zmld(ji,jj) ) 134 END_2D 135 END SELECT 136 ! ! convert density into buoyancy 137 zbm(:,:) = + grav * zbm(:,:) / MAX( e3t(:,:,1,Kmm), zmld(:,:) ) 138 ! 139 ! 140 ! !== Magnitude of the MLE stream function ==! 141 ! 142 ! di[bm] Ds 143 ! Psi = Ce H^2 ---------------- e2u mu(z) where fu Lf = MAX( fu*rn_fl , (Db H)^1/2 ) 144 ! e1u Lf fu and the e2u for the "transport" 145 ! (not *e3u as divided by e3u at the end) 146 ! 147 IF( nn_mle == 0 ) THEN ! Fox-Kemper et al. 2010 formulation 148 DO_2D( 1, 0, 1, 0 ) 149 zpsim_u(ji,jj) = rn_ce * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) & 150 & * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) & 151 & / ( 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) ) ) ) 152 207 ! 153 zpsim_v(ji,jj) = rn_ce * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) &154 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) &155 & / ( MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) ) )156 END_2D157 !158 ELSEIF( nn_mle == 1 ) THEN ! New formulation (Lf = 5km fo/ff with fo=Coriolis parameter at latitude rn_lat)159 DO_2D( 1, 0, 1, 0 )160 zpsim_u(ji,jj) = rc_f * zhu(ji,jj) * zhu(ji,jj) * e2_e1u(ji,jj) &161 & * ( 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) ) 162 217 ! 163 zpsim_v(ji,jj) = rc_f * zhv(ji,jj) * zhv(ji,jj) * e1_e2v(ji,jj) & 164 & * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) 165 END_2D 166 ENDIF 167 ! 168 IF( nn_conv == 1 ) THEN ! No MLE in case of convection 169 DO_2D( 1, 0, 1, 0 ) 170 IF( MIN( zn2(ji,jj) , zn2(ji+1,jj) ) < 0._wp ) zpsim_u(ji,jj) = 0._wp 171 IF( MIN( zn2(ji,jj) , zn2(ji,jj+1) ) < 0._wp ) zpsim_v(ji,jj) = 0._wp 172 END_2D 173 ENDIF 174 ! 175 ! !== structure function value at uw- and vw-points ==! 176 DO_2D( 1, 0, 1, 0 ) 177 zhu(ji,jj) = 1._wp / zhu(ji,jj) ! hu --> 1/hu 178 zhv(ji,jj) = 1._wp / zhv(ji,jj) 179 END_2D 180 ! 181 zpsi_uw(:,:,:) = 0._wp 182 zpsi_vw(:,:,:) = 0._wp 183 ! 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 ! 184 240 DO_3D( 1, 0, 1, 0, 2, ikmax ) ! start from 2 : surface value = 0 185 241 zcuw = 1._wp - ( gdepw(ji+1,jj,jk,Kmm) + gdepw(ji,jj,jk,Kmm) ) * zhu(ji,jj) … … 206 262 END DO 207 263 264 ! TEMP: [tiling] These changes not necessary if using XIOS (subdomain support) 208 265 IF( cdtype == 'TRA') THEN !== outputs ==! 209 ! 210 zLf_NH(:,:) = SQRT( rb_c * zmld(:,:) ) * r1_ft(:,:) ! Lf = N H / f 211 CALL iom_put( "Lf_NHpf" , zLf_NH ) ! Lf = N H / f 266 IF( ntile == 0 .OR. ntile == 1 ) THEN ! Do only on the first tile 267 ALLOCATE( zLf_NH(jpi,jpj), zpsiu_mle(jpi,jpj,jpk), zpsiv_mle(jpi,jpj,jpk) ) 268 zpsiu_mle(:,:,:) = 0._wp ; zpsiv_mle(:,:,:) = 0._wp 269 ENDIF 270 ! 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 212 280 ! 213 281 ! divide by cross distance to give streamfunction with dimensions m^2/s 214 DO jk = 1, ikmax+1 215 zpsi_uw(:,:,jk) = zpsi_uw(:,:,jk) * r1_e2u(:,:) 216 zpsi_vw(:,:,jk) = zpsi_vw(:,:,jk) * r1_e1v(:,:) 217 END DO 218 CALL iom_put( "psiu_mle", zpsi_uw ) ! i-mle streamfunction 219 CALL iom_put( "psiv_mle", zpsi_vw ) ! j-mle streamfunction 282 DO_3D( 0, 0, 0, 0, 1, ikmax+1 ) 283 zpsiu_mle(ji,jj,jk) = zpsi_uw(ji,jj,jk) * r1_e2u(ji,jj) 284 zpsiv_mle(ji,jj,jk) = zpsi_vw(ji,jj,jk) * r1_e1v(ji,jj) 285 END_3D 286 287 IF( ntile == 0 .OR. ntile == nijtile ) THEN ! Do only on the last tile 288 CALL iom_put( "Lf_NHpf" , zLf_NH ) ! Lf = N H / f 289 CALL iom_put( "psiu_mle", zpsiu_mle ) ! i-mle streamfunction 290 CALL iom_put( "psiv_mle", zpsiv_mle ) ! j-mle streamfunction 291 DEALLOCATE( zLf_NH, zpsiu_mle, zpsiv_mle ) 292 ENDIF 220 293 ENDIF 221 294 ! 222 295 END SUBROUTINE tra_mle_trp 223 224 296 225 297 SUBROUTINE tra_mle_init
Note: See TracChangeset
for help on using the changeset viewer.
