Changeset 1381
- Timestamp:
- 2009-04-06T15:12:22+02:00 (15 years ago)
- Location:
- branches/dev_004_VVL/NEMO/OPA_SRC
- Files:
-
- 4 edited
-
DYN/dynnxt.F90 (modified) (2 diffs)
-
DYN/dynspg_flt.F90 (modified) (1 diff)
-
DYN/wzvmod.F90 (modified) (5 diffs)
-
TRA/tranxt.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynnxt.F90
r1361 r1381 105 105 #else 106 106 IF( lk_vvl ) THEN ! Varying levels 107 !RB_vvl scale factors are wrong at this point108 107 DO jk = 1, jpkm1 109 ua(ji,jj,jk) = ( ub(:,:,jk) * fse3u (:,:,jk) &110 & + z2dt * ua(:,:,jk) * fse3u (:,:,jk) ) &111 & / fse3u (:,:,jk) * umask(:,:,jk)112 va(ji,jj,jk) = ( vb(:,:,jk) * fse3v (:,:,jk) &113 & + z2dt * va(:,:,jk) * fse3v (:,:,jk) ) &114 & / fse3v (:,:,jk) * vmask(:,:,jk)108 ua(ji,jj,jk) = ( ub(:,:,jk) * fse3u_b(:,:,jk) & 109 & + z2dt * ua(:,:,jk) * fse3u_n(:,:,jk) ) & 110 & / fse3u_a(:,:,jk) * umask(:,:,jk) 111 va(ji,jj,jk) = ( vb(:,:,jk) * fse3v_b(:,:,jk) & 112 & + z2dt * va(:,:,jk) * fse3v_n(:,:,jk) ) & 113 & / fse3v_a(:,:,jk) * vmask(:,:,jk) 115 114 END DO 116 115 ELSE … … 200 199 ELSE 201 200 IF( lk_vvl ) THEN ! Varying levels 202 ! Not done 201 DO jk = 1, jpkm1 ! filter applied on thickness weighted velocities 202 DO jj = 1, jpj 203 DO ji = 1, jpi 204 ze3u_a = fse3u_a(ji,jj,jk) 205 ze3v_a = fse3v_a(ji,jj,jk) 206 ze3u_n = fse3u_n(ji,jj,jk) 207 ze3v_n = fse3v_n(ji,jj,jk) 208 ze3u_b = fse3u_b(ji,jj,jk) 209 ze3v_b = fse3v_b(ji,jj,jk) 210 ! 211 zue3a = ua(ji,jj,jk) * ze3u_a 212 zve3a = va(ji,jj,jk) * ze3v_a 213 zue3n = un(ji,jj,jk) * ze3u_n 214 zve3n = vn(ji,jj,jk) * ze3v_n 215 zue3b = ub(ji,jj,jk) * ze3u_b 216 zve3b = vb(ji,jj,jk) * ze3v_b 217 ! 218 ub(ji,jj,jk) = ( atfp * ( zue3b + zue3a ) + atfp1 * zue3n ) & 219 & / ( atfp * ( ze3u_b + ze3u_a ) + atfp1 * ze3u_n ) * umask(ji,jj,jk) 220 vb(ji,jj,jk) = ( atfp * ( zve3b + zve3a ) + atfp1 * zve3n ) & 221 & / ( atfp * ( ze3v_b + ze3v_a ) + atfp1 * ze3v_n ) * vmask(ji,jj,jk) 222 un(ji,jj,jk) = ua(ji,jj,jk) * umask(ji,jj,jk) 223 vn(ji,jj,jk) = va(ji,jj,jk) * vmask(ji,jj,jk) 224 END DO 225 END DO 226 END DO 203 227 ELSE ! Fixed levels 204 228 !RB_vvl : should be done as in tranxt ? -
branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynspg_flt.F90
r1380 r1381 158 158 IF( lk_vvl ) THEN ! variable volume 159 159 160 DO jj = 1, jpjm1161 DO ji = 1,jpim1162 163 ! Sea Surface Height at u-point before164 zsshub(ji,jj) = 0.5 * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) ) &165 & * ( e1t(ji ,jj) * e2t(ji ,jj) * sshb(ji ,jj) &166 & + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) )167 168 ! Sea Surface Height at v-point before169 zsshvb(ji,jj) = 0.5 * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) ) &170 & * ( e1t(ji,jj ) * e2t(ji,jj ) * sshb(ji,jj ) &171 & + e1t(ji,jj+1) * e2t(ji,jj+1) * sshb(ji,jj+1) )172 173 ! Sea Surface Height at u-point after174 zsshua(ji,jj) = 0.5 * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) ) &175 & * ( e1t(ji ,jj) * e2t(ji ,jj) * ssha(ji ,jj) &176 & + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) )177 178 ! Sea Surface Height at v-point after179 zsshva(ji,jj) = 0.5 * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) ) &180 & * ( e1t(ji,jj ) * e2t(ji,jj ) * ssha(ji,jj ) &181 & + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) )182 183 END DO184 END DO185 186 ! Boundaries conditions187 CALL lbc_lnk( zsshub, 'U', 1. ) ; CALL lbc_lnk( zsshvb, 'V', 1. )188 CALL lbc_lnk( zsshua, 'U', 1. ) ; CALL lbc_lnk( zsshva, 'V', 1. )189 190 ! Thickness weighting191 ! -------------------192 DO jk = 1, jpkm1193 DO jj = 1, jpj194 DO ji = 1, jpi195 ua(ji,jj,jk) = ua(ji,jj,jk) * fse3u(ji,jj,jk)196 va(ji,jj,jk) = va(ji,jj,jk) * fse3v(ji,jj,jk)197 198 zub(ji,jj,jk) = ub(ji,jj,jk) * zfse3ub(ji,jj,jk)199 zvb(ji,jj,jk) = vb(ji,jj,jk) * zfse3vb(ji,jj,jk)200 END DO201 END DO202 END DO203 204 160 ! Evaluate the masked next velocity (effect of the additional force not included) 161 ! ------------------- (thickness weighted velocity, surface pressure gradient already included in dyn_hpg) 205 162 DO jk = 1, jpkm1 206 163 DO jj = 2, jpjm1 207 164 DO ji = fs_2, fs_jpim1 ! vector opt. 208 ua(ji,jj,jk) = ( zub(ji,jj,jk) + z2dt * ua(ji,jj,jk) ) /zfse3ua(ji,jj,jk) * umask(ji,jj,jk) 209 va(ji,jj,jk) = ( zvb(ji,jj,jk) + z2dt * va(ji,jj,jk) ) /zfse3va(ji,jj,jk) * vmask(ji,jj,jk) 165 ua(ji,jj,jk) = ( ub(ji,jj,jk) * fse3u_b(ji,jj,jk) & 166 & + z2dt * ua(ji,jj,jk) * fse3u_n(ji,jj,jk) ) & 167 & / fse3u_a(ji,jj,jk) * umask(ji,jj,jk) 168 va(ji,jj,jk) = ( vb(ji,jj,jk) * fse3v_b(ji,jj,jk) & 169 & + z2dt * va(ji,jj,jk) * fse3v_n(ji,jj,jk) ) & 170 & / fse3v_a(ji,jj,jk) * vmask(ji,jj,jk) 210 171 END DO 211 172 END DO -
branches/dev_004_VVL/NEMO/OPA_SRC/DYN/wzvmod.F90
r1368 r1381 114 114 ! ! Sea Surface Height at u-,v- and f-points (vvl case only) 115 115 IF( lk_vvl ) THEN ! (required only in key_vvl case) 116 !RB_vvl 117 ! Compute ssh at u,v, f points and update vertical coordinate 118 ! Currently done in dom_vvl 116 DO jj = 1, jpjm1 117 DO ji = 1, fs_jpim1 ! Vector Opt. 118 sshu_a(ji,jj) = 0.5 * umask(ji,jj,1) / ( e1u(ji ,jj) * e2u(ji ,jj) ) & 119 & * ( e1t(ji ,jj) * e2t(ji ,jj) * ssha(ji ,jj) & 120 & + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) ) 121 sshv_a(ji,jj) = 0.5 * vmask(ji,jj,1) / ( e1v(ji,jj ) * e2v(ji,jj ) ) & 122 & * ( e1t(ji,jj ) * e2t(ji,jj ) * ssha(ji,jj ) & 123 & + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) ) 124 !!gm bug used of fmask, even if thereafter multiplied by muf which is correctly masked) 125 sshf_a(ji,jj) = 0.25 * fmask(ji,jj,1) * ( ssha(ji ,jj) + ssha(ji ,jj+1) & 126 & + ssha(ji+1,jj) + ssha(ji+1,jj+1) ) 127 END DO 128 END DO 129 CALL lbc_lnk( sshu_a, 'U', 1. ) ! Boundaries conditions 130 CALL lbc_lnk( sshv_a, 'V', 1. ) 131 CALL lbc_lnk( sshf_a, 'F', 1. ) 119 132 ENDIF 120 133 … … 134 147 IF( lk_vvl ) THEN ! only in vvl case) 135 148 ! ! now local depth and scale factors (stored in fse3. arrays) 136 !RB_vvl to be done 137 138 149 DO jk = 1, jpkm1 150 fsdept(:,:,jk) = fsdept_n(:,:,jk) ! depths 151 fsdepw(:,:,jk) = fsdepw_n(:,:,jk) 152 fsde3w(:,:,jk) = fsde3w_n(:,:,jk) 153 ! 154 fse3t (:,:,jk) = fse3t_n (:,:,jk) ! vertical scale factors 155 fse3u (:,:,jk) = fse3u_n (:,:,jk) 156 fse3v (:,:,jk) = fse3v_n (:,:,jk) 157 fse3f (:,:,jk) = fse3f_n (:,:,jk) 158 fse3w (:,:,jk) = fse3w_n (:,:,jk) 159 fse3uw(:,:,jk) = fse3uw_n(:,:,jk) 160 fse3vw(:,:,jk) = fse3vw_n(:,:,jk) 161 END DO 162 ! ! ocean depth (at u- and v-points) 163 hu(:,:) = hu_0(:,:) + sshu_n(:,:) 164 hv(:,:) = hv_0(:,:) + sshv_n(:,:) 165 ! ! masked inverse of the ocean depth (at u- and v-points) 139 166 ENDIF 140 167 ! … … 171 198 ! ------------------------------- 172 199 173 !RB_vvl ssh at u, f,v point to be added174 200 175 201 IF( neuler == 0 .AND. kt == nit000 ) THEN ! Euler time stepping … … 178 204 ! before <-- now 179 205 sshb (ji,jj) = sshn(ji,jj) 206 sshu_b(ji,jj) = sshu_n(ji,jj) 207 sshv_b(ji,jj) = sshv_n(ji,jj) 208 sshf_b(ji,jj) = sshf_n(ji,jj) 180 209 ! now <-- after 181 210 sshn (ji,jj) = ssha (ji,jj) 211 sshu_n(ji,jj) = sshu_a(ji,jj) 212 sshv_n(ji,jj) = sshv_a(ji,jj) 213 sshf_n(ji,jj) = sshf_a(ji,jj) 182 214 END DO 183 215 END DO … … 187 219 ! before <-- now filtered 188 220 sshb (ji,jj) = sshn(ji,jj) + atfp * ( sshb (ji,jj) - 2 * sshn (ji,jj) + ssha (ji,jj) ) !& 221 sshu_b(ji,jj) = sshu_n(ji,jj) + atfp * ( sshu_b(ji,jj) - 2 * sshu_n(ji,jj) + sshu_a(ji,jj) ) !& 222 sshv_b(ji,jj) = sshv_n(ji,jj) + atfp * ( sshv_b(ji,jj) - 2 * sshv_n(ji,jj) + sshv_a(ji,jj) ) !& 223 sshf_b(ji,jj) = sshf_n(ji,jj) + atfp * ( sshf_b(ji,jj) - 2 * sshf_n(ji,jj) + sshf_a(ji,jj) ) !& 189 224 ! now <-- after 190 225 sshn (ji,jj) = ssha (ji,jj) 226 sshu_n(ji,jj) = sshu_a(ji,jj) 227 sshv_n(ji,jj) = sshv_a(ji,jj) 228 sshf_n(ji,jj) = sshf_a(ji,jj) 191 229 END DO 192 230 END DO -
branches/dev_004_VVL/NEMO/OPA_SRC/TRA/tranxt.F90
r1361 r1381 277 277 INTEGER, INTENT(in) :: kt ! ocean time-step index 278 278 !! 279 280 ! Not yet ready 281 WRITE(*,*) 'tra_next_vvl : you should not be there' 282 STOP 279 INTEGER :: ji, jj, jk ! dummy loop indices 280 REAL(wp) :: ztm , ztc_f , ztf , ztca, ztcn, ztcb ! temporary scalar 281 REAL(wp) :: zsm , zsc_f , zsf , zsca, zscn, zscb ! - - 282 REAL(wp) :: ze3mr, ze3fr ! - - 283 REAL(wp) :: ze3t_b, ze3t_n, ze3t_a, ze3t_f ! - - 284 !!---------------------------------------------------------------------- 285 286 IF( kt == nit000 ) THEN 287 IF(lwp) WRITE(numout,*) 288 IF(lwp) WRITE(numout,*) 'tra_nxt_vvl : time stepping' 289 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 290 ENDIF 291 292 ! ! ----------------------- ! 293 IF( ln_dynhpg_imp ) THEN ! semi-implicite hpg case ! 294 ! ! ----------------------- ! 295 ! 296 IF( neuler == 0 .AND. kt == nit000 ) THEN 297 DO jk = 1, jpkm1 298 DO jj = 1, jpj 299 DO ji = 1, jpi 300 ze3t_b = fse3t_b(ji,jj,jk) 301 ze3t_n = fse3t_n(ji,jj,jk) 302 ze3t_a = fse3t_a(ji,jj,jk) 303 ! ! tracer content at Before, now and after 304 ztcb = tb(ji,jj,jk) * ze3t_b ; zscb = sb(ji,jj,jk) * ze3t_b 305 ztcn = tn(ji,jj,jk) * ze3t_n ; zscn = sn(ji,jj,jk) * ze3t_n 306 ztca = ta(ji,jj,jk) * ze3t_a ; zsca = sa(ji,jj,jk) * ze3t_a 307 ! 308 ! ! mean thickness and tracer 309 ze3mr= 1.e0 / ( ze3t_a + 2.* ze3t_n + ze3t_b ) 310 ztm = ze3mr * ( ztca + 2.* ztcn + ztcb ) 311 zsm = ze3mr * ( zsca + 2.* zscn + zscb ) 312 !!gm mean e3t have to be saved and used in dynhpg or it can be recomputed in dynhpg !! 313 !!gm e3t_m(ji,jj,jk) = 0.25 / ze3mr 314 ! ! swap of arrays 315 tb(ji,jj,jk) = tn(ji,jj,jk) ! tb <-- tn 316 sb(ji,jj,jk) = sn(ji,jj,jk) 317 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 318 sn(ji,jj,jk) = sa(ji,jj,jk) 319 ta(ji,jj,jk) = ztm ! ta <-- mean t 320 sa(ji,jj,jk) = zsm 321 END DO 322 END DO 323 END DO 324 ELSE 325 DO jk = 1, jpkm1 326 DO jj = 1, jpj 327 DO ji = 1, jpi 328 ze3t_b = fse3t_b(ji,jj,jk) 329 ze3t_n = fse3t_n(ji,jj,jk) 330 ze3t_a = fse3t_a(ji,jj,jk) 331 ! ! tracer content at Before, now and after 332 ztcb = tb(ji,jj,jk) * ze3t_b ; zscb = sb(ji,jj,jk) * ze3t_b 333 ztcn = tn(ji,jj,jk) * ze3t_n ; zscn = sn(ji,jj,jk) * ze3t_n 334 ztca = ta(ji,jj,jk) * ze3t_a ; zsca = sa(ji,jj,jk) * ze3t_a 335 ! 336 ! ! Asselin filter on thickness and tracer content 337 ze3t_f = atfp * ( ze3t_a - 2.* ze3t_n + ze3t_b ) 338 ztc_f = atfp * ( ztca - 2.* ztcn + ztcb ) 339 zsc_f = atfp * ( zsca - 2.* zscn + zscb ) 340 ! 341 ! ! filtered tracer including the correction 342 ze3fr = 1.e0 / ( ze3t_n + ze3t_f ) 343 ztf = ze3fr * ( ztcn + ztc_f ) 344 zsf = ze3fr * ( zscn + zsc_f ) 345 ! ! mean thickness and tracer 346 ze3mr = 1.e0 / ( ze3t_a + 2.* ze3t_n + ze3t_b ) 347 ztm = ze3mr * ( ztca + 2.* ztcn + ztcb ) 348 zsm = ze3mr * ( zsca + 2.* zscn + zscb ) 349 !!gm mean e3t have to be saved and used in dynhpg or it can be recomputed in dynhpg !! 350 !!gm e3t_m(ji,jj,jk) = 0.25 / ze3mr 351 ! ! swap of arrays 352 tb(ji,jj,jk) = ztf ! tb <-- tn + filter 353 sb(ji,jj,jk) = zsf 354 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 355 sn(ji,jj,jk) = sa(ji,jj,jk) 356 ta(ji,jj,jk) = ztm ! ta <-- mean t 357 sa(ji,jj,jk) = zsm 358 END DO 359 END DO 360 END DO 361 ENDIF 362 ! ! ----------------------- ! 363 ELSE ! explicit hpg case ! 364 ! ! ----------------------- ! 365 ! 366 IF( neuler == 0 .AND. kt == nit000 ) THEN ! case of Euler time-stepping at first time-step 367 DO jk = 1, jpkm1 ! No filter nor thickness weighting computation required 368 DO jj = 1, jpj ! ONLY swap 369 DO ji = 1, jpi 370 tb(ji,jj,jk) = tn(ji,jj,jk) ! tb <-- tn 371 sb(ji,jj,jk) = sn(ji,jj,jk) 372 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 373 sn(ji,jj,jk) = sa(ji,jj,jk) 374 END DO 375 END DO 376 END DO 377 ! ! general case (Leapfrog + Asselin filter) 378 ELSE ! apply filter on thickness weighted tracer and swap 379 DO jk = 1, jpkm1 380 DO jj = 1, jpj 381 DO ji = 1, jpi 382 ze3t_b = fse3t_b(ji,jj,jk) 383 ze3t_n = fse3t_n(ji,jj,jk) 384 ze3t_a = fse3t_a(ji,jj,jk) 385 ! ! tracer content at Before, now and after 386 ztcb = tb(ji,jj,jk) * ze3t_b ; zscb = sb(ji,jj,jk) * ze3t_b 387 ztcn = tn(ji,jj,jk) * ze3t_n ; zscn = sn(ji,jj,jk) * ze3t_n 388 ztca = ta(ji,jj,jk) * ze3t_a ; zsca = sa(ji,jj,jk) * ze3t_a 389 ! 390 ! ! Asselin filter on thickness and tracer content 391 ze3t_f = atfp * ( ze3t_a - 2.* ze3t_n + ze3t_b ) 392 ztc_f = atfp * ( ztca - 2.* ztcn + ztcb ) 393 zsc_f = atfp * ( zsca - 2.* zscn + zscb ) 394 ! 395 !!gm tmask useless below 396 ! ! filtered tracer including the correction 397 ze3fr = tmask(ji,jj,jk) / ( ze3t_n + ze3t_f ) 398 ztf = ( ztcn + ztc_f ) * ze3fr 399 zsf = ( zscn + zsc_f ) * ze3fr 400 ! ! swap of arrays 401 tb(ji,jj,jk) = ztf ! tb <-- tn filtered 402 sb(ji,jj,jk) = zsf 403 tn(ji,jj,jk) = ta(ji,jj,jk) ! tn <-- ta 404 sn(ji,jj,jk) = sa(ji,jj,jk) 405 END DO 406 END DO 407 END DO 408 ENDIF 409 ENDIF 283 410 ! 284 411 END SUBROUTINE tra_nxt_vvl
Note: See TracChangeset
for help on using the changeset viewer.
