Ignore:
Timestamp:
2020-02-28T16:55:11+01:00 (14 months ago)
Author:
davestorkey
Message:

Preparation for new timestepping scheme #2390.
Main changes:

  1. Initial euler timestep now handled in stp and not in TRA/DYN routines.
  2. Renaming of all timestep parameters. In summary, the namelist parameter is now rn_Dt and the current timestep is rDt (and rDt_ice, rDt_trc etc).
  3. Renaming of a few miscellaneous parameters, eg. atfp → rn_atfp (namelist parameter used everywhere) and rau0 → rho0.

This version gives bit-comparable results to the previous version of the trunk.

Location:
NEMO/trunk/src/OCE/DYN
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/DYN/dynatf.F90

    r12377 r12489  
    8787      !!             arrays to start the next time step: 
    8888      !!                (puu(Kmm),pvv(Kmm)) = (puu(Kmm),pvv(Kmm))  
    89       !!                                    + atfp [ (puu(Kbb),pvv(Kbb)) + (puu(Kaa),pvv(Kaa)) - 2 (puu(Kmm),pvv(Kmm)) ] 
     89      !!                                    + rn_atfp [ (puu(Kbb),pvv(Kbb)) + (puu(Kaa),pvv(Kaa)) - 2 (puu(Kmm),pvv(Kmm)) ] 
    9090      !!             Note that with flux form advection and non linear free surface, 
    9191      !!             the time filter is applied on thickness weighted velocity. 
     
    157157      ! 
    158158      IF( l_trddyn ) THEN             ! prepare the atf trend computation + some diagnostics 
    159          z1_2dt = 1._wp / (2. * rdt)        ! Euler or leap-frog time step  
    160          IF( neuler == 0 .AND. kt == nit000 )   z1_2dt = 1._wp / rdt 
    161159         ! 
    162160         !                                  ! Kinetic energy and Conversion 
     
    164162         ! 
    165163         IF( ln_dyn_trd ) THEN              ! 3D output: total momentum trends 
    166             zua(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) * z1_2dt 
    167             zva(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) * z1_2dt 
     164            zua(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) * r1_Dt 
     165            zva(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) * r1_Dt 
    168166            CALL iom_put( "utrd_tot", zua )        ! total momentum trends, except the asselin time filter 
    169167            CALL iom_put( "vtrd_tot", zva ) 
     
    178176      ! ------------------------------------------ 
    179177          
    180       IF( .NOT.( neuler == 0 .AND. kt == nit000 ) ) THEN    !* Leap-Frog : Asselin time filter  
     178      IF( .NOT. l_1st_euler ) THEN    !* Leap-Frog : Asselin time filter  
    181179         !                                ! =============! 
    182180         IF( ln_linssh ) THEN             ! Fixed volume ! 
    183181            !                             ! =============! 
    184182            DO_3D_11_11( 1, jpkm1 ) 
    185                puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 
    186                pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 
     183               puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + rn_atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 
     184               pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + rn_atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 
    187185            END_3D 
    188186            !                             ! ================! 
     
    193191            ALLOCATE( ze3t_f(jpi,jpj,jpk), zwfld(jpi,jpj) ) 
    194192            DO jk = 1, jpkm1 
    195                ze3t_f(:,:,jk) = pe3t(:,:,jk,Kmm) + atfp * ( pe3t(:,:,jk,Kbb) - 2._wp * pe3t(:,:,jk,Kmm) + pe3t(:,:,jk,Kaa) ) 
     193               ze3t_f(:,:,jk) = pe3t(:,:,jk,Kmm) + rn_atfp * ( pe3t(:,:,jk,Kbb) - 2._wp * pe3t(:,:,jk,Kmm) + pe3t(:,:,jk,Kaa) ) 
    196194            END DO 
    197195            ! Add volume filter correction: compatibility with tracer advection scheme 
    198196            ! => time filter + conservation correction 
    199             zcoef = atfp * rdt * r1_rau0 
     197            zcoef = rn_atfp * rn_Dt * r1_rho0 
    200198            zwfld(:,:) = emp_b(:,:) - emp(:,:) 
    201199            IF ( ln_rnf ) zwfld(:,:) =  zwfld(:,:) - ( rnf_b(:,:) - rnf(:,:) ) 
     
    209207            !     to manage rnf, isf and possibly in the futur icb, tide water glacier (...) 
    210208            !     ...(kt, coef, ktop, kbot, hz, fwf_b, fwf) 
    211             IF ( ln_isf ) CALL isf_dynatf( kt, Kmm, ze3t_f, atfp * rdt ) 
     209            IF ( ln_isf ) CALL isf_dynatf( kt, Kmm, ze3t_f, rn_atfp * rn_Dt ) 
    212210            ! 
    213211            pe3t(:,:,1:jpkm1,Kmm) = ze3t_f(:,:,1:jpkm1)        ! filtered scale factor at T-points 
     
    218216               CALL dom_vvl_interpol( pe3t(:,:,:,Kmm), pe3v(:,:,:,Kmm), 'V' ) 
    219217               DO_3D_11_11( 1, jpkm1 ) 
    220                   puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 
    221                   pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 
     218                  puu(ji,jj,jk,Kmm) = puu(ji,jj,jk,Kmm) + rn_atfp * ( puu(ji,jj,jk,Kbb) - 2._wp * puu(ji,jj,jk,Kmm) + puu(ji,jj,jk,Kaa) ) 
     219                  pvv(ji,jj,jk,Kmm) = pvv(ji,jj,jk,Kmm) + rn_atfp * ( pvv(ji,jj,jk,Kbb) - 2._wp * pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk,Kaa) ) 
    222220               END_3D 
    223221               ! 
     
    236234                  zve3b = pe3v(ji,jj,jk,Kbb) * pvv(ji,jj,jk,Kbb) 
    237235                  ! 
    238                   puu(ji,jj,jk,Kmm) = ( zue3n + atfp * ( zue3b - 2._wp * zue3n  + zue3a ) ) / ze3u_f(ji,jj,jk) 
    239                   pvv(ji,jj,jk,Kmm) = ( zve3n + atfp * ( zve3b - 2._wp * zve3n  + zve3a ) ) / ze3v_f(ji,jj,jk) 
     236                  puu(ji,jj,jk,Kmm) = ( zue3n + rn_atfp * ( zue3b - 2._wp * zue3n  + zue3a ) ) / ze3u_f(ji,jj,jk) 
     237                  pvv(ji,jj,jk,Kmm) = ( zve3n + rn_atfp * ( zve3b - 2._wp * zve3n  + zve3a ) ) / ze3v_f(ji,jj,jk) 
    240238               END_3D 
    241239               pe3u(:,:,1:jpkm1,Kmm) = ze3u_f(:,:,1:jpkm1)   
     
    263261         ENDIF 
    264262         ! 
    265       ENDIF ! neuler /= 0 
     263      ENDIF ! .NOT. l_1st_euler 
    266264      ! 
    267265      ! Set "now" and "before" barotropic velocities for next time step: 
  • NEMO/trunk/src/OCE/DYN/dynspg.F90

    r12377 r12489  
    6767      !!              ln_apr_dyn=T : the atmospheric pressure forcing is applied  
    6868      !!             as the gradient of the inverse barometer ssh: 
    69       !!                apgu = - 1/rau0 di[apr] = 0.5*grav di[ssh_ib+ssh_ibb] 
    70       !!                apgv = - 1/rau0 dj[apr] = 0.5*grav dj[ssh_ib+ssh_ibb] 
    71       !!             Note that as all external forcing a time averaging over a two rdt 
     69      !!                apgu = - 1/rho0 di[apr] = 0.5*grav di[ssh_ib+ssh_ibb] 
     70      !!                apgv = - 1/rho0 dj[apr] = 0.5*grav dj[ssh_ib+ssh_ibb] 
     71      !!             Note that as all external forcing a time averaging over a two rn_Dt 
    7272      !!             period is used to prevent the divergence of odd and even time step. 
    7373      !!---------------------------------------------------------------------- 
     
    7878      ! 
    7979      INTEGER  ::   ji, jj, jk                   ! dummy loop indices 
    80       REAL(wp) ::   z2dt, zg_2, zintp, zgrau0r, zld   ! local scalars 
     80      REAL(wp) ::   z2dt, zg_2, zintp, zgrho0r, zld   ! local scalars 
    8181      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   zpice 
    8282      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdu, ztrdv 
     
    114114            ! 
    115115            ! Update tide potential at the beginning of current time step 
    116             zt0step = REAL(nsec_day, wp)-0.5_wp*rdt 
     116            zt0step = REAL(nsec_day, wp)-0.5_wp*rn_Dt 
    117117            CALL upd_tide(zt0step, Kmm) 
    118118            ! 
     
    134134            ALLOCATE( zpice(jpi,jpj) ) 
    135135            zintp = REAL( MOD( kt-1, nn_fsbc ) ) / REAL( nn_fsbc ) 
    136             zgrau0r     = - grav * r1_rau0 
    137             zpice(:,:) = (  zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:)  ) * zgrau0r 
     136            zgrho0r     = - grav * r1_rho0 
     137            zpice(:,:) = (  zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:)  ) * zgrho0r 
    138138            DO_2D_00_00 
    139139               spgu(ji,jj) = spgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj) 
     
    183183      NAMELIST/namdyn_spg/ ln_dynspg_exp       , ln_dynspg_ts,   & 
    184184      &                    ln_bt_fw, ln_bt_av  , ln_bt_auto  ,   & 
    185       &                    nn_baro , rn_bt_cmax, nn_bt_flt, rn_bt_alpha 
     185      &                    nn_e , rn_bt_cmax, nn_bt_flt, rn_bt_alpha 
    186186      !!---------------------------------------------------------------------- 
    187187      ! 
     
    222222      ! 
    223223      IF( nspg == np_TS ) THEN   ! split-explicit scheme initialisation 
    224          CALL dyn_spg_ts_init          ! do it first: set nn_baro used to allocate some arrays later on 
     224         CALL dyn_spg_ts_init          ! do it first: set nn_e used to allocate some arrays later on 
    225225      ENDIF 
    226226      ! 
  • NEMO/trunk/src/OCE/DYN/dynspg_exp.F90

    r12377 r12489  
    4949      !!              momentum trend the surface pressure gradient : 
    5050      !!                      (uu(rhs),vv(rhs)) = (uu(rhs),vv(rhs)) + (spgu,spgv) 
    51       !!              where spgu = -1/rau0 d/dx(ps) = -g/e1u di( ssh(now) ) 
    52       !!                    spgv = -1/rau0 d/dy(ps) = -g/e2v dj( ssh(now) ) 
     51      !!              where spgu = -1/rho0 d/dx(ps) = -g/e1u di( ssh(now) ) 
     52      !!                    spgv = -1/rho0 d/dy(ps) = -g/e2v dj( ssh(now) ) 
    5353      !! 
    5454      !! ** Action :   (puu(:,:,:,Krhs),pvv(:,:,:,Krhs))   trend of horizontal velocity increased by  
  • NEMO/trunk/src/OCE/DYN/dynspg_ts.F90

    r12377 r12489  
    7272   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   un_adv , vn_adv   !: Advection vel. at "now" barocl. step 
    7373   ! 
    74    INTEGER, SAVE :: icycle      ! Number of barotropic sub-steps for each internal step nn_baro <= 2.5 nn_baro 
    75    REAL(wp),SAVE :: rdtbt       ! Barotropic time step 
     74   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 
    7676   ! 
    7777   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)   ::   wgtbtp1, wgtbtp2   ! 1st & 2nd weights used in time filtering of barotropic fields 
     
    102102      ierr(:) = 0 
    103103      ! 
    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) ) 
    105105      IF( ln_dynvor_een .OR. ln_dynvor_eeT )   & 
    106106         &     ALLOCATE( ftnw(jpi,jpj) , ftne(jpi,jpj) , ftsw(jpi,jpj) , ftse(jpi,jpj), STAT=ierr(2)   ) 
     
    150150      LOGICAL  ::   ll_init               ! =T : special startup of 2d equations 
    151151      INTEGER  ::   noffset               ! local integers  : time offset for bdy update 
    152       REAL(wp) ::   r1_2dt_b, z1_hu, z1_hv          ! local scalars 
     152      REAL(wp) ::   r1_Dt_b, z1_hu, z1_hv          ! local scalars 
    153153      REAL(wp) ::   za0, za1, za2, za3              !   -      - 
    154154      REAL(wp) ::   zztmp, zldg               !   -      - 
     
    180180!     zwdramp = 1._wp / (rn_wdmin2 - rn_wdmin1) ! more general ramp 
    181181      !                                         ! inverse of baroclinic time step  
    182       IF( kt == nit000 .AND. neuler == 0 ) THEN   ;   r1_2dt_b = 1._wp / (         rdt ) 
    183       ELSE                                        ;   r1_2dt_b = 1._wp / ( 2._wp * rdt ) 
    184       ENDIF 
     182      r1_Dt_b = 1._wp / rDt  
    185183      ! 
    186184      ll_init     = ln_bt_av                    ! if no time averaging, then no specific restart  
    187185      ll_fw_start = .FALSE. 
    188186      !                                         ! time offset in steps for bdy data update 
    189       IF( .NOT.ln_bt_fw ) THEN   ;   noffset = - nn_baro 
     187      IF( .NOT.ln_bt_fw ) THEN   ;   noffset = - nn_e 
    190188      ELSE                       ;   noffset =   0  
    191189      ENDIF 
     
    198196         IF(lwp) WRITE(numout,*) 
    199197         ! 
    200          IF( neuler == 0 )   ll_init=.TRUE. 
    201          ! 
    202          IF( ln_bt_fw .OR. neuler == 0 ) THEN 
     198         IF( l_1st_euler )   ll_init=.TRUE. 
     199         ! 
     200         IF( ln_bt_fw .OR. l_1st_euler ) THEN 
    203201            ll_fw_start =.TRUE. 
    204202            noffset     = 0 
     
    209207         CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 
    210208         ! 
    211       ENDIF 
    212       ! 
    213       ! If forward start at previous time step, and centered integration,  
    214       ! then update averaging weights: 
    215       IF (.NOT.ln_bt_fw .AND.( neuler==0 .AND. kt==nit000+1 ) ) THEN 
    216          ll_fw_start=.FALSE. 
    217          CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 
    218       ENDIF 
    219       ! 
    220                            
     209      ELSEIF( kt == nit000 + 1 ) THEN           !* initialisation 2nd time-step 
     210         ! 
     211         IF( .NOT.ln_bt_fw ) THEN 
     212            ! If we did an Euler timestep on the first timestep we need to reset ll_fw_start 
     213            ! and the averaging weights. We don't have an easy way of telling whether we did 
     214            ! an Euler timestep on the first timestep (because l_1st_euler is reset to .false. 
     215            ! at the end of the first timestep) so just do this in all cases.  
     216            ll_fw_start = .FALSE. 
     217            CALL ts_wgt( ln_bt_av, ll_fw_start, icycle, wgtbtp1, wgtbtp2 ) 
     218         ENDIF 
     219         ! 
     220      ENDIF 
     221      ! 
    221222      ! ----------------------------------------------------------------------------- 
    222223      !  Phase 1 : Coupling between general trend and barotropic estimates (1st step) 
     
    302303      IF( ln_bt_fw ) THEN                        ! Add wind forcing 
    303304         DO_2D_00_00 
    304             zu_frc(ji,jj) =  zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) 
    305             zv_frc(ji,jj) =  zv_frc(ji,jj) + r1_rau0 * 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) 
    306307         END_2D 
    307308      ELSE 
    308          zztmp = r1_rau0 * r1_2 
     309         zztmp = r1_rho0 * r1_2 
    309310         DO_2D_00_00 
    310311            zu_frc(ji,jj) =  zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) 
     
    319320      !                                   ! ---------------------------------------------------  ! 
    320321      IF (ln_bt_fw) THEN                          ! FORWARD integration: use kt+1/2 fluxes (NOW+1/2) 
    321          zssh_frc(:,:) = r1_rau0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) 
     322         zssh_frc(:,:) = r1_rho0 * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) ) 
    322323      ELSE                                        ! CENTRED integration: use kt-1/2 + kt+1/2 fluxes (NOW) 
    323          zztmp = r1_rau0 * r1_2 
     324         zztmp = r1_rho0 * r1_2 
    324325         zssh_frc(:,:) = zztmp * (  emp(:,:)        + emp_b(:,:)                    & 
    325326                &                 - rnf(:,:)        - rnf_b(:,:)                    & 
     
    428429         ! Update tide potential at the beginning of current time substep 
    429430         IF( ln_tide_pot .AND. ln_tide ) THEN 
    430             zt0substep = REAL(nsec_day, wp) - 0.5_wp*rdt + (jn + noffset - 1) * rdt / REAL(nn_baro, wp) 
     431            zt0substep = REAL(nsec_day, wp) - 0.5_wp*rn_Dt + (jn + noffset - 1) * rn_Dt / REAL(nn_e, wp) 
    431432            CALL upd_tide(zt0substep, Kmm) 
    432433         END IF 
     
    494495         IF( .NOT.Agrif_Root() .AND. ln_bt_fw ) CALL agrif_dyn_ts_flux( jn, zhU, zhV ) 
    495496#endif 
    496          IF( ln_wd_il )   CALL wad_lmt_bt(zhU, zhV, sshn_e, zssh_frc, rdtbt)    !!gm wad_lmt_bt use of lbc_lnk on zhU, zhV 
     497         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 
    497498 
    498499         IF( ln_wd_dl ) THEN           ! un_e and vn_e are set to zero at faces where  
     
    509510         DO_2D_00_00 
    510511            zhdiv = (   zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1)   ) * r1_e1e2t(ji,jj) 
    511             ssha_e(ji,jj) = (  sshn_e(ji,jj) - rdtbt * ( 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) 
    512513         END_2D 
    513514         ! 
     
    599600            DO_2D_00_00 
    600601               ua_e(ji,jj) = (                                 un_e(ji,jj)   &  
    601                          &     + rdtbt * (                   zu_spg(ji,jj)   & 
     602                         &     + rDt_e * (                   zu_spg(ji,jj)   & 
    602603                         &                                 + zu_trd(ji,jj)   & 
    603604                         &                                 + zu_frc(ji,jj) ) &  
     
    605606 
    606607               va_e(ji,jj) = (                                 vn_e(ji,jj)   & 
    607                          &     + rdtbt * (                   zv_spg(ji,jj)   & 
     608                         &     + rDt_e * (                   zv_spg(ji,jj)   & 
    608609                         &                                 + zv_trd(ji,jj)   & 
    609610                         &                                 + zv_frc(ji,jj) ) & 
     
    624625               ! 
    625626               ua_e(ji,jj) = (               hu_e  (ji,jj) *   un_e (ji,jj)      &  
    626                     &            + rdtbt * (  zhu_bck        * zu_spg (ji,jj)  &   ! 
     627                    &            + rDt_e * (  zhu_bck        * zu_spg (ji,jj)  &   ! 
    627628                    &                       + zhup2_e(ji,jj) * zu_trd (ji,jj)  &   ! 
    628629                    &                       +  hu(ji,jj,Kmm) * zu_frc (ji,jj)  )   ) * z1_hu 
    629630               ! 
    630631               va_e(ji,jj) = (               hv_e  (ji,jj) *   vn_e (ji,jj)      & 
    631                     &            + rdtbt * (  zhv_bck        * zv_spg (ji,jj)  &   ! 
     632                    &            + rDt_e * (  zhv_bck        * zv_spg (ji,jj)  &   ! 
    632633                    &                       + zhvp2_e(ji,jj) * zv_trd (ji,jj)  &   ! 
    633634                    &                       +  hv(ji,jj,Kmm) * zv_frc (ji,jj)  )   ) * z1_hv 
     
    637638         IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 
    638639            DO_2D_00_00 
    639                   ua_e(ji,jj) =  ua_e(ji,jj) /(1.0 -   rdtbt * zCdU_u(ji,jj) * hur_e(ji,jj)) 
    640                   va_e(ji,jj) =  va_e(ji,jj) /(1.0 -   rdtbt * 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)) 
    641642            END_2D 
    642643         ENDIF 
     
    701702      ! Set advection velocity correction: 
    702703      IF (ln_bt_fw) THEN 
    703          IF( .NOT.( kt == nit000 .AND. neuler==0 ) ) THEN 
     704         IF( .NOT.( kt == nit000 .AND. l_1st_euler ) ) THEN 
    704705            DO_2D_11_11 
    705706               zun_save = un_adv(ji,jj) 
    706707               zvn_save = vn_adv(ji,jj) 
    707708               !                          ! apply the previously computed correction  
    708                un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - atfp * un_bf(ji,jj) ) 
    709                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) ) 
    710711               !                          ! Update corrective fluxes for next time step 
    711                un_bf(ji,jj)  = atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) ) 
    712                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) ) 
    713714               !                          ! Save integrated transport for next computation 
    714715               ub2_b(ji,jj) = zun_save 
     
    728729      IF( ln_dynadv_vec .OR. ln_linssh ) THEN 
    729730         DO jk=1,jpkm1 
    730             puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) ) * r1_2dt_b 
    731             pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) ) * r1_2dt_b 
     731            puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) ) * r1_Dt_b 
     732            pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) ) * r1_Dt_b 
    732733         END DO 
    733734      ELSE 
     
    744745         ! 
    745746         DO jk=1,jpkm1 
    746             puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_2dt_b 
    747             pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_2dt_b 
     747            puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_Dt_b 
     748            pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_Dt_b 
    748749         END DO 
    749750         ! Save barotropic velocities not transport: 
     
    808809      LOGICAL, INTENT(in) ::   ll_fw      ! forward time splitting =.true. 
    809810      INTEGER, INTENT(inout) :: jpit      ! cycle length     
    810       REAL(wp), DIMENSION(3*nn_baro), INTENT(inout) ::   zwgt1, & ! Primary weights 
     811      REAL(wp), DIMENSION(3*nn_e), INTENT(inout) ::   zwgt1, & ! Primary weights 
    811812                                                         zwgt2    ! Secondary weights 
    812813       
     
    820821      ! Set time index when averaged value is requested 
    821822      IF (ll_fw) THEN  
    822          jic = nn_baro 
     823         jic = nn_e 
    823824      ELSE 
    824          jic = 2 * nn_baro 
     825         jic = 2 * nn_e 
    825826      ENDIF 
    826827 
     
    828829      IF (ll_av) THEN 
    829830           ! Define simple boxcar window for primary weights  
    830            ! (width = nn_baro, centered around jic)      
     831           ! (width = nn_e, centered around jic)      
    831832         SELECT CASE ( nn_bt_flt ) 
    832833              CASE( 0 )  ! No averaging 
     
    834835                 jpit = jic 
    835836 
    836               CASE( 1 )  ! Boxcar, width = nn_baro 
    837                  DO jn = 1, 3*nn_baro 
    838                     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)  
    839840                    IF (za1 < 0.5_wp) THEN 
    840841                      zwgt1(jn) = 1._wp 
     
    843844                 ENDDO 
    844845 
    845               CASE( 2 )  ! Boxcar, width = 2 * nn_baro 
    846                  DO jn = 1, 3*nn_baro 
    847                     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)  
    848849                    IF (za1 < 1._wp) THEN 
    849850                      zwgt1(jn) = 1._wp 
     
    889890      IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise  
    890891         !                                   ! --------------- 
    891          IF( ln_rstart .AND. ln_bt_fw .AND. (neuler/=0) ) THEN    !* Read the restart file 
     892         IF( ln_rstart .AND. ln_bt_fw .AND. (.NOT.l_1st_euler) ) THEN    !* Read the restart file 
    892893            CALL iom_get( numror, jpdom_autoglo, 'ub2_b'  , ub2_b  (:,:), ldxios = lrxios )    
    893894            CALL iom_get( numror, jpdom_autoglo, 'vb2_b'  , vb2_b  (:,:), ldxios = lrxios )  
     
    975976 
    976977      ! Estimate number of iterations to satisfy a max courant number= rn_bt_cmax 
    977       IF( ln_bt_auto )   nn_baro = CEILING( rdt / rn_bt_cmax * zcmax) 
     978      IF( ln_bt_auto )   nn_e = CEILING( rn_Dt / rn_bt_cmax * zcmax) 
    978979       
    979       rdtbt = rdt / REAL( nn_baro , wp ) 
    980       zcmax = zcmax * rdtbt 
     980      rDt_e = rn_Dt / REAL( nn_e , wp ) 
     981      zcmax = zcmax * rDt_e 
    981982      ! Print results 
    982983      IF(lwp) WRITE(numout,*) 
     
    984985      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    985986      IF( ln_bt_auto ) THEN 
    986          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 ' 
    987988         IF(lwp) WRITE(numout,*) '     Max. courant number allowed: ', rn_bt_cmax 
    988989      ELSE 
    989          IF(lwp) WRITE(numout,*) '     ln_ts_auto=.false.: Use nn_baro in namelist   nn_baro = ', nn_baro 
     990         IF(lwp) WRITE(numout,*) '     ln_ts_auto=.false.: Use nn_e in namelist   nn_e = ', nn_e 
    990991      ENDIF 
    991992 
    992993      IF(ln_bt_av) THEN 
    993          IF(lwp) WRITE(numout,*) '     ln_bt_av =.true.  ==> Time averaging over nn_baro time steps is on ' 
     994         IF(lwp) WRITE(numout,*) '     ln_bt_av =.true.  ==> Time averaging over nn_e time steps is on ' 
    994995      ELSE 
    995996         IF(lwp) WRITE(numout,*) '     ln_bt_av =.false. => No time averaging of barotropic variables ' 
     
    10111012      SELECT CASE ( nn_bt_flt ) 
    10121013         CASE( 0 )      ;   IF(lwp) WRITE(numout,*) '           Dirac' 
    1013          CASE( 1 )      ;   IF(lwp) WRITE(numout,*) '           Boxcar: width = nn_baro' 
    1014          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'  
    10151016         CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for nn_bt_flt: should 0,1, or 2' ) 
    10161017      END SELECT 
    10171018      ! 
    10181019      IF(lwp) WRITE(numout,*) ' ' 
    1019       IF(lwp) WRITE(numout,*) '     nn_baro = ', nn_baro 
    1020       IF(lwp) WRITE(numout,*) '     Barotropic time step [s] is :', rdtbt 
     1020      IF(lwp) WRITE(numout,*) '     nn_e = ', nn_e 
     1021      IF(lwp) WRITE(numout,*) '     Barotropic time step [s] is :', rDt_e 
    10211022      IF(lwp) WRITE(numout,*) '     Maximum Courant number is   :', zcmax 
    10221023      ! 
     
    10301031      ENDIF 
    10311032      IF( zcmax>0.9_wp ) THEN 
    1032          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 !' )           
    10331034      ENDIF 
    10341035      ! 
     
    14291430      ! 
    14301431      IF( ln_wd_il ) THEN      ! W/D : use the "clipped" bottom friction   !!gm   explain WHY, please ! 
    1431          zztmp = -1._wp / rdtbt 
     1432         zztmp = -1._wp / rDt_e 
    14321433         DO_2D_00_00 
    14331434            pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) *  wdrampu(ji,jj) * MAX(                                 &  
  • NEMO/trunk/src/OCE/DYN/dynzdf.F90

    r12377 r12489  
    9292         ENDIF 
    9393      ENDIF 
    94       !                             !* set time step 
    95       IF( neuler == 0 .AND. kt == nit000     ) THEN   ;   r2dt =      rdt   ! = rdt (restart with Euler time stepping) 
    96       ELSEIF(               kt <= nit000 + 1 ) THEN   ;   r2dt = 2. * rdt   ! = 2 rdt (leapfrog) 
    97       ENDIF 
    98       ! 
    9994      !                             !* explicit top/bottom drag case 
    10095      IF( .NOT.ln_drgimp )   CALL zdf_drg_exp( kt, Kmm, puu(:,:,:,Kbb), pvv(:,:,:,Kbb), puu(:,:,:,Krhs), pvv(:,:,:,Krhs) )  ! add top/bottom friction trend to (puu(Kaa),pvv(Kaa)) 
     
    112107      IF( ln_dynadv_vec .OR. ln_linssh ) THEN   ! applied on velocity 
    113108         DO jk = 1, jpkm1 
    114             puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kbb) + r2dt * puu(:,:,jk,Krhs) ) * umask(:,:,jk) 
    115             pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kbb) + r2dt * pvv(:,:,jk,Krhs) ) * vmask(:,:,jk) 
     109            puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kbb) + rDt * puu(:,:,jk,Krhs) ) * umask(:,:,jk) 
     110            pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kbb) + rDt * pvv(:,:,jk,Krhs) ) * vmask(:,:,jk) 
    116111         END DO 
    117112      ELSE                                      ! applied on thickness weighted velocity 
    118113         DO jk = 1, jpkm1 
    119114            puu(:,:,jk,Kaa) = (         e3u(:,:,jk,Kbb) * puu(:,:,jk,Kbb)  & 
    120                &          + r2dt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs)  ) / e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     115               &          + rDt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs)  ) / e3u(:,:,jk,Kaa) * umask(:,:,jk) 
    121116            pvv(:,:,jk,Kaa) = (         e3v(:,:,jk,Kbb) * pvv(:,:,jk,Kbb)  & 
    122                &          + r2dt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs)  ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
     117               &          + rDt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs)  ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    123118         END DO 
    124119      ENDIF 
     
    138133            ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 
    139134            ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 
    140             puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + r2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
    141             pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + r2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
     135            puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
     136            pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
    142137         END_2D 
    143138         IF( ln_isfcav ) THEN    ! Ocean cavities (ISF) 
     
    147142               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 
    148143               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 
    149                puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
    150                pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + r2dt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
     144               puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
     145               pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
    151146            END_2D 
    152147         END IF 
     
    156151      ! 
    157152      !                    !* Matrix construction 
    158       zdt = r2dt * 0.5 
     153      zdt = rDt * 0.5 
    159154      IF( ln_zad_Aimp ) THEN   !! 
    160155         SELECT CASE( nldf_dyn ) 
     
    232227            iku = mbku(ji,jj)       ! ocean bottom level at u- and v-points 
    233228            ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
    234             zwd(ji,jj,iku) = zwd(ji,jj,iku) - r2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 
     229            zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 
    235230         END_2D 
    236231         IF ( ln_isfcav ) THEN   ! top friction (always implicit) 
     
    239234               iku = miku(ji,jj)       ! ocean top level at u- and v-points  
    240235               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
    241                zwd(ji,jj,iku) = zwd(ji,jj,iku) - r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 
     236               zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 
    242237            END_2D 
    243238         END IF 
     
    265260      DO_2D_00_00 
    266261         ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa)  
    267          puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + r2dt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) )   & 
    268             &                                      / ( ze3ua * rau0 ) * umask(ji,jj,1)  
     262         puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) )   & 
     263            &                                      / ( ze3ua * rho0 ) * umask(ji,jj,1)  
    269264      END_2D 
    270265      DO_3D_00_00( 2, jpkm1 ) 
     
    282277      ! 
    283278      !                       !* Matrix construction 
    284       zdt = r2dt * 0.5 
     279      zdt = rDt * 0.5 
    285280      IF( ln_zad_Aimp ) THEN   !! 
    286281         SELECT CASE( nldf_dyn ) 
     
    357352            ikv = mbkv(ji,jj)       ! (deepest ocean u- and v-points) 
    358353            ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
    359             zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - r2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va            
     354            zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va            
    360355         END_2D 
    361356         IF ( ln_isfcav ) THEN 
     
    363358               ikv = mikv(ji,jj)       ! (first wet ocean u- and v-points) 
    364359               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
    365                zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - r2dt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 
     360               zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 
    366361            END_2D 
    367362         ENDIF 
     
    389384      DO_2D_00_00 
    390385         ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa)  
    391          pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + r2dt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) )   & 
    392             &                                      / ( ze3va * rau0 ) * vmask(ji,jj,1)  
     386         pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) )   & 
     387            &                                      / ( ze3va * rho0 ) * vmask(ji,jj,1)  
    393388      END_2D 
    394389      DO_3D_00_00( 2, jpkm1 ) 
     
    404399      ! 
    405400      IF( l_trddyn )   THEN                      ! save the vertical diffusive trends for further diagnostics 
    406          ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / r2dt - ztrdu(:,:,:) 
    407          ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / r2dt - ztrdv(:,:,:) 
     401         ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / rDt - ztrdu(:,:,:) 
     402         ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / rDt - ztrdv(:,:,:) 
    408403         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt, Kmm ) 
    409404         DEALLOCATE( ztrdu, ztrdv )  
  • NEMO/trunk/src/OCE/DYN/sshwzv.F90

    r12377 r12489  
    7575      REAL(wp), DIMENSION(jpi,jpj,jpt), INTENT(inout) ::   pssh           ! sea-surface height 
    7676      !  
    77       INTEGER  ::   jk            ! dummy loop indice 
    78       REAL(wp) ::   z2dt, zcoef   ! local scalars 
     77      INTEGER  ::   jk      ! dummy loop index 
     78      REAL(wp) ::   zcoef   ! local scalar 
    7979      REAL(wp), DIMENSION(jpi,jpj) ::   zhdiv   ! 2D workspace 
    8080      !!---------------------------------------------------------------------- 
     
    8888      ENDIF 
    8989      ! 
    90       z2dt = 2._wp * rdt                          ! set time step size (Euler/Leapfrog) 
    91       IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt 
    92       zcoef = 0.5_wp * r1_rau0 
     90      zcoef = 0.5_wp * r1_rho0 
    9391 
    9492      !                                           !------------------------------! 
     
    9694      !                                           !------------------------------! 
    9795      IF(ln_wd_il) THEN 
    98          CALL wad_lmt(pssh(:,:,Kbb), zcoef * (emp_b(:,:) + emp(:,:)), z2dt, Kmm, uu, vv ) 
     96         CALL wad_lmt(pssh(:,:,Kbb), zcoef * (emp_b(:,:) + emp(:,:)), rDt, Kmm, uu, vv ) 
    9997      ENDIF 
    10098 
     
    109107      ! compute the vertical velocity which can be used to compute the non-linear terms of the momentum equations. 
    110108      !  
    111       pssh(:,:,Kaa) = (  pssh(:,:,Kbb) - z2dt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * ssmask(:,:) 
     109      pssh(:,:,Kaa) = (  pssh(:,:,Kbb) - rDt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * ssmask(:,:) 
    112110      ! 
    113111#if defined key_agrif 
     
    152150      ! 
    153151      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    154       REAL(wp) ::   z1_2dt       ! local scalars 
    155152      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zhdiv 
    156153      !!---------------------------------------------------------------------- 
     
    168165      !                                           !     Now Vertical Velocity    ! 
    169166      !                                           !------------------------------! 
    170       z1_2dt = 1. / ( 2. * rdt )                         ! set time step size (Euler/Leapfrog) 
    171       IF( neuler == 0 .AND. kt == nit000 )   z1_2dt = 1. / rdt 
    172167      ! 
    173168      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN      ! z_tilde and layer cases 
     
    187182            ! computation of w 
    188183            pww(:,:,jk) = pww(:,:,jk+1) - (  e3t(:,:,jk,Kmm) * hdiv(:,:,jk) + zhdiv(:,:,jk)    & 
    189                &                         + z1_2dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )     ) * tmask(:,:,jk) 
     184               &                         + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )     ) * tmask(:,:,jk) 
    190185         END DO 
    191186         !          IF( ln_vvl_layer ) pww(:,:,:) = 0.e0 
     
    195190            ! computation of w 
    196191            pww(:,:,jk) = pww(:,:,jk+1) - (  e3t(:,:,jk,Kmm) * hdiv(:,:,jk)                 & 
    197                &                         + z1_2dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )  ) * tmask(:,:,jk) 
     192               &                         + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )  ) * tmask(:,:,jk) 
    198193         END DO 
    199194      ENDIF 
     
    227222      !! ** Method  : - apply Asselin time fiter to now ssh (excluding the forcing 
    228223      !!              from the filter, see Leclair and Madec 2010) and swap : 
    229       !!                pssh(:,:,Kmm) = pssh(:,:,Kaa) + atfp * ( pssh(:,:,Kbb) -2 pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
    230       !!                            - atfp * rdt * ( emp_b - emp ) / rau0 
     224      !!                pssh(:,:,Kmm) = pssh(:,:,Kaa) + rn_atfp * ( pssh(:,:,Kbb) -2 pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
     225      !!                            - rn_atfp * rn_Dt * ( emp_b - emp ) / rho0 
    231226      !! 
    232227      !! ** action  : - pssh(:,:,Kmm) time filtered 
     
    249244      ENDIF 
    250245      !              !==  Euler time-stepping: no filter, just swap  ==! 
    251       IF ( .NOT.( neuler == 0 .AND. kt == nit000 ) ) THEN   ! Only do time filtering for leapfrog timesteps 
     246      IF ( .NOT.( l_1st_euler ) ) THEN   ! Only do time filtering for leapfrog timesteps 
    252247         !                                                  ! filtered "now" field 
    253          pssh(:,:,Kmm) = pssh(:,:,Kmm) + atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
     248         pssh(:,:,Kmm) = pssh(:,:,Kmm) + rn_atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
    254249         IF( .NOT.ln_linssh ) THEN                          ! "now" <-- with forcing removed 
    255             zcoef = atfp * rdt * r1_rau0 
     250            zcoef = rn_atfp * rn_Dt * r1_rho0 
    256251            pssh(:,:,Kmm) = pssh(:,:,Kmm) - zcoef * (     emp_b(:,:) - emp   (:,:)   & 
    257252               &                             - rnf_b(:,:)        + rnf   (:,:)       & 
     
    260255 
    261256            ! ice sheet coupling 
    262             IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) pssh(:,:,Kbb) = pssh(:,:,Kbb) - atfp * rdt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 
     257            IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) pssh(:,:,Kbb) = pssh(:,:,Kbb) - rn_atfp * rn_Dt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 
    263258 
    264259         ENDIF 
     
    311306         DO_3D_00_00( 1, jpkm1 ) 
    312307            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    313             ! 2*rdt and not r2dt (for restartability) 
    314             Cu_adv(ji,jj,jk) = 2._wp * rdt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )                       &   
     308            ! 2*rn_Dt and not rDt (for restartability) 
     309            Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )                       &   
    315310               &                             + ( MAX( e2u(ji  ,jj)*e3u(ji  ,jj,jk,Kmm)*uu(ji  ,jj,jk,Kmm) + un_td(ji  ,jj,jk), 0._wp ) -   & 
    316311               &                                 MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm) + un_td(ji-1,jj,jk), 0._wp ) )   & 
     
    324319         DO_3D_00_00( 1, jpkm1 ) 
    325320            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    326             ! 2*rdt and not r2dt (for restartability) 
    327             Cu_adv(ji,jj,jk) = 2._wp * rdt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )   &  
     321            ! 2*rn_Dt and not rDt (for restartability) 
     322            Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )   &  
    328323               &                             + ( MAX( e2u(ji  ,jj)*e3u(ji  ,jj,jk,Kmm)*uu(ji  ,jj,jk,Kmm), 0._wp ) -   & 
    329324               &                                 MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm), 0._wp ) )   & 
  • NEMO/trunk/src/OCE/DYN/wet_dry.F90

    r12377 r12489  
    270270 
    271271 
    272    SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rdtbt ) 
     272   SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rDt_e ) 
    273273      !!---------------------------------------------------------------------- 
    274274      !!                  ***  ROUTINE wad_lmt  *** 
     
    280280      !! ** Action  : - calculate flux limiter and W/D flag 
    281281      !!---------------------------------------------------------------------- 
    282       REAL(wp)                , INTENT(in   ) ::   rdtbt    ! ocean time-step index 
     282      REAL(wp)                , INTENT(in   ) ::   rDt_e    ! ocean time-step index 
    283283      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zflxu,  zflxv, sshn_e, zssh_frc   
    284284      ! 
     
    299299      zdepwd = 50._wp   ! maximum depth that ocean cells can have W/D processes 
    300300      ! 
    301       z2dt = rdtbt    
     301      z2dt = rDt_e    
    302302      ! 
    303303      zflxp(:,:)   = 0._wp 
Note: See TracChangeset for help on using the changeset viewer.