- Timestamp:
- 2019-11-22T15:29:17+01:00 (4 years ago)
- Location:
- NEMO/branches/2019/dev_r11943_MERGE_2019/tests
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
CANAL/MY_SRC/trazdf.F90 (modified) (13 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11943_MERGE_2019/tests
- Property svn:mergeinfo deleted
-
NEMO/branches/2019/dev_r11943_MERGE_2019/tests/CANAL/MY_SRC/trazdf.F90
r10572 r11949 44 44 CONTAINS 45 45 46 SUBROUTINE tra_zdf( kt )46 SUBROUTINE tra_zdf( kt, Kbb, Kmm, Krhs, pts, Kaa ) 47 47 !!---------------------------------------------------------------------- 48 48 !! *** ROUTINE tra_zdf *** … … 50 50 !! ** Purpose : compute the vertical ocean tracer physics. 51 51 !!--------------------------------------------------------------------- 52 INTEGER, INTENT(in) :: kt ! ocean time-step index 52 INTEGER , INTENT(in) :: kt ! ocean time-step index 53 INTEGER , INTENT(in) :: Kbb, Kmm, Krhs, Kaa ! time level indices 54 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts ! active tracers and RHS of tracer equation 53 55 ! 54 56 INTEGER :: jk ! Dummy loop indices … … 70 72 IF( l_trdtra ) THEN !* Save ta and sa trends 71 73 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 72 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)73 ztrds(:,:,:) = tsa(:,:,:,jp_sal)74 ztrdt(:,:,:) = pts(:,:,:,jp_tem,Kaa) 75 ztrds(:,:,:) = pts(:,:,:,jp_sal,Kaa) 74 76 ENDIF 75 77 ! 76 78 ! !* compute lateral mixing trend and add it to the general trend 77 CALL tra_zdf_imp( kt, nit000, 'TRA', r2dt, tsb, tsa, jpts )79 CALL tra_zdf_imp( kt, nit000, 'TRA', r2dt, Kbb, Kmm, Krhs, pts, Kaa, jpts ) 78 80 79 81 !!gm WHY here ! and I don't like that ! … … 81 83 ! JMM avoid negative salinities near river outlet ! Ugly fix 82 84 ! JMM : restore negative salinities to small salinities: 83 !!$ WHERE( tsa(:,:,:,jp_sal) < 0._wp ) tsa(:,:,:,jp_sal) = 0.1_wp85 !!$ WHERE( pts(:,:,:,jp_sal,Kaa) < 0._wp ) pts(:,:,:,jp_sal,Kaa) = 0.1_wp 84 86 !!gm 85 87 86 88 IF( l_trdtra ) THEN ! save the vertical diffusive trends for further diagnostics 87 89 DO jk = 1, jpkm1 88 ztrdt(:,:,jk) = ( ( tsa(:,:,jk,jp_tem)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_tem)*e3t_b(:,:,jk) ) &89 & / (e3t _n(:,:,jk)*r2dt) ) - ztrdt(:,:,jk)90 ztrds(:,:,jk) = ( ( tsa(:,:,jk,jp_sal)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_sal)*e3t_b(:,:,jk) ) &91 & / (e3t _n(:,:,jk)*r2dt) ) - ztrds(:,:,jk)90 ztrdt(:,:,jk) = ( ( pts(:,:,jk,jp_tem,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_tem,Kbb)*e3t(:,:,jk,Kbb) ) & 91 & / (e3t(:,:,jk,Kmm)*r2dt) ) - ztrdt(:,:,jk) 92 ztrds(:,:,jk) = ( ( pts(:,:,jk,jp_sal,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_sal,Kbb)*e3t(:,:,jk,Kbb) ) & 93 & / (e3t(:,:,jk,Kmm)*r2dt) ) - ztrds(:,:,jk) 92 94 END DO 93 95 !!gm this should be moved in trdtra.F90 and done on all trends 94 96 CALL lbc_lnk_multi( 'trazdf', ztrdt, 'T', 1. , ztrds, 'T', 1. ) 95 97 !!gm 96 CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdf, ztrdt )97 CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdf, ztrds )98 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_zdf, ztrdt ) 99 CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_zdf, ztrds ) 98 100 DEALLOCATE( ztrdt , ztrds ) 99 101 ENDIF 100 102 ! ! print mean trends (used for debugging) 101 IF(ln_ctl) CALL prt_ctl( tab3d_1= tsa(:,:,:,jp_tem), clinfo1=' zdf - Ta: ', mask1=tmask, &102 & tab3d_2= tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )103 IF(ln_ctl) CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Kaa), clinfo1=' zdf - Ta: ', mask1=tmask, & 104 & tab3d_2=pts(:,:,:,jp_sal,Kaa), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' ) 103 105 ! 104 106 IF( ln_timing ) CALL timing_stop('tra_zdf') … … 107 109 108 110 109 SUBROUTINE tra_zdf_imp( kt, kit000, cdtype, p2dt, ptb, pta, kjpt )111 SUBROUTINE tra_zdf_imp( kt, kit000, cdtype, p2dt, Kbb, Kmm, Krhs, pt, Kaa, kjpt ) 110 112 !!---------------------------------------------------------------------- 111 113 !! *** ROUTINE tra_zdf_imp *** … … 125 127 !! If iso-neutral mixing, add to avt the contribution due to lateral mixing. 126 128 !! 127 !! ** Action : - pt abecomes the after tracer128 !!--------------------------------------------------------------------- 129 INTEGER , INTENT(in ) :: kt ! ocean time-step index130 INTEGER , INTENT(in ) :: kit000 ! first time step index131 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)132 INTEGER , INTENT(in ) :: kjpt ! number of tracers133 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step134 REAL(wp) , DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields135 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt ), INTENT(inout) :: pta ! in: tracer trend ; out: after tracer field129 !! ** Action : - pt(:,:,:,:,Kaa) becomes the after tracer 130 !!--------------------------------------------------------------------- 131 INTEGER , INTENT(in ) :: kt ! ocean time-step index 132 INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs, Kaa ! ocean time level indices 133 INTEGER , INTENT(in ) :: kit000 ! first time step index 134 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 135 INTEGER , INTENT(in ) :: kjpt ! number of tracers 136 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 137 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) :: pt ! tracers and RHS of tracer equation 136 138 ! 137 139 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 181 183 DO jj = 2, jpjm1 182 184 DO ji = fs_2, fs_jpim1 ! vector opt. (ensure same order of calculation as below if wi=0.) 183 zzwi = - p2dt * zwt(ji,jj,jk ) / e3w _n(ji,jj,jk)184 zzws = - p2dt * zwt(ji,jj,jk+1) / e3w _n(ji,jj,jk+1)185 zwd(ji,jj,jk) = e3t _a(ji,jj,jk) - zzwi - zzws &185 zzwi = - p2dt * zwt(ji,jj,jk ) / e3w(ji,jj,jk ,Kmm) 186 zzws = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) 187 zwd(ji,jj,jk) = e3t(ji,jj,jk,Kaa) - zzwi - zzws & 186 188 & + p2dt * ( MAX( wi(ji,jj,jk ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) 187 189 zwi(ji,jj,jk) = zzwi + p2dt * MIN( wi(ji,jj,jk ) , 0._wp ) … … 194 196 DO jj = 2, jpjm1 195 197 DO ji = fs_2, fs_jpim1 ! vector opt. 196 zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk ) / e3w _n(ji,jj,jk)197 zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w _n(ji,jj,jk+1)198 zwd(ji,jj,jk) = e3t _a(ji,jj,jk) - zwi(ji,jj,jk) - zws(ji,jj,jk)198 zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk ) / e3w(ji,jj,jk,Kmm) 199 zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) 200 zwd(ji,jj,jk) = e3t(ji,jj,jk,Kaa) - zwi(ji,jj,jk) - zws(ji,jj,jk) 199 201 END DO 200 202 END DO … … 218 220 ! The solution will be in the 4d array pta. 219 221 ! The 3d array zwt is used as a work space array. 220 ! En route to the solution pt ais used a to evaluate the rhs and then222 ! En route to the solution pt(:,:,:,:,Kaa) is used a to evaluate the rhs and then 221 223 ! used as a work space array: its value is modified. 222 224 ! … … 238 240 DO jj = 2, jpjm1 !* 2nd recurrence: Zk = Yk - Ik / Tk-1 Zk-1 239 241 DO ji = fs_2, fs_jpim1 240 pt a(ji,jj,1,jn) = e3t_b(ji,jj,1) * ptb(ji,jj,1,jn) + p2dt * e3t_n(ji,jj,1) * pta(ji,jj,1,jn)242 pt(ji,jj,1,jn,Kaa) = e3t(ji,jj,1,Kbb) * pt(ji,jj,1,jn,Kbb) + p2dt * e3t(ji,jj,1,Kmm) * pt(ji,jj,1,jn,Krhs) 241 243 END DO 242 244 END DO … … 244 246 DO jj = 2, jpjm1 245 247 DO ji = fs_2, fs_jpim1 246 zrhs = e3t _b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * e3t_n(ji,jj,jk) * pta(ji,jj,jk,jn) ! zrhs=right hand side247 pt a(ji,jj,jk,jn) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pta(ji,jj,jk-1,jn)248 zrhs = e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * e3t(ji,jj,jk,Kmm) * pt(ji,jj,jk,jn,Krhs) ! zrhs=right hand side 249 pt(ji,jj,jk,jn,Kaa) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pt(ji,jj,jk-1,jn,Kaa) 248 250 END DO 249 251 END DO … … 252 254 DO jj = 2, jpjm1 !* 3d recurrence: Xk = (Zk - Sk Xk+1 ) / Tk (result is the after tracer) 253 255 DO ji = fs_2, fs_jpim1 254 pt a(ji,jj,jpkm1,jn) = pta(ji,jj,jpkm1,jn) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1)256 pt(ji,jj,jpkm1,jn,Kaa) = pt(ji,jj,jpkm1,jn,Kaa) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1) 255 257 END DO 256 258 END DO … … 258 260 DO jj = 2, jpjm1 259 261 DO ji = fs_2, fs_jpim1 260 pt a(ji,jj,jk,jn) = ( pta(ji,jj,jk,jn) - zws(ji,jj,jk) * pta(ji,jj,jk+1,jn) ) &262 pt(ji,jj,jk,jn,Kaa) = ( pt(ji,jj,jk,jn,Kaa) - zws(ji,jj,jk) * pt(ji,jj,jk+1,jn,Kaa) ) & 261 263 & / zwt(ji,jj,jk) * tmask(ji,jj,jk) 262 264 END DO
Note: See TracChangeset
for help on using the changeset viewer.
