Changeset 1361 for branches/dev_004_VVL/NEMO/OPA_SRC/TRA/tranxt.F90
- Timestamp:
- 2009-03-31T18:26:41+02:00 (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_004_VVL/NEMO/OPA_SRC/TRA/tranxt.F90
r1146 r1361 39 39 PUBLIC tra_nxt ! routine called by step.F90 40 40 41 REAL(wp), DIMENSION(jpk) :: r2dt_t ! vertical profile time step, =2*rdttra (leapfrog) or =rdttra (Euler) 42 41 43 !! * Substitutions 42 44 # include "domzgr_substitute.h90" … … 67 69 !! at the AGRIF zoom boundaries (lk_agrif=T) 68 70 !! 69 !! - Apply the Asselin time filter on now fields, 70 !! save in (ta,sa) an average over the three time levels 71 !! which will be used to compute rdn and thus the semi-implicit 72 !! hydrostatic pressure gradient (ln_dynhpg_imp = T), and 73 !! swap tracer fields to prepare the next time_step. 74 !! This can be summurized for tempearture as: 75 !! zt = (ta+2tn+tb)/4 ln_dynhpg_imp = T 76 !! zt = 0 otherwise 77 !! tb = tn + atfp*[ tb - 2 tn + ta ] 78 !! tn = ta 79 !! ta = zt (NB: reset to 0 after eos_bn2 call) 71 !! - Update lateral boundary conditions on AGRIF children 72 !! domains (lk_agrif=T) 80 73 !! 81 74 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step … … 87 80 INTEGER, INTENT(in) :: kt ! ocean time-step index 88 81 !! 89 INTEGER :: j i, jj, jk ! dummy loop indices90 REAL(wp) :: z t, zs, zfact ! temporary scalars82 INTEGER :: jk ! dummy loop indices 83 REAL(wp) :: zfact ! temporary scalars 91 84 !!---------------------------------------------------------------------- 92 85 … … 111 104 CALL Agrif_tra ! AGRIF zoom boundaries 112 105 #endif 106 107 ! set time step size (Euler/Leapfrog) 108 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt_t(:) = rdttra(:) ! at nit000 (Euler) 109 ELSEIF( kt <= nit000 + 1 ) THEN ; r2dt_t(:) = 2.* rdttra(:) ! at nit000 or nit000+1 (Leapfrog) 110 ENDIF 113 111 114 112 ! trends computation initialisation … … 118 116 ENDIF 119 117 120 ! Asselin time filter and swap of arrays 121 ! 122 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler 1st time step : swap only 123 DO jk = 1, jpkm1 124 tb(:,:,jk) = tn(:,:,jk) ! ta, sa remain at their values which 125 sb(:,:,jk) = sn(:,:,jk) ! correspond to tn, sn after the sawp 126 tn(:,:,jk) = ta(:,:,jk) 127 sn(:,:,jk) = sa(:,:,jk) 128 END DO 129 ! 130 ELSE ! Leap-Frog : filter + swap 131 ! 132 IF( ln_dynhpg_imp ) THEN ! semi-implicite hpg case 133 DO jk = 1, jpkm1 ! (save the averaged of the 3 time steps 134 DO jj = 1, jpj ! in the after fields) 135 DO ji = 1, jpi 136 zt = ( ta(ji,jj,jk) + 2. * tn(ji,jj,jk) + tb(ji,jj,jk) ) * 0.25 137 zs = ( sa(ji,jj,jk) + 2. * sn(ji,jj,jk) + sb(ji,jj,jk) ) * 0.25 138 tb(ji,jj,jk) = atfp * ( tb(ji,jj,jk) + ta(ji,jj,jk) ) + atfp1 * tn(ji,jj,jk) 139 sb(ji,jj,jk) = atfp * ( sb(ji,jj,jk) + sa(ji,jj,jk) ) + atfp1 * sn(ji,jj,jk) 140 tn(ji,jj,jk) = ta(ji,jj,jk) 141 sn(ji,jj,jk) = sa(ji,jj,jk) 142 ta(ji,jj,jk) = zt 143 sa(ji,jj,jk) = zs 144 END DO 145 END DO 146 END DO 147 ELSE ! explicit hpg case 148 DO jk = 1, jpkm1 149 DO jj = 1, jpj 150 DO ji = 1, jpi 151 tb(ji,jj,jk) = atfp * ( tb(ji,jj,jk) + ta(ji,jj,jk) ) + atfp1 * tn(ji,jj,jk) 152 sb(ji,jj,jk) = atfp * ( sb(ji,jj,jk) + sa(ji,jj,jk) ) + atfp1 * sn(ji,jj,jk) 153 tn(ji,jj,jk) = ta(ji,jj,jk) 154 sn(ji,jj,jk) = sa(ji,jj,jk) 155 END DO 156 END DO 157 END DO 158 ENDIF 159 ! 118 ! Leap-Frog + Asselin filter time stepping 119 IF( lk_vvl ) THEN ; CALL tra_nxt_vvl( kt ) ! variable volume level (vvl) 120 ELSE ; CALL tra_nxt_fix( kt ) ! fixed volume level 160 121 ENDIF 161 122 … … 168 129 IF( l_trdtra ) THEN 169 130 DO jk = 1, jpkm1 170 zfact = 1.e0 / ( 2.*rdttra(jk) ) ! NB: euler case, (tb filtered - tb)=0 so 2dt always OK131 zfact = 1.e0 / r2dt_t(jk) 171 132 ztrdt(:,:,jk) = ( tb(:,:,jk) - ztrdt(:,:,jk) ) * zfact 172 133 ztrds(:,:,jk) = ( sb(:,:,jk) - ztrds(:,:,jk) ) * zfact … … 180 141 END SUBROUTINE tra_nxt 181 142 143 144 SUBROUTINE tra_nxt_fix( kt ) 145 !!---------------------------------------------------------------------- 146 !! *** ROUTINE tra_nxt_fix *** 147 !! 148 !! ** Purpose : fixed volume: apply the Asselin time filter and 149 !! swap the tracer fields. 150 !! 151 !! ** Method : - Apply a Asselin time filter on now fields. 152 !! - save in (ta,sa) an average over the three time levels 153 !! which will be used to compute rdn and thus the semi-implicit 154 !! hydrostatic pressure gradient (ln_dynhpg_imp = T) 155 !! - swap tracer fields to prepare the next time_step. 156 !! This can be summurized for tempearture as: 157 !! ztm = (ta+2tn+tb)/4 ln_dynhpg_imp = T 158 !! ztm = 0 otherwise 159 !! tb = tn + atfp*[ tb - 2 tn + ta ] 160 !! tn = ta 161 !! ta = ztm (NB: reset to 0 after eos_bn2 call) 162 !! 163 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step 164 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 165 !!---------------------------------------------------------------------- 166 INTEGER, INTENT(in) :: kt ! ocean time-step index 167 !! 168 INTEGER :: ji, jj, jk ! dummy loop indices 169 REAL(wp) :: ztm, ztf ! temporary scalars 170 REAL(wp) :: zsm, zsf ! - - 171 !!---------------------------------------------------------------------- 172 173 IF( kt == nit000 ) THEN 174 IF(lwp) WRITE(numout,*) 175 IF(lwp) WRITE(numout,*) 'tra_nxt_fix : time stepping' 176 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 177 ENDIF 178 ! 179 ! ! ----------------------- ! 180 IF( ln_dynhpg_imp ) THEN ! semi-implicite hpg case ! 181 ! ! ----------------------- ! 182 ! 183 IF( neuler == 0 .AND. kt == nit000 ) THEN ! case of Euler time-stepping at first time-step 184 DO jk = 1, jpkm1 185 DO jj = 1, jpj 186 DO ji = 1, jpi 187 ztm = 0.25 * ( ta(ji,jj,jk) + 2. * tn(ji,jj,jk) + tb(ji,jj,jk) ) ! mean t 188 zsm = 0.25 * ( sa(ji,jj,jk) + 2. * sn(ji,jj,jk) + sb(ji,jj,jk) ) 189 tb(ji,jj,jk) = tn(ji,jj,jk) ! tb <-- tn 190 sb(ji,jj,jk) = sn(ji,jj,jk) 191 tn(ji,jj,jk) = ta(ji,jj,jk) ! tb <-- tn 192 sn(ji,jj,jk) = sa(ji,jj,jk) 193 ta(ji,jj,jk) = ztm ! ta <-- mean t 194 sa(ji,jj,jk) = zsm 195 END DO 196 END DO 197 END DO 198 ELSE ! general case (Leapfrog + Asselin filter 199 DO jk = 1, jpkm1 200 DO jj = 1, jpj 201 DO ji = 1, jpi 202 ztm = 0.25 * ( ta(ji,jj,jk) + 2.* tn(ji,jj,jk) + tb(ji,jj,jk) ) ! mean t 203 zsm = 0.25 * ( sa(ji,jj,jk) + 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 204 ztf = atfp * ( ta(ji,jj,jk) - 2.* tn(ji,jj,jk) + tb(ji,jj,jk) ) ! Asselin filter on t 205 zsf = atfp * ( sa(ji,jj,jk) - 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 206 tb(ji,jj,jk) = tn(ji,jj,jk) + ztf ! tb <-- filtered tn 207 sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 208 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 209 sn(ji,jj,jk) = sa(ji,jj,jk) 210 ta(ji,jj,jk) = ztm ! ta <-- mean t 211 sa(ji,jj,jk) = zsm 212 END DO 213 END DO 214 END DO 215 ENDIF 216 ! ! ----------------------- ! 217 ELSE ! explicit hpg case ! 218 ! ! ----------------------- ! 219 ! 220 IF( neuler == 0 .AND. kt == nit000 ) THEN ! case of Euler time-stepping at first time-step 221 DO jk = 1, jpkm1 222 DO jj = 1, jpj 223 DO ji = 1, jpi 224 tb(ji,jj,jk) = tn(ji,jj,jk) ! tb <-- tn 225 sb(ji,jj,jk) = sn(ji,jj,jk) 226 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 227 sn(ji,jj,jk) = sa(ji,jj,jk) 228 END DO 229 END DO 230 END DO 231 ELSE ! general case (Leapfrog + Asselin filter 232 DO jk = 1, jpkm1 233 DO jj = 1, jpj 234 DO ji = 1, jpi 235 !RBvvl for reproducibility 236 ! ztf = atfp * ( ta(ji,jj,jk) - 2.* tn(ji,jj,jk) + tb(ji,jj,jk) ) ! Asselin filter on t 237 ! zsf = atfp * ( sa(ji,jj,jk) - 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 238 ! tb(ji,jj,jk) = tn(ji,jj,jk) + ztf ! tb <-- filtered tn 239 ! sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 240 tb(ji,jj,jk) = atfp * ( tb(ji,jj,jk) + ta(ji,jj,jk) ) + atfp1 * tn(ji,jj,jk) 241 sb(ji,jj,jk) = atfp * ( sb(ji,jj,jk) + sa(ji,jj,jk) ) + atfp1 * sn(ji,jj,jk) 242 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 243 sn(ji,jj,jk) = sa(ji,jj,jk) 244 END DO 245 END DO 246 END DO 247 ENDIF 248 ENDIF 249 ! 250 END SUBROUTINE tra_nxt_fix 251 252 253 SUBROUTINE tra_nxt_vvl( kt ) 254 !!---------------------------------------------------------------------- 255 !! *** ROUTINE tra_nxt_vvl *** 256 !! 257 !! ** Purpose : Time varying volume: apply the Asselin time filter 258 !! and swap the tracer fields. 259 !! 260 !! ** Method : - Apply a thickness weighted Asselin time filter on now fields. 261 !! - save in (ta,sa) a thickness weighted average over the three 262 !! time levels which will be used to compute rdn and thus the semi- 263 !! implicit hydrostatic pressure gradient (ln_dynhpg_imp = T) 264 !! - swap tracer fields to prepare the next time_step. 265 !! This can be summurized for tempearture as: 266 !! ztm = (e3t_a*ta+2*e3t_n*tn+e3t_b*tb) ln_dynhpg_imp = T 267 !! /(e3t_a +2*e3t_n +e3t_b ) 268 !! ztm = 0 otherwise 269 !! tb = ( e3t_n*tn + atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) 270 !! /( e3t_n + atfp*[ e3t_b - 2 e3t_n + e3t_a ] ) 271 !! tn = ta 272 !! ta = zt (NB: reset to 0 after eos_bn2 call) 273 !! 274 !! ** Action : - (tb,sb) and (tn,sn) ready for the next time step 275 !! - (ta,sa) time averaged (t,s) (ln_dynhpg_imp = T) 276 !!---------------------------------------------------------------------- 277 INTEGER, INTENT(in) :: kt ! ocean time-step index 278 !! 279 280 ! Not yet ready 281 WRITE(*,*) 'tra_next_vvl : you should not be there' 282 STOP 283 ! 284 END SUBROUTINE tra_nxt_vvl 285 182 286 !!====================================================================== 183 287 END MODULE tranxt
Note: See TracChangeset
for help on using the changeset viewer.