- Timestamp:
- 2017-06-06T15:55:44+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r8093 r8143 16 16 !! 3.7 ! 2015-11 (J. Chanut) free surface simplification 17 17 !! - ! 2016-12 (G. Madec, E. Clementi) update for Stoke-Drift divergence 18 !! 4.0 ! 2017-05 (G. Madec) drag coef. defined at t-point (zdfdrg.F90) 18 19 !!--------------------------------------------------------------------- 19 20 … … 27 28 USE dom_oce ! ocean space and time domain 28 29 USE sbc_oce ! surface boundary condition: ocean 29 USE zdf_oce ! Bottom friction coefts 30 !!gm new 31 USE zdfdrg ! vertical physics: top/bottom drag coef. 32 !!gm 30 USE zdf_oce ! vertical physics: variables 31 USE zdfdrg ! vertical physics: top/bottom drag coef. 33 32 USE sbcisf ! ice shelf variable (fwfisf) 34 33 USE sbcapr ! surface boundary condition: atmospheric pressure … … 43 42 USE updtide ! tide potential 44 43 USE sbcwave ! surface wave 44 USE diatmb ! Top,middle,bottom output 45 #if defined key_agrif 46 USE agrif_opa_interp ! agrif 47 #endif 48 #if defined key_asminc 49 USE asminc ! Assimilation increment 50 #endif 45 51 ! 46 52 USE in_out_manager ! I/O manager … … 50 56 USE iom ! IOM library 51 57 USE restart ! only for lrst_oce 52 USE wrk_nemo ! Memory Allocation53 58 USE timing ! Timing 54 USE diatmb ! Top,middle,bottom output55 #if defined key_agrif56 USE agrif_opa_interp ! agrif57 #endif58 #if defined key_asminc59 USE asminc ! Assimilation increment60 #endif61 62 59 63 60 IMPLICIT NONE … … 69 66 PUBLIC ts_rst ! " " " " 70 67 71 INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_baro <= 2.5 nn_baro72 REAL(wp),SAVE :: rdtbt ! Barotropic time step73 74 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: wgtbtp1, wgtbtp2 !: 1st & 2nd weights used in time filtering of barotropic fields75 76 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zwz !: ff_f/h at F points77 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ftnw, ftne !: triad of coriolis parameter78 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ftsw, ftse !: (only used with een vorticity scheme)79 80 68 !! Time filtered arrays at baroclinic time step: 81 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv , vn_adv !: Advection vel. at "now" barocl. step 69 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv , vn_adv !: Advection vel. at "now" barocl. step 70 71 INTEGER , SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_baro <= 2.5 nn_baro 72 REAL(wp), SAVE :: rdtbt ! Barotropic time step 73 ! 74 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: wgtbtp1, wgtbtp2 ! 1st & 2nd weights used in time filtering of barotropic fields 75 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zwz ! ff_f/h at F points 76 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ftnw, ftne ! triad of coriolis parameter 77 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ftsw, ftse ! (only used with een vorticity scheme) 78 79 REAL(wp) :: r1_12 = 1._wp / 12._wp ! local ratios 80 REAL(wp) :: r1_8 = 0.125_wp ! 81 REAL(wp) :: r1_4 = 0.25_wp ! 82 REAL(wp) :: r1_2 = 0.5_wp ! 82 83 83 84 !! * Substitutions 84 85 # include "vectopt_loop_substitute.h90" 85 86 !!---------------------------------------------------------------------- 86 !! NEMO/OPA 3.5 , NEMO Consortium (2013)87 !! NEMO/OPA 4.0 , NEMO Consortium (2017) 87 88 !! $Id$ 88 89 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 140 141 INTEGER, INTENT(in) :: kt ! ocean time-step index 141 142 ! 142 LOGICAL :: ll_fw_start ! if true, forward integration143 LOGICAL :: ll_init ! if true, special startup of 2d equations144 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables used in W/D145 143 INTEGER :: ji, jj, jk, jn ! dummy loop indices 146 INTEGER :: ikbu, ikbv, noffset ! local integers 147 INTEGER :: iktu, iktv ! local integers 148 REAL(wp) :: zmdi 149 REAL(wp) :: zraur, z1_2dt_b, z2dt_bf ! local scalars 150 REAL(wp) :: zx1, zy1, zx2, zy2 ! - - 151 REAL(wp) :: z1_12, z1_8, z1_4, z1_2 ! - - 152 REAL(wp) :: zu_spg, zv_spg ! - - 153 REAL(wp) :: zhura, zhvra ! - - 154 REAL(wp) :: za0, za1, za2, za3 ! - - 155 REAL(wp) :: zztmp ! - - 156 REAL(wp), DIMENSION(jpi,jpj) :: zsshp2_e 157 REAL(wp), DIMENSION(jpi,jpj) :: zu_trd, zv_trd, zu_frc, zv_frc, zssh_frc 158 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zwy, zhdiv 159 REAL(wp), DIMENSION(jpi,jpj) :: zhup2_e, zhvp2_e, zhust_e, zhvst_e 160 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_a, zsshv_a 161 REAL(wp), DIMENSION(jpi,jpj) :: zhf 144 LOGICAL :: ll_fw_start ! =T : forward integration 145 LOGICAL :: ll_init ! =T : special startup of 2d equations 146 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables used in W/D 147 INTEGER :: ikbu, iktu, noffset ! local integers 148 INTEGER :: ikbv, iktv ! - - 149 REAL(wp) :: z1_2dt_b, z2dt_bf ! local scalars 150 REAL(wp) :: zx1, zx2, zu_spg, zhura ! - - 151 REAL(wp) :: zy1, zy2, zv_spg, zhvra ! - - 152 REAL(wp) :: za0, za1, za2, za3 ! - - 153 REAL(wp) :: zmdi, zztmp ! - - 154 REAL(wp), DIMENSION(jpi,jpj) :: zsshp2_e, zhf 155 REAL(wp), DIMENSION(jpi,jpj) :: zwx, zu_trd, zu_frc, zssh_frc 156 REAL(wp), DIMENSION(jpi,jpj) :: zwy, zv_trd, zv_frc, zhdiv 157 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_a, zhup2_e, zhust_e 158 REAL(wp), DIMENSION(jpi,jpj) :: zsshv_a, zhvp2_e, zhvst_e 162 159 REAL(wp), DIMENSION(jpi,jpj) :: zCdU_u, zCdU_v ! top/bottom stress at u- & v-points 163 160 ! … … 170 167 ! 171 168 zmdi=1.e+20 ! missing data indicator for masking 172 ! !* Local constant initialization 173 z1_12 = 1._wp / 12._wp 174 z1_8 = 0.125_wp 175 z1_4 = 0.25_wp 176 z1_2 = 0.5_wp 177 zraur = 1._wp / rau0 169 ! 178 170 ! ! reciprocal of baroclinic time step 179 171 IF( kt == nit000 .AND. neuler == 0 ) THEN ; z2dt_bf = rdt … … 210 202 ENDIF 211 203 ! 212 !!gm old/new213 204 IF( ln_isfcav ) THEN ! top+bottom friction (ocean cavities) 214 zCdU_u(:,:) = bfrua(:,:) + tfrua(:,:) 215 zCdU_v(:,:) = bfrva(:,:) + tfrva(:,:) 205 DO jj = 2, jpjm1 206 DO ji = fs_2, fs_jpim1 ! vector opt. 207 zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 208 zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 209 END DO 210 END DO 216 211 ELSE ! bottom friction only 217 zCdU_u(:,:) = bfrua(:,:) 218 zCdU_v(:,:) = bfrva(:,:) 219 ENDIF 220 !!gm new 221 ! IF( ln_isfcav ) THEN ! top+bottom friction (ocean cavities) 222 ! DO jj = 2, jpjm1 223 ! DO ji = fs_2, fs_jpim1 ! vector opt. 224 ! zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 225 ! zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 226 ! END DO 227 ! END DO 228 ! ELSE ! bottom friction only 229 ! DO jj = 2, jpjm1 230 ! DO ji = fs_2, fs_jpim1 ! vector opt. 231 ! zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 232 ! zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 233 ! END DO 234 ! END DO 235 ! ENDIF 236 !!gm 212 DO jj = 2, jpjm1 213 DO ji = fs_2, fs_jpim1 ! vector opt. 214 zCdU_u(ji,jj) = 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 215 zCdU_v(ji,jj) = 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 216 END DO 217 END DO 218 ENDIF 237 219 ! 238 220 ! Set arrays to remove/compute coriolis trend. … … 287 269 !!gm 288 270 !! 289 IF ( .not.ln_sco ) THEN271 IF( .NOT.ln_sco ) THEN 290 272 291 273 !!gm agree the JC comment : this should be done in a much clear way … … 338 320 IF (.NOT.ln_bt_fw .AND.( neuler==0 .AND. kt==nit000+1 ) ) THEN 339 321 ll_fw_start=.FALSE. 340 CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2)322 CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 341 323 ENDIF 342 324 … … 387 369 zx2 = ( zwx(ji ,jj) + zwx(ji ,jj+1) ) * r1_e2v(ji,jj) 388 370 ! energy conserving formulation for planetary vorticity term 389 zu_trd(ji,jj) = z1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 )390 zv_trd(ji,jj) = -z1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 )371 zu_trd(ji,jj) = r1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 372 zv_trd(ji,jj) = - r1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 ) 391 373 END DO 392 374 END DO … … 395 377 DO jj = 2, jpjm1 396 378 DO ji = fs_2, fs_jpim1 ! vector opt. 397 zy1 = z1_8 * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) &379 zy1 = r1_8 * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) & 398 380 & + zwy(ji ,jj ) + zwy(ji+1,jj ) ) * r1_e1u(ji,jj) 399 zx1 = - z1_8 * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) &381 zx1 = - r1_8 * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) & 400 382 & + zwx(ji ,jj ) + zwx(ji ,jj+1) ) * r1_e2v(ji,jj) 401 383 zu_trd(ji,jj) = zy1 * ( zwz(ji ,jj-1) + zwz(ji,jj) ) … … 407 389 DO jj = 2, jpjm1 408 390 DO ji = fs_2, fs_jpim1 ! vector opt. 409 zu_trd(ji,jj) = + z1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zwy(ji ,jj ) &391 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zwy(ji ,jj ) & 410 392 & + ftnw(ji+1,jj) * zwy(ji+1,jj ) & 411 393 & + ftse(ji,jj ) * zwy(ji ,jj-1) & 412 394 & + ftsw(ji+1,jj) * zwy(ji+1,jj-1) ) 413 zv_trd(ji,jj) = - z1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zwx(ji-1,jj+1) &395 zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zwx(ji-1,jj+1) & 414 396 & + ftse(ji,jj+1) * zwx(ji ,jj+1) & 415 397 & + ftnw(ji,jj ) * zwx(ji-1,jj ) & … … 423 405 ! ! ---------------------------------------------------- 424 406 IF( .NOT.ln_linssh ) THEN ! Variable volume : remove surface pressure gradient 425 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters426 DO jj = 2, jpjm1427 DO ji = 2, jpim1428 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > &429 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. &430 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) &407 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters 408 DO jj = 2, jpjm1 409 DO ji = 2, jpim1 410 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 411 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 412 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 431 413 & > rn_wdmin1 + rn_wdmin2 432 ll_tmp2 = ( ABS( sshn(ji+1,jj) - sshn(ji ,jj)) > 1.E-12 ).AND.( & 433 & MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 434 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 435 436 IF(ll_tmp1) THEN 437 zcpx(ji,jj) = 1.0_wp 438 ELSE IF(ll_tmp2) THEN 439 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 440 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 441 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 442 ELSE 443 zcpx(ji,jj) = 0._wp 444 END IF 445 446 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 414 ll_tmp2 = ( ABS( sshn(ji+1,jj) - sshn(ji ,jj)) > 1.E-12 ).AND.( & 415 & MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 416 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 417 ! 418 IF(ll_tmp1) THEN 419 zcpx(ji,jj) = 1.0_wp 420 ELSE IF(ll_tmp2) THEN ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 421 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 422 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 423 ELSE 424 zcpx(ji,jj) = 0._wp 425 ENDIF 426 ! 427 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 447 428 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 448 429 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 449 430 & > rn_wdmin1 + rn_wdmin2 450 ll_tmp2 = ( ABS( sshn(ji,jj) - sshn(ji,jj+1)) > 1.E-12 ).AND.( &431 ll_tmp2 = ( ABS( sshn(ji,jj) - sshn(ji,jj+1)) > 1.E-12 ).AND.( & 451 432 & MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 452 433 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 ) 453 454 IF(ll_tmp1) THEN455 zcpy(ji,jj) = 1.0_wp456 ELSE IF(ll_tmp2) THEN457 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here458 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) &459 &/ (sshn(ji,jj+1) - sshn(ji,jj )) )460 ELSE461 zcpy(ji,jj) = 0._wp462 ENDIF463 END DO464 END DO465 466 DO jj = 2, jpjm1467 DO ji = 2, jpim1468 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) &469 &* r1_e1u(ji,jj) * zcpx(ji,jj)470 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) &471 &* r1_e2v(ji,jj) * zcpy(ji,jj)472 END DO473 END DO474 434 ! 435 IF(ll_tmp1) THEN 436 zcpy(ji,jj) = 1.0_wp 437 ELSE IF(ll_tmp2) THEN 438 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 439 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 440 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 441 ELSE 442 zcpy(ji,jj) = 0._wp 443 ENDIF 444 END DO 445 END DO 446 ! 447 DO jj = 2, jpjm1 448 DO ji = 2, jpim1 449 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 450 & * r1_e1u(ji,jj) * zcpx(ji,jj) 451 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) & 452 & * r1_e2v(ji,jj) * zcpy(ji,jj) 453 END DO 454 END DO 455 ! 475 456 ELSE 476 477 DO jj = 2, jpjm1478 DO ji = fs_2, fs_jpim1 ! vector opt.479 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) * r1_e1u(ji,jj)480 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) * r1_e2v(ji,jj)481 END DO482 END DO483 ENDIF484 485 ENDIF 486 457 ! 458 DO jj = 2, jpjm1 459 DO ji = fs_2, fs_jpim1 ! vector opt. 460 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) * r1_e1u(ji,jj) 461 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) * r1_e2v(ji,jj) 462 END DO 463 END DO 464 ENDIF 465 ! 466 ENDIF 467 ! 487 468 DO jj = 2, jpjm1 ! Remove coriolis term (and possibly spg) from barotropic trend 488 469 DO ji = fs_2, fs_jpim1 … … 492 473 END DO 493 474 ! 494 ! ! Add bottomstress contribution from baroclinic velocities:495 IF (ln_bt_fw) THEN475 ! ! Add BOTTOM stress contribution from baroclinic velocities: 476 IF( ln_bt_fw ) THEN 496 477 DO jj = 2, jpjm1 497 478 DO ji = fs_2, fs_jpim1 ! vector opt. … … 518 499 DO jj = 2, jpjm1 519 500 DO ji = fs_2, fs_jpim1 ! vector opt. 520 !!gm old 521 zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * bfrua(ji,jj) , zztmp ) * zwx(ji,jj) 522 zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * bfrva(ji,jj) , zztmp ) * zwy(ji,jj) 523 !!gm new 524 ! zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) * zwx(ji,jj) 525 ! zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) * zwy(ji,jj) 526 !!gm 501 zu_frc(ji,jj) = zu_frc(ji,jj) + MAX( r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) * zwx(ji,jj) 502 zv_frc(ji,jj) = zv_frc(ji,jj) + MAX( r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) * zwy(ji,jj) 527 503 END DO 528 504 END DO … … 530 506 DO jj = 2, jpjm1 531 507 DO ji = fs_2, fs_jpim1 ! vector opt. 532 !!gm old 533 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * bfrua(ji,jj) * zwx(ji,jj) 534 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * bfrva(ji,jj) * zwy(ji,jj) 535 !!gm new 536 ! zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 537 ! zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 538 !!gm 508 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 509 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 539 510 END DO 540 511 END DO 541 512 END IF 542 513 ! 543 ! ! Add top stress contribution from baroclinic velocities: 544 IF( ln_bt_fw ) THEN 514 IF( ln_isfcav ) THEN ! Add TOP stress contribution from baroclinic velocities: 515 IF( ln_bt_fw ) THEN 516 DO jj = 2, jpjm1 517 DO ji = fs_2, fs_jpim1 ! vector opt. 518 iktu = miku(ji,jj) 519 iktv = mikv(ji,jj) 520 zwx(ji,jj) = un(ji,jj,iktu) - un_b(ji,jj) ! NOW top baroclinic velocities 521 zwy(ji,jj) = vn(ji,jj,iktv) - vn_b(ji,jj) 522 END DO 523 END DO 524 ELSE 525 DO jj = 2, jpjm1 526 DO ji = fs_2, fs_jpim1 ! vector opt. 527 iktu = miku(ji,jj) 528 iktv = mikv(ji,jj) 529 zwx(ji,jj) = ub(ji,jj,iktu) - ub_b(ji,jj) ! BEFORE top baroclinic velocities 530 zwy(ji,jj) = vb(ji,jj,iktv) - vb_b(ji,jj) 531 END DO 532 END DO 533 ENDIF 534 ! 535 ! Note that the "unclipped" top friction parameter is used even with explicit drag 536 DO jj = 2, jpjm1 537 DO ji = fs_2, fs_jpim1 ! vector opt. 538 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zwx(ji,jj) 539 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * zwy(ji,jj) 540 END DO 541 END DO 542 ENDIF 543 ! 544 IF( ln_bt_fw ) THEN ! Add wind forcing 545 545 DO jj = 2, jpjm1 546 546 DO ji = fs_2, fs_jpim1 ! vector opt. 547 iktu = miku(ji,jj) 548 iktv = mikv(ji,jj) 549 zwx(ji,jj) = un(ji,jj,iktu) - un_b(ji,jj) ! NOW top baroclinic velocities 550 zwy(ji,jj) = vn(ji,jj,iktv) - vn_b(ji,jj) 547 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 548 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) 551 549 END DO 552 550 END DO 553 551 ELSE 552 zztmp = r1_rau0 * r1_2 554 553 DO jj = 2, jpjm1 555 554 DO ji = fs_2, fs_jpim1 ! vector opt. 556 iktu = miku(ji,jj) 557 iktv = mikv(ji,jj) 558 zwx(ji,jj) = ub(ji,jj,iktu) - ub_b(ji,jj) ! BEFORE top baroclinic velocities 559 zwy(ji,jj) = vb(ji,jj,iktv) - vb_b(ji,jj) 560 END DO 561 END DO 562 ENDIF 563 ! 564 ! Note that the "unclipped" top friction parameter is used even with explicit drag 565 DO jj = 2, jpjm1 566 DO ji = fs_2, fs_jpim1 ! vector opt. 567 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * tfrua(ji,jj) * zwx(ji,jj) 568 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * tfrva(ji,jj) * zwy(ji,jj) 569 END DO 570 END DO 571 ! 572 IF (ln_bt_fw) THEN ! Add wind forcing 573 DO jj = 2, jpjm1 574 DO ji = fs_2, fs_jpim1 ! vector opt. 575 zu_frc(ji,jj) = zu_frc(ji,jj) + zraur * utau(ji,jj) * r1_hu_n(ji,jj) 576 zv_frc(ji,jj) = zv_frc(ji,jj) + zraur * vtau(ji,jj) * r1_hv_n(ji,jj) 577 END DO 578 END DO 579 ELSE 580 DO jj = 2, jpjm1 581 DO ji = fs_2, fs_jpim1 ! vector opt. 582 zu_frc(ji,jj) = zu_frc(ji,jj) + zraur * z1_2 * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 583 zv_frc(ji,jj) = zv_frc(ji,jj) + zraur * z1_2 * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 555 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 556 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 584 557 END DO 585 558 END DO 586 559 ENDIF 587 560 ! 588 IF ( ln_apr_dyn ) THEN! Add atm pressure forcing589 IF (ln_bt_fw) THEN561 IF( ln_apr_dyn ) THEN ! Add atm pressure forcing 562 IF( ln_bt_fw ) THEN 590 563 DO jj = 2, jpjm1 591 564 DO ji = fs_2, fs_jpim1 ! vector opt. … … 597 570 END DO 598 571 ELSE 572 zztmp = grav * r1_2 599 573 DO jj = 2, jpjm1 600 574 DO ji = fs_2, fs_jpim1 ! vector opt. 601 zu_spg = grav * z1_2* ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) &602 & 603 zv_spg = grav * z1_2* ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) &604 & 575 zu_spg = zztmp * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) & 576 & + ssh_ibb(ji+1,jj ) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) 577 zv_spg = zztmp * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) & 578 & + ssh_ibb(ji ,jj+1) - ssh_ibb(ji,jj) ) * r1_e2v(ji,jj) 605 579 zu_frc(ji,jj) = zu_frc(ji,jj) + zu_spg 606 580 zv_frc(ji,jj) = zv_frc(ji,jj) + zv_spg … … 613 587 ! ! Surface net water flux and rivers 614 588 IF (ln_bt_fw) THEN 615 zssh_frc(:,:) = zraur* ( emp(:,:) - rnf(:,:) + fwfisf(:,:) )589 zssh_frc(:,:) = r1_rau0 * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) 616 590 ELSE 617 zssh_frc(:,:) = zraur * z1_2 * ( emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:) & 618 & + fwfisf(:,:) + fwfisf_b(:,:) ) 591 zztmp = r1_rau0 * r1_2 592 zssh_frc(:,:) = zztmp * ( emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:) & 593 & + fwfisf(:,:) + fwfisf_b(:,:) ) 619 594 ENDIF 620 595 ! … … 712 687 DO jj = 2, jpjm1 ! Sea Surface Height at u- & v-points 713 688 DO ji = 2, fs_jpim1 ! Vector opt. 714 zwx(ji,jj) = z1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) &689 zwx(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 715 690 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 716 691 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) 717 zwy(ji,jj) = z1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) &692 zwy(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 718 693 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 719 694 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) … … 789 764 DO jj = 2, jpjm1 790 765 DO ji = 2, jpim1 ! NO Vector Opt. 791 zsshu_a(ji,jj) = z1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) &766 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 792 767 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 793 768 & + e1e2t(ji+1,jj ) * ssha_e(ji+1,jj ) ) 794 zsshv_a(ji,jj) = z1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) &769 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 795 770 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 796 771 & + e1e2t(ji ,jj+1) * ssha_e(ji ,jj+1) ) … … 868 843 DO jj = 2, jpjm1 869 844 DO ji = 2, jpim1 870 zx1 = z1_2 * ssumask(ji ,jj) * r1_e1e2u(ji ,jj) &845 zx1 = r1_2 * ssumask(ji ,jj) * r1_e1e2u(ji ,jj) & 871 846 & * ( e1e2t(ji ,jj ) * zsshp2_e(ji ,jj) & 872 847 & + e1e2t(ji+1,jj ) * zsshp2_e(ji+1,jj ) ) 873 zy1 = z1_2 * ssvmask(ji ,jj) * r1_e1e2v(ji ,jj ) &848 zy1 = r1_2 * ssvmask(ji ,jj) * r1_e1e2v(ji ,jj ) & 874 849 & * ( e1e2t(ji ,jj ) * zsshp2_e(ji ,jj ) & 875 850 & + e1e2t(ji ,jj+1) * zsshp2_e(ji ,jj+1) ) … … 895 870 zx1 = ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) ) * r1_e2v(ji,jj) 896 871 zx2 = ( zwx(ji ,jj ) + zwx(ji ,jj+1) ) * r1_e2v(ji,jj) 897 zu_trd(ji,jj) = z1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 )898 zv_trd(ji,jj) =- z1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 )872 zu_trd(ji,jj) = r1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 873 zv_trd(ji,jj) =-r1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 ) 899 874 END DO 900 875 END DO … … 903 878 DO jj = 2, jpjm1 904 879 DO ji = fs_2, fs_jpim1 ! vector opt. 905 zy1 = z1_8 * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) &880 zy1 = r1_8 * ( zwy(ji ,jj-1) + zwy(ji+1,jj-1) & 906 881 & + zwy(ji ,jj ) + zwy(ji+1,jj ) ) * r1_e1u(ji,jj) 907 zx1 = - z1_8 * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) &882 zx1 = - r1_8 * ( zwx(ji-1,jj ) + zwx(ji-1,jj+1) & 908 883 & + zwx(ji ,jj ) + zwx(ji ,jj+1) ) * r1_e2v(ji,jj) 909 884 zu_trd(ji,jj) = zy1 * ( zwz(ji ,jj-1) + zwz(ji,jj) ) … … 915 890 DO jj = 2, jpjm1 916 891 DO ji = fs_2, fs_jpim1 ! vector opt. 917 zu_trd(ji,jj) = + z1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zwy(ji ,jj ) &892 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zwy(ji ,jj ) & 918 893 & + ftnw(ji+1,jj) * zwy(ji+1,jj ) & 919 894 & + ftse(ji,jj ) * zwy(ji ,jj-1) & 920 895 & + ftsw(ji+1,jj) * zwy(ji+1,jj-1) ) 921 zv_trd(ji,jj) = - z1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zwx(ji-1,jj+1) &896 zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zwx(ji-1,jj+1) & 922 897 & + ftse(ji,jj+1) * zwx(ji ,jj+1) & 923 898 & + ftnw(ji,jj ) * zwx(ji-1,jj ) & … … 1082 1057 vn_adv(:,:) = zwy(:,:) * r1_hv_n(:,:) 1083 1058 ELSE 1084 un_adv(:,:) = z1_2 * ( ub2_b(:,:) + zwx(:,:) ) * r1_hu_n(:,:)1085 vn_adv(:,:) = z1_2 * ( vb2_b(:,:) + zwy(:,:) ) * r1_hv_n(:,:)1059 un_adv(:,:) = r1_2 * ( ub2_b(:,:) + zwx(:,:) ) * r1_hu_n(:,:) 1060 vn_adv(:,:) = r1_2 * ( vb2_b(:,:) + zwy(:,:) ) * r1_hv_n(:,:) 1086 1061 END IF 1087 1062 … … 1101 1076 DO jj = 1, jpjm1 1102 1077 DO ji = 1, jpim1 ! NO Vector Opt. 1103 zsshu_a(ji,jj) = z1_2 * umask(ji,jj,1) * r1_e1e2u(ji,jj) &1078 zsshu_a(ji,jj) = r1_2 * umask(ji,jj,1) * r1_e1e2u(ji,jj) & 1104 1079 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 1105 1080 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 1106 zsshv_a(ji,jj) = z1_2 * vmask(ji,jj,1) * r1_e1e2v(ji,jj) &1081 zsshv_a(ji,jj) = r1_2 * vmask(ji,jj,1) * r1_e1e2v(ji,jj) & 1107 1082 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 1108 1083 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) … … 1299 1274 INTEGER :: ji ,jj ! dummy loop indices 1300 1275 REAL(wp) :: zxr2, zyr2, zcmax ! local scalar 1301 REAL(wp), POINTER, DIMENSION(:,:) :: zcu1276 REAL(wp), DIMENSION(jpi,jpj) :: zcu 1302 1277 !!---------------------------------------------------------------------- 1303 1278 ! 1304 1279 ! Max courant number for ext. grav. waves 1305 !1306 CALL wrk_alloc( jpi,jpj, zcu )1307 1280 ! 1308 1281 DO jj = 1, jpj … … 1371 1344 ENDIF 1372 1345 ! 1373 CALL wrk_dealloc( jpi,jpj, zcu )1374 !1375 1346 END SUBROUTINE dyn_spg_ts_init 1376 1347
Note: See TracChangeset
for help on using the changeset viewer.
