Changeset 12443 for NEMO/branches/2020/KERNEL-03_Storkey_Coward_RK3_stage2/src/OCE/DYN/dynspg_ts.F90
- Timestamp:
- 2020-02-24T14:00:21+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/KERNEL-03_Storkey_Coward_RK3_stage2/src/OCE/DYN/dynspg_ts.F90
r12424 r12443 72 72 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: un_adv , vn_adv !: Advection vel. at "now" barocl. step 73 73 ! 74 INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_ baro <= 2.5 nn_baro75 REAL(wp),SAVE :: r dtbt! Barotropic time step74 INTEGER, SAVE :: icycle ! Number of barotropic sub-steps for each internal step nn_e <= 2.5 nn_e 75 REAL(wp),SAVE :: rDt_e ! Barotropic time step 76 76 ! 77 77 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: wgtbtp1, wgtbtp2 ! 1st & 2nd weights used in time filtering of barotropic fields … … 102 102 ierr(:) = 0 103 103 ! 104 ALLOCATE( wgtbtp1(3*nn_ baro), wgtbtp2(3*nn_baro), zwz(jpi,jpj), STAT=ierr(1) )104 ALLOCATE( wgtbtp1(3*nn_e), wgtbtp2(3*nn_e), zwz(jpi,jpj), STAT=ierr(1) ) 105 105 IF( ln_dynvor_een .OR. ln_dynvor_eeT ) & 106 106 & ALLOCATE( ftnw(jpi,jpj) , ftne(jpi,jpj) , ftsw(jpi,jpj) , ftse(jpi,jpj), STAT=ierr(2) ) … … 185 185 ll_fw_start = .FALSE. 186 186 ! ! time offset in steps for bdy data update 187 IF( .NOT.ln_bt_fw ) THEN ; noffset = - nn_ baro187 IF( .NOT.ln_bt_fw ) THEN ; noffset = - nn_e 188 188 ELSE ; noffset = 0 189 189 ENDIF … … 303 303 IF( ln_bt_fw ) THEN ! Add wind forcing 304 304 DO_2D_00_00 305 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_r au0 * utau(ji,jj) * r1_hu(ji,jj,Kmm)306 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_r au0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm)305 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) 306 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) 307 307 END_2D 308 308 ELSE 309 zztmp = r1_r au0 * r1_2309 zztmp = r1_rho0 * r1_2 310 310 DO_2D_00_00 311 311 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) … … 320 320 ! ! --------------------------------------------------- ! 321 321 IF (ln_bt_fw) THEN ! FORWARD integration: use kt+1/2 fluxes (NOW+1/2) 322 zssh_frc(:,:) = r1_r au0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) )322 zssh_frc(:,:) = r1_rho0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) 323 323 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 fluxes (NOW) 324 zztmp = r1_r au0 * r1_2324 zztmp = r1_rho0 * r1_2 325 325 zssh_frc(:,:) = zztmp * ( emp(:,:) + emp_b(:,:) & 326 326 & - rnf(:,:) - rnf_b(:,:) & … … 429 429 ! Update tide potential at the beginning of current time substep 430 430 IF( ln_tide_pot .AND. ln_tide ) THEN 431 zt0substep = REAL(nsec_day, wp) - 0.5_wp*rn_Dt + (jn + noffset - 1) * rn_Dt / REAL(nn_ baro, wp)431 zt0substep = REAL(nsec_day, wp) - 0.5_wp*rn_Dt + (jn + noffset - 1) * rn_Dt / REAL(nn_e, wp) 432 432 CALL upd_tide(zt0substep, Kmm) 433 433 END IF … … 495 495 IF( .NOT.Agrif_Root() .AND. ln_bt_fw ) CALL agrif_dyn_ts_flux( jn, zhU, zhV ) 496 496 #endif 497 IF( ln_wd_il ) CALL wad_lmt_bt(zhU, zhV, sshn_e, zssh_frc, r dtbt) !!gm wad_lmt_bt use of lbc_lnk on zhU, zhV497 IF( ln_wd_il ) CALL wad_lmt_bt(zhU, zhV, sshn_e, zssh_frc, rDt_e) !!gm wad_lmt_bt use of lbc_lnk on zhU, zhV 498 498 499 499 IF( ln_wd_dl ) THEN ! un_e and vn_e are set to zero at faces where … … 510 510 DO_2D_00_00 511 511 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 512 ssha_e(ji,jj) = ( sshn_e(ji,jj) - r dtbt* ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj)512 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rDt_e * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) 513 513 END_2D 514 514 ! … … 552 552 ! 553 553 ! ! Surface pressure gradient 554 zldg = ( 1._wp - rn_ scal_load ) * grav ! local factor554 zldg = ( 1._wp - rn_load ) * grav ! local factor 555 555 DO_2D_00_00 556 556 zu_spg(ji,jj) = - zldg * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) … … 600 600 DO_2D_00_00 601 601 ua_e(ji,jj) = ( un_e(ji,jj) & 602 & + r dtbt* ( zu_spg(ji,jj) &602 & + rDt_e * ( zu_spg(ji,jj) & 603 603 & + zu_trd(ji,jj) & 604 604 & + zu_frc(ji,jj) ) & … … 606 606 607 607 va_e(ji,jj) = ( vn_e(ji,jj) & 608 & + r dtbt* ( zv_spg(ji,jj) &608 & + rDt_e * ( zv_spg(ji,jj) & 609 609 & + zv_trd(ji,jj) & 610 610 & + zv_frc(ji,jj) ) & … … 625 625 ! 626 626 ua_e(ji,jj) = ( hu_e (ji,jj) * un_e (ji,jj) & 627 & + r dtbt* ( zhu_bck * zu_spg (ji,jj) & !627 & + rDt_e * ( zhu_bck * zu_spg (ji,jj) & ! 628 628 & + zhup2_e(ji,jj) * zu_trd (ji,jj) & ! 629 629 & + hu(ji,jj,Kmm) * zu_frc (ji,jj) ) ) * z1_hu 630 630 ! 631 631 va_e(ji,jj) = ( hv_e (ji,jj) * vn_e (ji,jj) & 632 & + r dtbt* ( zhv_bck * zv_spg (ji,jj) & !632 & + rDt_e * ( zhv_bck * zv_spg (ji,jj) & ! 633 633 & + zhvp2_e(ji,jj) * zv_trd (ji,jj) & ! 634 634 & + hv(ji,jj,Kmm) * zv_frc (ji,jj) ) ) * z1_hv … … 638 638 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 639 639 DO_2D_00_00 640 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - r dtbt* zCdU_u(ji,jj) * hur_e(ji,jj))641 va_e(ji,jj) = va_e(ji,jj) /(1.0 - r dtbt* zCdU_v(ji,jj) * hvr_e(ji,jj))640 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rDt_e * zCdU_u(ji,jj) * hur_e(ji,jj)) 641 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rDt_e * zCdU_v(ji,jj) * hvr_e(ji,jj)) 642 642 END_2D 643 643 ENDIF … … 707 707 zvn_save = vn_adv(ji,jj) 708 708 ! ! apply the previously computed correction 709 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - atfp * un_bf(ji,jj) )710 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - atfp * vn_bf(ji,jj) )709 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - rn_atfp * un_bf(ji,jj) ) 710 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - rn_atfp * vn_bf(ji,jj) ) 711 711 ! ! Update corrective fluxes for next time step 712 un_bf(ji,jj) = atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) )713 vn_bf(ji,jj) = atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) )712 un_bf(ji,jj) = rn_atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) ) 713 vn_bf(ji,jj) = rn_atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) ) 714 714 ! ! Save integrated transport for next computation 715 715 ub2_b(ji,jj) = zun_save … … 809 809 LOGICAL, INTENT(in) :: ll_fw ! forward time splitting =.true. 810 810 INTEGER, INTENT(inout) :: jpit ! cycle length 811 REAL(wp), DIMENSION(3*nn_ baro), INTENT(inout) :: zwgt1, & ! Primary weights811 REAL(wp), DIMENSION(3*nn_e), INTENT(inout) :: zwgt1, & ! Primary weights 812 812 zwgt2 ! Secondary weights 813 813 … … 821 821 ! Set time index when averaged value is requested 822 822 IF (ll_fw) THEN 823 jic = nn_ baro823 jic = nn_e 824 824 ELSE 825 jic = 2 * nn_ baro825 jic = 2 * nn_e 826 826 ENDIF 827 827 … … 829 829 IF (ll_av) THEN 830 830 ! Define simple boxcar window for primary weights 831 ! (width = nn_ baro, centered around jic)831 ! (width = nn_e, centered around jic) 832 832 SELECT CASE ( nn_bt_flt ) 833 833 CASE( 0 ) ! No averaging … … 835 835 jpit = jic 836 836 837 CASE( 1 ) ! Boxcar, width = nn_ baro838 DO jn = 1, 3*nn_ baro839 za1 = ABS(float(jn-jic))/float(nn_ baro)837 CASE( 1 ) ! Boxcar, width = nn_e 838 DO jn = 1, 3*nn_e 839 za1 = ABS(float(jn-jic))/float(nn_e) 840 840 IF (za1 < 0.5_wp) THEN 841 841 zwgt1(jn) = 1._wp … … 844 844 ENDDO 845 845 846 CASE( 2 ) ! Boxcar, width = 2 * nn_ baro847 DO jn = 1, 3*nn_ baro848 za1 = ABS(float(jn-jic))/float(nn_ baro)846 CASE( 2 ) ! Boxcar, width = 2 * nn_e 847 DO jn = 1, 3*nn_e 848 za1 = ABS(float(jn-jic))/float(nn_e) 849 849 IF (za1 < 1._wp) THEN 850 850 zwgt1(jn) = 1._wp … … 976 976 977 977 ! Estimate number of iterations to satisfy a max courant number= rn_bt_cmax 978 IF( ln_bt_auto ) nn_ baro= CEILING( rn_Dt / rn_bt_cmax * zcmax)978 IF( ln_bt_auto ) nn_e = CEILING( rn_Dt / rn_bt_cmax * zcmax) 979 979 980 r dtbt = rn_Dt / REAL( nn_baro, wp )981 zcmax = zcmax * r dtbt980 rDt_e = rn_Dt / REAL( nn_e , wp ) 981 zcmax = zcmax * rDt_e 982 982 ! Print results 983 983 IF(lwp) WRITE(numout,*) … … 985 985 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' 986 986 IF( ln_bt_auto ) THEN 987 IF(lwp) WRITE(numout,*) ' ln_ts_auto =.true. Automatically set nn_ baro'987 IF(lwp) WRITE(numout,*) ' ln_ts_auto =.true. Automatically set nn_e ' 988 988 IF(lwp) WRITE(numout,*) ' Max. courant number allowed: ', rn_bt_cmax 989 989 ELSE 990 IF(lwp) WRITE(numout,*) ' ln_ts_auto=.false.: Use nn_ baro in namelist nn_baro = ', nn_baro990 IF(lwp) WRITE(numout,*) ' ln_ts_auto=.false.: Use nn_e in namelist nn_e = ', nn_e 991 991 ENDIF 992 992 993 993 IF(ln_bt_av) THEN 994 IF(lwp) WRITE(numout,*) ' ln_bt_av =.true. ==> Time averaging over nn_ barotime steps is on '994 IF(lwp) WRITE(numout,*) ' ln_bt_av =.true. ==> Time averaging over nn_e time steps is on ' 995 995 ELSE 996 996 IF(lwp) WRITE(numout,*) ' ln_bt_av =.false. => No time averaging of barotropic variables ' … … 1012 1012 SELECT CASE ( nn_bt_flt ) 1013 1013 CASE( 0 ) ; IF(lwp) WRITE(numout,*) ' Dirac' 1014 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = nn_ baro'1015 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = 2*nn_ baro'1014 CASE( 1 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = nn_e' 1015 CASE( 2 ) ; IF(lwp) WRITE(numout,*) ' Boxcar: width = 2*nn_e' 1016 1016 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for nn_bt_flt: should 0,1, or 2' ) 1017 1017 END SELECT 1018 1018 ! 1019 1019 IF(lwp) WRITE(numout,*) ' ' 1020 IF(lwp) WRITE(numout,*) ' nn_ baro = ', nn_baro1021 IF(lwp) WRITE(numout,*) ' Barotropic time step [s] is :', r dtbt1020 IF(lwp) WRITE(numout,*) ' nn_e = ', nn_e 1021 IF(lwp) WRITE(numout,*) ' Barotropic time step [s] is :', rDt_e 1022 1022 IF(lwp) WRITE(numout,*) ' Maximum Courant number is :', zcmax 1023 1023 ! … … 1031 1031 ENDIF 1032 1032 IF( zcmax>0.9_wp ) THEN 1033 CALL ctl_stop( 'dynspg_ts ERROR: Maximum Courant number is greater than 0.9: Inc. nn_ baro!' )1033 CALL ctl_stop( 'dynspg_ts ERROR: Maximum Courant number is greater than 0.9: Inc. nn_e !' ) 1034 1034 ENDIF 1035 1035 ! … … 1430 1430 ! 1431 1431 IF( ln_wd_il ) THEN ! W/D : use the "clipped" bottom friction !!gm explain WHY, please ! 1432 zztmp = -1._wp / r dtbt1432 zztmp = -1._wp / rDt_e 1433 1433 DO_2D_00_00 1434 1434 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( &
Note: See TracChangeset
for help on using the changeset viewer.
