- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90
r4328 r6225 27 27 USE dom_oce ! ocean space and time domain variables 28 28 USE sbc_oce ! surface boundary condition: ocean 29 USE zdf_oce ! ??? 29 USE sbcrnf ! river runoffs 30 USE sbcisf ! ice shelf melting 31 USE zdf_oce ! ocean vertical mixing 30 32 USE domvvl ! variable volume 31 USE dynspg_oce ! surface pressure gradient variables 32 USE dynhpg ! hydrostatic pressure gradient 33 USE trdmod_oce ! ocean space and time domain variables 34 USE trdtra ! ocean active tracers trends 35 USE phycst 36 USE bdy_oce 33 USE trd_oce ! trends: ocean variables 34 USE trdtra ! trends manager: tracers 35 USE traqsr ! penetrative solar radiation (needed for nksr) 36 USE phycst ! physical constant 37 USE ldftra ! lateral physics on tracers 38 USE ldfslp 39 USE bdy_oce ! BDY open boundary condition variables 37 40 USE bdytra ! open boundary condition (bdy_tra routine) 41 ! 38 42 USE in_out_manager ! I/O manager 39 43 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 40 44 USE prtctl ! Print control 41 USE traqsr ! penetrative solar radiation (needed for nksr) 45 USE wrk_nemo ! Memory allocation 46 USE timing ! Timing 42 47 #if defined key_agrif 43 USE agrif_opa_update44 48 USE agrif_opa_interp 45 49 #endif 46 USE wrk_nemo ! Memory allocation47 USE timing ! Timing48 50 49 51 IMPLICIT NONE … … 54 56 PUBLIC tra_nxt_vvl ! to be used in trcnxt 55 57 56 REAL(wp) :: rbcp ! Brown & Campana parameters for semi-implicit hpg57 58 58 !! * Substitutions 59 # include " domzgr_substitute.h90"59 # include "vectopt_loop_substitute.h90" 60 60 !!---------------------------------------------------------------------- 61 61 !! NEMO/OPA 3.3 , NEMO-Consortium (2010) … … 80 80 !! at the local domain boundaries through lbc_lnk call, 81 81 !! at the one-way open boundaries (lk_bdy=T), 82 !! at the AGRIF zoom 82 !! at the AGRIF zoom boundaries (lk_agrif=T) 83 83 !! 84 84 !! - Update lateral boundary conditions on AGRIF children 85 85 !! domains (lk_agrif=T) 86 86 !! 87 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step 88 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 87 !! ** Action : - tsb & tsn ready for the next time step 89 88 !!---------------------------------------------------------------------- 90 89 INTEGER, INTENT(in) :: kt ! ocean time-step index 91 90 !! 92 INTEGER :: j k, jn! dummy loop indices93 REAL(wp) :: zfact ! local scalars91 INTEGER :: ji, jj, jk, jn ! dummy loop indices 92 REAL(wp) :: zfact ! local scalars 94 93 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdt, ztrds 95 94 !!---------------------------------------------------------------------- … … 101 100 IF(lwp) WRITE(numout,*) 'tra_nxt : achieve the time stepping by Asselin filter and array swap' 102 101 IF(lwp) WRITE(numout,*) '~~~~~~~' 103 !104 rbcp = 0.25_wp * (1._wp + atfp) * (1._wp + atfp) * ( 1._wp - atfp) ! Brown & Campana parameter for semi-implicit hpg105 102 ENDIF 106 103 107 104 ! Update after tracer on domain lateral boundaries 108 105 ! 106 #if defined key_agrif 107 CALL Agrif_tra ! AGRIF zoom boundaries 108 #endif 109 ! 109 110 CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1._wp ) ! local domain boundaries (T-point, unchanged sign) 110 111 CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1._wp ) … … 113 114 IF( lk_bdy ) CALL bdy_tra( kt ) ! BDY open boundaries 114 115 #endif 115 #if defined key_agrif116 CALL Agrif_tra ! AGRIF zoom boundaries117 #endif118 116 119 117 ! set time step size (Euler/Leapfrog) 120 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt ra(:) = rdttra(:)! at nit000 (Euler)121 ELSEIF( kt <= nit000 + 1 ) THEN ; r2dt ra(:) = 2._wp* rdttra(:)! at nit000 or nit000+1 (Leapfrog)118 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt = rdt ! at nit000 (Euler) 119 ELSEIF( kt <= nit000 + 1 ) THEN ; r2dt = 2._wp* rdt ! at nit000 or nit000+1 (Leapfrog) 122 120 ENDIF 123 121 … … 127 125 ztrdt(:,:,:) = tsn(:,:,:,jp_tem) 128 126 ztrds(:,:,:) = tsn(:,:,:,jp_sal) 127 IF( ln_traldf_iso ) THEN ! diagnose the "pure" Kz diffusive trend 128 CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdfp, ztrdt ) 129 CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdfp, ztrds ) 130 ENDIF 129 131 ENDIF 130 132 … … 135 137 END DO 136 138 END DO 139 ! 137 140 ELSE ! Leap-Frog + Asselin filter time stepping 138 141 ! 139 IF( lk_vvl ) THEN ; CALL tra_nxt_vvl( kt, nit000, 'TRA', tsb, tsn, tsa, jpts ) ! variable volume level (vvl) 140 ELSE ; CALL tra_nxt_fix( kt, nit000, 'TRA', tsb, tsn, tsa, jpts ) ! fixed volume level 142 IF( ln_linssh ) THEN ; CALL tra_nxt_fix( kt, nit000, 'TRA', tsb, tsn, tsa, jpts ) ! linear free surface 143 ELSE ; CALL tra_nxt_vvl( kt, nit000, rdt, 'TRA', tsb, tsn, tsa, & 144 & sbc_tsc, sbc_tsc_b, jpts ) ! non-linear free surface 141 145 ENDIF 142 ENDIF143 !144 #if defined key_agrif 145 ! Update tracer at AGRIF zoom boundaries146 IF( .NOT.Agrif_Root() ) CALL Agrif_Update_Tra( kt ) ! children only147 #endif 148 !149 ! trends computation146 ! 147 DO jn = 1, jpts 148 CALL lbc_lnk( tsb(:,:,:,jn), 'T', 1._wp ) 149 CALL lbc_lnk( tsn(:,:,:,jn), 'T', 1._wp ) 150 CALL lbc_lnk( tsa(:,:,:,jn), 'T', 1._wp ) 151 END DO 152 ENDIF 153 ! 150 154 IF( l_trdtra ) THEN ! trend of the Asselin filter (tb filtered - tb)/dt 151 155 DO jk = 1, jpkm1 152 zfact = 1. e0_wp / r2dtra(jk)156 zfact = 1._wp / r2dt 153 157 ztrdt(:,:,jk) = ( tsb(:,:,jk,jp_tem) - ztrdt(:,:,jk) ) * zfact 154 158 ztrds(:,:,jk) = ( tsb(:,:,jk,jp_sal) - ztrds(:,:,jk) ) * zfact 155 159 END DO 156 CALL trd_tra( kt, 'TRA', jp_tem, jptra_ trd_atf, ztrdt )157 CALL trd_tra( kt, 'TRA', jp_sal, jptra_ trd_atf, ztrds )160 CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt ) 161 CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds ) 158 162 CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 159 163 END IF … … 163 167 & tab3d_2=tsn(:,:,:,jp_sal), clinfo2= ' Sn: ', mask2=tmask ) 164 168 ! 165 ! 166 IF( nn_timing == 1 ) CALL timing_stop('tra_nxt') 169 IF( nn_timing == 1 ) CALL timing_stop('tra_nxt') 167 170 ! 168 171 END SUBROUTINE tra_nxt … … 177 180 !! 178 181 !! ** Method : - Apply a Asselin time filter on now fields. 179 !! - save in (ta,sa) an average over the three time levels180 !! which will be used to compute rdn and thus the semi-implicit181 !! hydrostatic pressure gradient (ln_dynhpg_imp = T)182 182 !! - swap tracer fields to prepare the next time_step. 183 !! This can be summurized for tempearture as: 184 !! ztm = tn + rbcp * [ta -2 tn + tb ] ln_dynhpg_imp = T 185 !! ztm = 0 otherwise 186 !! with rbcp=1/4 * (1-atfp^4) / (1-atfp) 187 !! tb = tn + atfp*[ tb - 2 tn + ta ] 188 !! tn = ta 189 !! ta = ztm (NB: reset to 0 after eos_bn2 call) 190 !! 191 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step 192 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 193 !!---------------------------------------------------------------------- 194 INTEGER , INTENT(in ) :: kt ! ocean time-step index 195 INTEGER , INTENT(in ) :: kit000 ! first time step index 196 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 197 INTEGER , INTENT(in ) :: kjpt ! number of tracers 198 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: ptb ! before tracer fields 199 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: ptn ! now tracer fields 200 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: pta ! tracer trend 183 !! 184 !! ** Action : - tsb & tsn ready for the next time step 185 !!---------------------------------------------------------------------- 186 INTEGER , INTENT(in ) :: kt ! ocean time-step index 187 INTEGER , INTENT(in ) :: kit000 ! first time step index 188 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 189 INTEGER , INTENT(in ) :: kjpt ! number of tracers 190 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptb ! before tracer fields 191 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptn ! now tracer fields 192 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 201 193 ! 202 194 INTEGER :: ji, jj, jk, jn ! dummy loop indices 203 LOGICAL :: ll_tra_hpg ! local logical204 195 REAL(wp) :: ztn, ztd ! local scalars 205 196 !!---------------------------------------------------------------------- 206 197 ! 207 198 IF( kt == kit000 ) THEN 208 199 IF(lwp) WRITE(numout,*) … … 211 202 ENDIF 212 203 ! 213 IF( cdtype == 'TRA' ) THEN ; ll_tra_hpg = ln_dynhpg_imp ! active tracers case and semi-implicit hpg214 ELSE ; ll_tra_hpg = .FALSE. ! passive tracers case or NO semi-implicit hpg215 ENDIF216 !217 204 DO jn = 1, kjpt 218 205 ! 219 206 DO jk = 1, jpkm1 220 DO jj = 1, jpj221 DO ji = 1, jpi207 DO jj = 2, jpjm1 208 DO ji = fs_2, fs_jpim1 222 209 ztn = ptn(ji,jj,jk,jn) 223 ztd = pta(ji,jj,jk,jn) - 2. * ztn + ptb(ji,jj,jk,jn) ! time laplacian on tracers 224 ! 225 ptb(ji,jj,jk,jn) = ztn + atfp * ztd ! ptb <-- filtered ptn 226 ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn) ! ptn <-- pta 227 ! 228 IF( ll_tra_hpg ) pta(ji,jj,jk,jn) = ztn + rbcp * ztd ! pta <-- Brown & Campana average 210 ztd = pta(ji,jj,jk,jn) - 2._wp * ztn + ptb(ji,jj,jk,jn) ! time laplacian on tracers 211 ! 212 ptb(ji,jj,jk,jn) = ztn + atfp * ztd ! ptb <-- filtered ptn 213 ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn) ! ptn <-- pta 229 214 END DO 230 215 END DO … … 236 221 237 222 238 SUBROUTINE tra_nxt_vvl( kt, kit000, cdtype, ptb, ptn, pta, kjpt )223 SUBROUTINE tra_nxt_vvl( kt, kit000, p2dt, cdtype, ptb, ptn, pta, psbc_tc, psbc_tc_b, kjpt ) 239 224 !!---------------------------------------------------------------------- 240 225 !! *** ROUTINE tra_nxt_vvl *** … … 244 229 !! 245 230 !! ** Method : - Apply a thickness weighted Asselin time filter on now fields. 246 !! - save in (ta,sa) a thickness weighted average over the three247 !! time levels which will be used to compute rdn and thus the semi-248 !! implicit hydrostatic pressure gradient (ln_dynhpg_imp = T)249 231 !! - swap tracer fields to prepare the next time_step. 250 !! This can be summurized for tempearture as:251 !! ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) ln_dynhpg_imp = T252 !! /( e3t_n + rbcp*[ e3t_b - 2 e3t_n + e3t_a ] )253 !! ztm = 0 otherwise254 232 !! tb = ( e3t_n*tn + atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) 255 233 !! /( e3t_n + atfp*[ e3t_b - 2 e3t_n + e3t_a ] ) 256 234 !! tn = ta 257 !! ta = zt (NB: reset to 0 after eos_bn2 call) 258 !! 259 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step 260 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 261 !!---------------------------------------------------------------------- 262 INTEGER , INTENT(in ) :: kt ! ocean time-step index 263 INTEGER , INTENT(in ) :: kit000 ! first time step index 264 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 265 INTEGER , INTENT(in ) :: kjpt ! number of tracers 266 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: ptb ! before tracer fields 267 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: ptn ! now tracer fields 268 REAL(wp) , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt) :: pta ! tracer trend 269 !! 270 LOGICAL :: ll_tra, ll_tra_hpg, ll_traqsr ! local logical 235 !! 236 !! ** Action : - tsb & tsn ready for the next time step 237 !!---------------------------------------------------------------------- 238 INTEGER , INTENT(in ) :: kt ! ocean time-step index 239 INTEGER , INTENT(in ) :: kit000 ! first time step index 240 REAL(wp) , INTENT(in ) :: p2dt ! time-step 241 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 242 INTEGER , INTENT(in ) :: kjpt ! number of tracers 243 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptb ! before tracer fields 244 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptn ! now tracer fields 245 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 246 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: psbc_tc ! surface tracer content 247 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: psbc_tc_b ! before surface tracer content 248 ! 249 LOGICAL :: ll_traqsr, ll_rnf, ll_isf ! local logical 271 250 INTEGER :: ji, jj, jk, jn ! dummy loop indices 272 251 REAL(wp) :: zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d ! local scalar … … 281 260 ! 282 261 IF( cdtype == 'TRA' ) THEN 283 ll_tra = .TRUE. ! active tracers case284 ll_tra_hpg = ln_dynhpg_imp ! active tracers case and semi-implicit hpg285 262 ll_traqsr = ln_traqsr ! active tracers case and solar penetration 286 ELSE 287 ll_tra = .FALSE. ! passive tracers case 288 ll_tra_hpg = .FALSE. ! passive tracers case or NO semi-implicit hpg 289 ll_traqsr = .FALSE. ! active tracers case and NO solar penetration 263 ll_rnf = ln_rnf ! active tracers case and river runoffs 264 ll_isf = ln_isf ! active tracers case and ice shelf melting 265 ELSE ! passive tracers case 266 ll_traqsr = .FALSE. ! NO solar penetration 267 ll_rnf = .FALSE. ! NO river runoffs ???? !!gm BUG ? 268 ll_isf = .FALSE. ! NO ice shelf melting/freezing !!gm BUG ?? 290 269 ENDIF 291 270 ! 292 271 DO jn = 1, kjpt 293 272 DO jk = 1, jpkm1 294 zfact1 = atfp * rdttra(jk)295 zfact2 = zfact1 /rau0296 DO jj = 1, jpj297 DO ji = 1, jpi298 ze3t_b = fse3t_b(ji,jj,jk)299 ze3t_n = fse3t_n(ji,jj,jk)300 ze3t_a = fse3t_a(ji,jj,jk)273 zfact1 = atfp * p2dt 274 zfact2 = zfact1 * r1_rau0 275 DO jj = 2, jpjm1 276 DO ji = fs_2, fs_jpim1 277 ze3t_b = e3t_b(ji,jj,jk) 278 ze3t_n = e3t_n(ji,jj,jk) 279 ze3t_a = e3t_a(ji,jj,jk) 301 280 ! ! tracer content at Before, now and after 302 281 ztc_b = ptb(ji,jj,jk,jn) * ze3t_b … … 310 289 ztc_f = ztc_n + atfp * ztc_d 311 290 ! 312 IF( ll_tra .AND. jk == 1 ) THEN ! first level only for T & S 313 ze3t_f = ze3t_f - zfact2 * ( emp_b(ji,jj) - emp(ji,jj) ) 314 ztc_f = ztc_f - zfact1 * ( sbc_tsc(ji,jj,jn) - sbc_tsc_b(ji,jj,jn) ) 291 IF( jk == mikt(ji,jj) ) THEN ! first level 292 ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj) - emp(ji,jj) ) & 293 & - (rnf_b(ji,jj) - rnf(ji,jj) ) & 294 & + (fwfisf_b(ji,jj) - fwfisf(ji,jj)) ) 295 ztc_f = ztc_f - zfact1 * ( psbc_tc(ji,jj,jn) - psbc_tc_b(ji,jj,jn) ) 315 296 ENDIF 316 IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr ) & ! solar penetration (temperature only) 297 ! 298 ! solar penetration (temperature only) 299 IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr ) & 317 300 & ztc_f = ztc_f - zfact1 * ( qsr_hc(ji,jj,jk) - qsr_hc_b(ji,jj,jk) ) 318 319 ze3t_f = 1.e0 / ze3t_f 320 ptb(ji,jj,jk,jn) = ztc_f * ze3t_f ! ptb <-- ptn filtered 321 ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn) ! ptn <-- pta 322 ! 323 IF( ll_tra_hpg ) THEN ! semi-implicit hpg (T & S only) 324 ze3t_d = 1.e0 / ( ze3t_n + rbcp * ze3t_d ) 325 pta(ji,jj,jk,jn) = ze3t_d * ( ztc_n + rbcp * ztc_d ) ! ta <-- Brown & Campana average 326 ENDIF 301 ! 302 ! river runoff 303 IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) ) & 304 & ztc_f = ztc_f - zfact1 * ( rnf_tsc(ji,jj,jn) - rnf_tsc_b(ji,jj,jn) ) & 305 & * e3t_n(ji,jj,jk) / h_rnf(ji,jj) 306 ! 307 ! ice shelf 308 IF( ll_isf ) THEN 309 ! level fully include in the Losch_2008 ice shelf boundary layer 310 IF ( jk >= misfkt(ji,jj) .AND. jk < misfkb(ji,jj) ) & 311 ztc_f = ztc_f - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) ) & 312 & * e3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) 313 ! level partially include in Losch_2008 ice shelf boundary layer 314 IF ( jk == misfkb(ji,jj) ) & 315 ztc_f = ztc_f - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) ) & 316 & * e3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) * ralpha(ji,jj) 317 END IF 318 ! 319 ze3t_f = 1.e0 / ze3t_f 320 ptb(ji,jj,jk,jn) = ztc_f * ze3t_f ! ptb <-- ptn filtered 321 ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn) ! ptn <-- pta 322 ! 327 323 END DO 328 324 END DO
Note: See TracChangeset
for help on using the changeset viewer.
