New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 11480 for NEMO/branches/2019 – NEMO

Changeset 11480 for NEMO/branches/2019


Ignore:
Timestamp:
2019-08-29T11:23:25+02:00 (5 years ago)
Author:
davestorkey
Message:

2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Merge in changes from branch of branch.
Main changes:

  1. "nxt" modules renamed as "atf" and now just do Asselin time filtering. The time level swapping is achieved by swapping indices.
  2. Some additional prognostic grid variables changed to use a time dimension.

Notes:

  1. This merged branch passes SETTE tests but does not identical results to the SETTE tests with the trunk@10721 unless minor bugs to do with Euler timestepping and the OFF timestepping are fixed in the trunk (NEMO tickets #2310 and #2311).
  2. The nn_dttrc > 1 option for TOP (TOP has a different timestep to OCE) doesn't work. But it doesn't work in the trunk or NEMO 4.0 release either.
Location:
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps
Files:
3 deleted
38 edited
3 copied

Legend:

Unmodified
Added
Removed

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/ICE/icedyn.F90

    r10564 r11480  
    6060CONTAINS 
    6161 
    62    SUBROUTINE ice_dyn( kt ) 
     62   SUBROUTINE ice_dyn( kt, Kmm ) 
    6363      !!------------------------------------------------------------------- 
    6464      !!               ***  ROUTINE ice_dyn  *** 
    … …  
    7373      !!-------------------------------------------------------------------- 
    7474      INTEGER, INTENT(in) ::   kt     ! ice time step 
     75      INTEGER, INTENT(in) ::   Kmm    ! ocean time level index 
    7576      !! 
    7677      INTEGER  ::   ji, jj, jl        ! dummy loop indices 
    … …  
    9293         tau_icebfr(:,:) = 0._wp 
    9394         DO jl = 1, jpl 
    94             WHERE( h_i_b(:,:,jl) > ht_n(:,:) * rn_depfra )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr 
     95            WHERE( h_i_b(:,:,jl) > ht(:,:) * rn_depfra )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr 
    9596         END DO 
    9697      ENDIF 
    … …  
    121122 
    122123      CASE ( np_dynALL )           !==  all dynamical processes  ==! 
    123          CALL ice_dyn_rhg   ( kt )                                                 ! -- rheology   
     124         CALL ice_dyn_rhg   ( kt, Kmm )                                            ! -- rheology   
    124125         CALL ice_dyn_adv   ( kt )   ;   CALL Hbig( zhi_max, zhs_max, zhip_max )   ! -- advection of ice + correction on ice thickness 
    125126         CALL ice_dyn_rdgrft( kt )                                                 ! -- ridging/rafting  
    … …  
    127128 
    128129      CASE ( np_dynRHGADV  )       !==  no ridge/raft & no corrections ==! 
    129          CALL ice_dyn_rhg   ( kt )                                                 ! -- rheology   
     130         CALL ice_dyn_rhg   ( kt, Kmm )                                            ! -- rheology   
    130131         CALL ice_dyn_adv   ( kt )   ;   CALL Hbig( zhi_max, zhs_max, zhip_max )   ! -- advection of ice + correction on ice thickness 
    131132         CALL Hpiling                                                              ! -- simple pile-up (replaces ridging/rafting) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/ICE/icedyn_rhg.F90

    r10413 r11480  
    4747CONTAINS 
    4848 
    49    SUBROUTINE ice_dyn_rhg( kt ) 
     49   SUBROUTINE ice_dyn_rhg( kt, Kmm ) 
    5050      !!------------------------------------------------------------------- 
    5151      !!               ***  ROUTINE ice_dyn_rhg  *** 
    … …  
    5858      !!-------------------------------------------------------------------- 
    5959      INTEGER, INTENT(in) ::   kt     ! ice time step 
     60      INTEGER, INTENT(in) ::   Kmm    ! ocean time level index 
    6061      !!-------------------------------------------------------------------- 
    6162      ! controls 
    … …  
    7677      CASE( np_rhgEVP )                ! Elasto-Viscous-Plastic ! 
    7778         !                             !------------------------! 
    78          CALL ice_dyn_rhg_evp( kt, stress1_i, stress2_i, stress12_i, shear_i, divu_i, delta_i ) 
     79         CALL ice_dyn_rhg_evp( kt, Kmm, stress1_i, stress2_i, stress12_i, shear_i, divu_i, delta_i ) 
    7980         !          
    8081      END SELECT 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/ICE/icedyn_rhg_evp.F90

    r10555 r11480  
    5656CONTAINS 
    5757 
    58    SUBROUTINE ice_dyn_rhg_evp( kt, pstress1_i, pstress2_i, pstress12_i, pshear_i, pdivu_i, pdelta_i ) 
     58   SUBROUTINE ice_dyn_rhg_evp( kt, Kmm, pstress1_i, pstress2_i, pstress12_i, pshear_i, pdivu_i, pdelta_i ) 
    5959      !!------------------------------------------------------------------- 
    6060      !!                 ***  SUBROUTINE ice_dyn_rhg_evp  *** 
    … …  
    109109      !!------------------------------------------------------------------- 
    110110      INTEGER                 , INTENT(in   ) ::   kt                                    ! time step 
     111      INTEGER                 , INTENT(in   ) ::   Kmm                                   ! ocean time level index 
    111112      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   pstress1_i, pstress2_i, pstress12_i   ! 
    112113      REAL(wp), DIMENSION(:,:), INTENT(  out) ::   pshear_i  , pdivu_i   , pdelta_i      ! 
    … …  
    335336               zvV = 0.5_wp * ( vt_i(ji,jj) * e1e2t(ji,jj) + vt_i(ji,jj+1) * e1e2t(ji,jj+1) ) * r1_e1e2v(ji,jj) * vmask(ji,jj,1) 
    336337               ! ice-bottom stress at U points 
    337                zvCr = zaU(ji,jj) * rn_depfra * hu_n(ji,jj) 
     338               zvCr = zaU(ji,jj) * rn_depfra * hu(ji,jj,Kmm) 
    338339               zTauU_ib(ji,jj)   = rn_icebfr * MAX( 0._wp, zvU - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaU(ji,jj) ) ) 
    339340               ! ice-bottom stress at V points 
    340                zvCr = zaV(ji,jj) * rn_depfra * hv_n(ji,jj) 
     341               zvCr = zaV(ji,jj) * rn_depfra * hv(ji,jj,Kmm) 
    341342               zTauV_ib(ji,jj)   = rn_icebfr * MAX( 0._wp, zvV - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaV(ji,jj) ) ) 
    342343               ! ice_bottom stress at T points 
    343                zvCr = at_i(ji,jj) * rn_depfra * ht_n(ji,jj) 
     344               zvCr = at_i(ji,jj) * rn_depfra * ht(ji,jj) 
    344345               tau_icebfr(ji,jj) = rn_icebfr * MAX( 0._wp, vt_i(ji,jj) - zvCr ) * EXP( -rn_crhg * ( 1._wp - at_i(ji,jj) ) ) 
    345346            END DO 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/ICE/icestp.F90

    r10998 r11480  
    9595CONTAINS 
    9696 
    97    SUBROUTINE ice_stp( kt, Kbb, ksbc ) 
     97   SUBROUTINE ice_stp( kt, Kbb, Kmm, ksbc ) 
    9898      !!--------------------------------------------------------------------- 
    9999      !!                  ***  ROUTINE ice_stp  *** 
    … …  
    115115      !!                utau, vtau, taum, wndm, qns , qsr, emp , sfx 
    116116      !!--------------------------------------------------------------------- 
    117       INTEGER, INTENT(in) ::   kt     ! ocean time step 
    118       INTEGER, INTENT(in) ::   Kbb    ! ocean time level index 
    119       INTEGER, INTENT(in) ::   ksbc   ! flux formulation (user defined, bulk, or Pure Coupled) 
     117      INTEGER, INTENT(in) ::   kt       ! ocean time step 
     118      INTEGER, INTENT(in) ::   Kbb, Kmm ! ocean time level indices 
     119      INTEGER, INTENT(in) ::   ksbc     ! flux formulation (user defined, bulk, or Pure Coupled) 
    120120      ! 
    121121      INTEGER ::   jl   ! dummy loop index 
    … …  
    161161         ! 
    162162         IF( ln_icedyn .AND. .NOT.lk_c1d )   & 
    163             &                           CALL ice_dyn( kt )            ! -- Ice dynamics 
     163            &                           CALL ice_dyn( kt, Kmm )       ! -- Ice dynamics 
    164164         ! 
    165165         !                          !==  lateral boundary conditions  ==! 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST/agrif_oce_interp.F90

    r11428 r11480  
    116116            END DO 
    117117            DO jj = 1, jpj 
    118                uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu_a(ibdy1:ibdy2,jj) 
     118               uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    119119            END DO 
    120120         ENDIF 
    … …  
    137137         END DO 
    138138         DO jj=1,jpj 
    139             zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 
     139            zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    140140         END DO 
    141141             
    … …  
    158158            END DO 
    159159            DO jj = 1, jpj 
    160                zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj) 
     160               zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv(ibdy1:ibdy2,jj,Krhs_a) 
    161161            END DO 
    162162            DO jk = 1, jpkm1 
    … …  
    194194            END DO 
    195195            DO jj = 1, jpj 
    196                uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu_a(ibdy1:ibdy2,jj) 
     196               uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    197197            END DO 
    198198         ENDIF 
    … …  
    217217         END DO 
    218218         DO jj=1,jpj 
    219             zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 
     219            zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    220220         END DO 
    221221             
    … …  
    242242            END DO 
    243243            DO jj = 1, jpj 
    244                zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj) 
     244               zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv(ibdy1:ibdy2,jj,Krhs_a) 
    245245            END DO 
    246246            DO jk = 1, jpkm1 
    … …  
    278278            END DO 
    279279            DO ji=1,jpi 
    280                vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv_a(ji,jbdy1:jbdy2) 
     280               vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    281281            END DO 
    282282         ENDIF 
    … …  
    302302         END DO 
    303303         DO ji = 1, jpi 
    304             zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 
     304            zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    305305         END DO 
    306306 
    … …  
    325325            END DO 
    326326            DO ji = 1, jpi 
    327                zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2) 
     327               zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu(ji,jbdy1:jbdy2,Krhs_a) 
    328328            END DO 
    329329                
    … …  
    362362            END DO 
    363363            DO ji=1,jpi 
    364                vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv_a(ji,jbdy1:jbdy2) 
     364               vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    365365            END DO 
    366366         ENDIF 
    … …  
    386386         END DO 
    387387         DO ji = 1, jpi 
    388             zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 
     388            zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    389389         END DO 
    390390 
    … …  
    411411            END DO 
    412412            DO ji = 1, jpi 
    413                zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2) 
     413               zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu(ji,jbdy1:jbdy2,Krhs_a) 
    414414            END DO 
    415415                
    … …  
    11291129      ! 
    11301130      IF( before ) THEN  
    1131          ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * uu_b(i1:i2,j1:j2,Kmm_a) 
     1131         ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu(i1:i2,j1:j2,Kmm_a) * uu_b(i1:i2,j1:j2,Kmm_a) 
    11321132      ELSE 
    11331133         western_side  = (nb == 1).AND.(ndir == 1) 
    … …  
    11821182      !  
    11831183      IF( before ) THEN  
    1184          ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vv_b(i1:i2,j1:j2,Kmm_a) 
     1184         ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv(i1:i2,j1:j2,Kmm_a) * vv_b(i1:i2,j1:j2,Kmm_a) 
    11851185      ELSE 
    11861186         western_side  = (nb == 1).AND.(ndir == 1) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/NST/agrif_oce_update.F90

    r11428 r11480  
    234234!      uu(:,:,:,Krhs_a) = e3u(:,:,:,Kbb_a) 
    235235!      vv(:,:,:,Krhs_a) = e3v(:,:,:,Kbb_a) 
    236       hu_a(:,:) = hu_n(:,:) 
    237       hv_a(:,:) = hv_n(:,:) 
     236      hu(:,:,Krhs_a) = hu(:,:,Kmm_a) 
     237      hv(:,:,Krhs_a) = hv(:,:,Kmm_a) 
    238238 
    239239      ! 1) NOW fields 
    … …  
    251251         ! Update total depths: 
    252252         ! -------------------- 
    253       hu_n(:,:) = 0._wp                        ! Ocean depth at U-points 
    254       hv_n(:,:) = 0._wp                        ! Ocean depth at V-points 
     253      hu(:,:,Kmm_a) = 0._wp                        ! Ocean depth at U-points 
     254      hv(:,:,Kmm_a) = 0._wp                        ! Ocean depth at V-points 
    255255      DO jk = 1, jpkm1 
    256          hu_n(:,:) = hu_n(:,:) + e3u(:,:,jk,Kmm_a) * umask(:,:,jk) 
    257          hv_n(:,:) = hv_n(:,:) + e3v(:,:,jk,Kmm_a) * vmask(:,:,jk) 
     256         hu(:,:,Kmm_a) = hu(:,:,Kmm_a) + e3u(:,:,jk,Kmm_a) * umask(:,:,jk) 
     257         hv(:,:,Kmm_a) = hv(:,:,Kmm_a) + e3v(:,:,jk,Kmm_a) * vmask(:,:,jk) 
    258258      END DO 
    259259      !                                        ! Inverse of the local depth 
    260       r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
    261       r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     260      r1_hu(:,:,Kmm_a) = ssumask(:,:) / ( hu(:,:,Kmm_a) + 1._wp - ssumask(:,:) ) 
     261      r1_hv(:,:,Kmm_a) = ssvmask(:,:) / ( hv(:,:,Kmm_a) + 1._wp - ssvmask(:,:) ) 
    262262 
    263263 
    … …  
    276276         ! Update total depths: 
    277277         ! -------------------- 
    278          hu_b(:,:) = 0._wp                     ! Ocean depth at U-points 
    279          hv_b(:,:) = 0._wp                     ! Ocean depth at V-points 
     278         hu(:,:,Kbb_a) = 0._wp                     ! Ocean depth at U-points 
     279         hv(:,:,Kbb_a) = 0._wp                     ! Ocean depth at V-points 
    280280         DO jk = 1, jpkm1 
    281             hu_b(:,:) = hu_b(:,:) + e3u(:,:,jk,Kbb_a) * umask(:,:,jk) 
    282             hv_b(:,:) = hv_b(:,:) + e3v(:,:,jk,Kbb_a) * vmask(:,:,jk) 
     281            hu(:,:,Kbb_a) = hu(:,:,Kbb_a) + e3u(:,:,jk,Kbb_a) * umask(:,:,jk) 
     282            hv(:,:,Kbb_a) = hv(:,:,Kbb_a) + e3v(:,:,jk,Kbb_a) * vmask(:,:,jk) 
    283283         END DO 
    284284         !                                     ! Inverse of the local depth 
    285          r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) ) 
    286          r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
     285         r1_hu(:,:,Kbb_a) = ssumask(:,:) / ( hu(:,:,Kbb_a) + 1._wp - ssumask(:,:) ) 
     286         r1_hv(:,:,Kbb_a) = ssvmask(:,:) / ( hv(:,:,Kbb_a) + 1._wp - ssvmask(:,:) ) 
    287287      ENDIF 
    288288      ! 
    … …  
    636636         IF (western_side) THEN 
    637637            DO jj=j1,j2 
    638                zcor = uu_b(i1-1,jj,Kmm_a) * hu_a(i1-1,jj) * r1_hu_n(i1-1,jj) - uu_b(i1-1,jj,Kmm_a) 
     638               zcor = uu_b(i1-1,jj,Kmm_a) * hu(i1-1,jj,Krhs_a) * r1_hu(i1-1,jj,Kmm_a) - uu_b(i1-1,jj,Kmm_a) 
    639639               uu_b(i1-1,jj,Kmm_a) = uu_b(i1-1,jj,Kmm_a) + zcor 
    640640               DO jk=1,jpkm1 
    … …  
    646646         IF (eastern_side) THEN 
    647647            DO jj=j1,j2 
    648                zcor = uu_b(i2+1,jj,Kmm_a) * hu_a(i2+1,jj) * r1_hu_n(i2+1,jj) - uu_b(i2+1,jj,Kmm_a) 
     648               zcor = uu_b(i2+1,jj,Kmm_a) * hu(i2+1,jj,Krhs_a) * r1_hu(i2+1,jj,Kmm_a) - uu_b(i2+1,jj,Kmm_a) 
    649649               uu_b(i2+1,jj,Kmm_a) = uu_b(i2+1,jj,Kmm_a) + zcor 
    650650               DO jk=1,jpkm1 
    … …  
    829829         IF (southern_side) THEN 
    830830            DO ji=i1,i2 
    831                zcor = vv_b(ji,j1-1,Kmm_a) * hv_a(ji,j1-1) * r1_hv_n(ji,j1-1) - vv_b(ji,j1-1,Kmm_a) 
     831               zcor = vv_b(ji,j1-1,Kmm_a) * hv(ji,j1-1,Krhs_a) * r1_hv(ji,j1-1,Kmm_a) - vv_b(ji,j1-1,Kmm_a) 
    832832               vv_b(ji,j1-1,Kmm_a) = vv_b(ji,j1-1,Kmm_a) + zcor 
    833833               DO jk=1,jpkm1 
    … …  
    839839         IF (northern_side) THEN 
    840840            DO ji=i1,i2 
    841                zcor = vv_b(ji,j2+1,Kmm_a) * hv_a(ji,j2+1) * r1_hv_n(ji,j2+1) - vv_b(ji,j2+1,Kmm_a) 
     841               zcor = vv_b(ji,j2+1,Kmm_a) * hv(ji,j2+1,Krhs_a) * r1_hv(ji,j2+1,Kmm_a) - vv_b(ji,j2+1,Kmm_a) 
    842842               vv_b(ji,j2+1,Kmm_a) = vv_b(ji,j2+1,Kmm_a) + zcor 
    843843               DO jk=1,jpkm1 
    … …  
    869869         DO jj=j1,j2 
    870870            DO ji=i1,i2 
    871                tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu_n(ji,jj) * e2u(ji,jj) 
     871               tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu(ji,jj,Kmm_a) * e2u(ji,jj) 
    872872            END DO 
    873873         END DO 
    … …  
    883883               END DO 
    884884               ! 
    885                zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj) 
     885               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a) 
    886886               DO jk=1,jpkm1               
    887887                  uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk)            
    … …  
    891891               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 
    892892                  IF ( .NOT.( lk_agrif_fstep .AND. (neuler==0) ) ) THEN                          ! Add asselin part 
    893                      zcorr = ( tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu_a(ji,jj) ) * r1_hu_b(ji,jj) 
     893                     zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a) 
    894894                     uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + atfp * zcorr * umask(ji,jj,1) 
    895895                  END IF 
    896896               ENDIF     
    897                uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1) 
     897               uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1) 
    898898               !        
    899899               ! Correct "before" velocities to hold correct bt component: 
    … …  
    903903               END DO 
    904904               ! 
    905                zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu_b(ji,jj) 
     905               zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a) 
    906906               DO jk=1,jpkm1               
    907907                  uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk)            
    … …  
    935935         DO jj=j1,j2 
    936936            DO ji=i1,i2 
    937                tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv_n(ji,jj) * e1v(ji,jj)  
     937               tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv(ji,jj,Kmm_a) * e1v(ji,jj)  
    938938            END DO 
    939939         END DO 
    … …  
    949949               END DO 
    950950               ! 
    951                zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj) 
     951               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a) 
    952952               DO jk=1,jpkm1               
    953953                  vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk)            
    … …  
    957957               IF ( .NOT.ln_dynspg_ts .OR. ( ln_dynspg_ts .AND. ( .NOT.ln_bt_fw ) ) ) THEN 
    958958                  IF ( .NOT. ( lk_agrif_fstep .AND. ( neuler==0 ) ) ) THEN                       ! Add asselin part 
    959                      zcorr = ( tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv_a(ji,jj) ) * r1_hv_b(ji,jj) 
     959                     zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a) 
    960960                     vv_b(ji,jj,Kbb_a) =       vv_b(ji,jj,Kbb_a) + atfp * zcorr  *   vmask(ji,jj,1) 
    961961                  END IF 
    962962               ENDIF               
    963                vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1) 
     963               vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1) 
    964964               !        
    965965               ! Correct "before" velocities to hold correct bt component: 
    … …  
    969969               END DO 
    970970               ! 
    971                zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv_b(ji,jj) 
     971               zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a) 
    972972               DO jk=1,jpkm1               
    973973                  vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk)            
    … …  
    13871387         ! 
    13881388         ! Update total depth: 
    1389          ht_n(i1:i2,j1:j2) = 0._wp 
     1389         ht(i1:i2,j1:j2) = 0._wp 
    13901390         DO jk = 1, jpkm1 
    1391             ht_n(i1:i2,j1:j2) = ht_n(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm_a) * tmask(i1:i2,j1:j2,jk) 
     1391            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm_a) * tmask(i1:i2,j1:j2,jk) 
    13921392         END DO 
    13931393         ! 
    … …  
    13961396         gdept(i1:i2,j1:j2,1,Kmm_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kmm_a) 
    13971397         gdepw(i1:i2,j1:j2,1,Kmm_a) = 0.0_wp 
    1398          gde3w(i1:i2,j1:j2,1      ) = gdept(i1:i2,j1:j2,1,Kmm_a) - ( ht_n(i1:i2,j1:j2) - ht_0(i1:i2,j1:j2) ) ! Last term in the rhs is ssh 
     1398         gde3w(i1:i2,j1:j2,1      ) = gdept(i1:i2,j1:j2,1,Kmm_a) - ( ht(i1:i2,j1:j2) - ht_0(i1:i2,j1:j2) ) ! Last term in the rhs is ssh 
    13991399         ! 
    14001400         DO jk = 2, jpk 
    … …  
    14091409               gdept(ji,jj,jk,Kmm_a) =             zcoef   * ( gdepw(ji,jj,jk  ,Kmm_a) + 0.5 * e3w(ji,jj,jk  ,Kmm_a) )  & 
    14101410                 &                     + ( 1._wp - zcoef ) * ( gdept(ji,jj,jk-1,Kmm_a) +       e3w(ji,jj,jk  ,Kmm_a) )  
    1411                gde3w(ji,jj,jk      ) =                         gdept(ji,jj,jk  ,Kmm_a) - ( ht_n(ji,jj)-ht_0(ji,jj) )    ! Last term in the rhs is ssh 
     1411               gde3w(ji,jj,jk      ) =                         gdept(ji,jj,jk  ,Kmm_a) - ( ht(ji,jj)-ht_0(ji,jj) )    ! Last term in the rhs is ssh 
    14121412               END DO 
    14131413            END DO 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ASM/asminc.F90

    r10969 r11480  
    805805      ELSE  
    806806         ALLOCATE( ztim(jpi,jpj) ) 
    807          ztim(:,:) = ssh_iau(:,:) / ( ht_n(:,:) + 1.0 - ssmask(:,:) ) 
     807         ztim(:,:) = ssh_iau(:,:) / ( ht(:,:) + 1.0 - ssmask(:,:) ) 
    808808         DO jk = 1, jpkm1                                  
    809809            phdivn(:,:,jk) = phdivn(:,:,jk) - ztim(:,:) * tmask(:,:,jk)  

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdydta.F90

    r10957 r11480  
    255255                  CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 
    256256                       & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy,   & 
    257                        & fvl=ln_full_vel_array(jbdy) ) 
     257                       & fvl=ln_full_vel_array(jbdy), Kmm=Kmm ) 
    258258               ENDIF 
    259259               ! If full velocities in boundary data then split into barotropic and baroclinic data 
    … …  
    270270                             &                       + e3u(ii,ij,ik,Kmm) * umask(ii,ij,ik) * dta%u3d(ib,ik) 
    271271                     END DO 
    272                      dta%u2d(ib) =  dta%u2d(ib) * r1_hu_n(ii,ij) 
     272                     dta%u2d(ib) =  dta%u2d(ib) * r1_hu(ii,ij,Kmm) 
    273273                     DO ik = 1, jpkm1 
    274274                        dta%u3d(ib,ik) = dta%u3d(ib,ik) - dta%u2d(ib) 
    … …  
    284284                             &                       + e3v(ii,ij,ik,Kmm) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 
    285285                     END DO 
    286                      dta%v2d(ib) =  dta%v2d(ib) * r1_hv_n(ii,ij) 
     286                     dta%v2d(ib) =  dta%v2d(ib) * r1_hv(ii,ij,Kmm) 
    287287                     DO ik = 1, jpkm1 
    288288                        dta%v3d(ib,ik) = dta%v3d(ib,ik) - dta%v2d(ib) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdydyn.F90

    r10957 r11480  
    7878         zva2d(:,:) = zva2d(:,:) + e3v(:,:,jk,Kaa) * pvv(:,:,jk,Kaa) * vmask(:,:,jk) 
    7979      END DO 
    80       zua2d(:,:) = zua2d(:,:) * r1_hu_a(:,:) 
    81       zva2d(:,:) = zva2d(:,:) * r1_hv_a(:,:) 
     80      zua2d(:,:) = zua2d(:,:) * r1_hu(:,:,Kaa) 
     81      zva2d(:,:) = zva2d(:,:) * r1_hv(:,:,Kaa) 
    8282 
    8383      DO jk = 1 , jpkm1 
    … …  
    9999      !------------------------------------------------------- 
    100100 
    101       IF( ll_dyn2d )   CALL bdy_dyn2d( kt, zua2d, zva2d, uu_b(:,:,Kbb), vv_b(:,:,Kbb), r1_hu_a(:,:), r1_hv_a(:,:), ssh(:,:,Kaa) ) 
     101      IF( ll_dyn2d )   CALL bdy_dyn2d( kt, zua2d, zva2d, uu_b(:,:,Kbb), vv_b(:,:,Kbb), r1_hu(:,:,Kaa), r1_hv(:,:,Kaa), ssh(:,:,Kaa) ) 
    102102 
    103103      IF( ll_dyn3d )   CALL bdy_dyn3d( kt, Kbb, puu, pvv, Kaa ) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/C1D/step_c1d.F90

    r11001 r11480  
    107107      IF(.NOT.ln_linssh)CALL tra_adv( kstp, Nbb, Nnn, ts, Nrhs  )  ! horizontal & vertical advection 
    108108      IF( ln_zdfosm  )  CALL tra_osm( kstp, Nnn     , ts, Nrhs  )  ! OSMOSIS non-local tracer fluxes 
    109                         CALL tra_zdf( kstp, Nbb, Nnn, Nrhs, ts, Naa   ) ! vertical mixing 
     109                        CALL tra_zdf( kstp, Nbb, Nnn, Nrhs, ts, Naa   )         ! vertical mixing 
    110110                        CALL eos( ts(:,:,:,:,Nnn), rhd, rhop, gdept_0(:,:,:) )  ! now potential density for zdfmxl 
    111       IF( ln_zdfnpc )   CALL tra_npc( kstp,      Nnn, Nrhs, ts, Naa   ) ! applied non penetrative convective adjustment on (t,s) 
    112                         CALL tra_nxt( kstp, Nbb, Nnn, Nrhs,     Naa   ) ! tracer fields at next time step 
     111      IF( ln_zdfnpc )   CALL tra_npc( kstp,      Nnn, Nrhs, ts, Naa   )         ! applied non penetrative convective adjustment on (t,s) 
     112                        CALL tra_atf( kstp, Nbb, Nnn, Nrhs,     Naa, ts   )     ! time filtering of "now" tracer fields 
    113113 
    114114 
    … …  
    124124      IF( ln_zdfosm  )  CALL dyn_osm    ( kstp,      Nnn      , uu, vv, Nrhs )  ! OSMOSIS non-local velocity fluxes 
    125125                        CALL dyn_zdf    ( kstp, Nbb, Nnn, Nrhs, uu, vv, Naa  )  ! vertical diffusion 
    126                         CALL dyn_nxt    ( kstp, Nbb, Nnn, Naa )                 ! lateral velocity at next time step 
    127       IF(.NOT.ln_linssh)CALL ssh_swp    ( kstp, Nbb, Nnn, Naa )                 ! swap of sea surface height 
    128  
    129       IF(.NOT.ln_linssh)CALL dom_vvl_sf_swp( kstp, Nbb, Nnn, Naa )              ! swap of vertical scale factors 
     126                        CALL dyn_atf    ( kstp, Nbb, Nnn, Naa , uu, vv, e3t, e3u, e3v )  ! time filtering of "now" fields 
     127      IF(.NOT.ln_linssh)CALL ssh_atf    ( kstp, Nbb, Nnn, Naa , ssh )                    ! time filtering of "now" sea surface height 
     128      ! 
     129      ! Swap time levels 
     130      Nrhs = Nbb 
     131      Nbb = Nnn 
     132      Nnn = Naa 
     133      Naa = Nrhs 
     134      ! 
     135      IF(.NOT.ln_linssh)CALL dom_vvl_sf_swp( kstp, Nbb, Nnn, Naa )                       ! swap of vertical scale factors 
    130136 
    131137      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/dom_oce.F90

    r11036 r11480  
    1212   !!            3.7  ! 2015-11  (G. Madec) introduce surface and scale factor ratio 
    1313   !!             -   ! 2015-11  (G. Madec, A. Coward)  time varying zgr by default 
     14   !!            4.1  ! 2019-08  (A. Coward, D. Storkey) rename prognostic variables in preparation for new time scheme. 
    1415   !!---------------------------------------------------------------------- 
    1516 
    … …  
    121122   REAL(wp), PUBLIC, ALLOCATABLE, SAVE        , DIMENSION(:,:) ::   e1e2f , r1_e1e2f                !: associated metrics at f-point 
    122123   ! 
    123    REAL(wp), PUBLIC, ALLOCATABLE, SAVE        , DIMENSION(:,:) ::   ff_f  , ff_t                    !: Coriolis factor at f- & t-points  [1/s] 
     124   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ff_f  , ff_t                    !: Coriolis factor at f- & t-points  [1/s] 
    124125   !!---------------------------------------------------------------------- 
    125126   !! vertical coordinate and scale factors 
    … …  
    149150   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  gde3w   
    150151    
    151    !                                                      !  ref. ! before  !   now   !  after  ! 
    152    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_0            ,    ht_n             !: t-depth              [m] 
    153    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu_0  ,    hu_b ,    hu_n ,    hu_a   !: u-depth              [m] 
    154    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hv_0  ,    hv_b ,    hv_n ,    hv_a   !: v-depth              [m] 
    155    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::           r1_hu_b , r1_hu_n , r1_hu_a   !: inverse of u-depth [1/m] 
    156    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::           r1_hv_b , r1_hv_n , r1_hv_a   !: inverse of v-depth [1/m] 
     152   !                                                      !  reference heights of water column 
     153   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_0  !: t-depth              [m] 
     154   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu_0  !: u-depth              [m] 
     155   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hv_0  !: v-depth              [m] 
     156                                                          ! time-dependent heights of water column 
     157   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ht                     !: height of water column at T-points [m] 
     158   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hu, hv, r1_hu, r1_hv   !: height of water column [m] and reciprocal [1/m] 
    157159 
    158160   INTEGER, PUBLIC ::   nla10              !: deepest    W level Above  ~10m (nlb10 - 1) 
    … …  
    268270         &      e3t  (jpi,jpj,jpk,jpt) , e3u  (jpi,jpj,jpk,jpt) , e3v  (jpi,jpj,jpk,jpt) , e3f  (jpi,jpj,jpk) , e3w  (jpi,jpj,jpk,jpt) ,   &  
    269271         &      e3uw_0(jpi,jpj,jpk)     , e3vw_0(jpi,jpj,jpk)     ,         & 
    270          &      e3uw  (jpi,jpj,jpk,jpt) , e3vw  (jpi,jpj,jpk,jpt) ,     STAT=ierr(5) )                        
    271          ! 
    272       ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj) , hv_0(jpi,jpj) ,                                           & 
    273          &                      hu_b(jpi,jpj) , hv_b(jpi,jpj) , r1_hu_b(jpi,jpj) , r1_hv_b(jpi,jpj) ,     & 
    274          &      ht_n(jpi,jpj) , hu_n(jpi,jpj) , hv_n(jpi,jpj) , r1_hu_n(jpi,jpj) , r1_hv_n(jpi,jpj) ,     & 
    275          &                      hu_a(jpi,jpj) , hv_a(jpi,jpj) , r1_hu_a(jpi,jpj) , r1_hv_a(jpi,jpj) , STAT=ierr(6)  ) 
     272         &      e3uw  (jpi,jpj,jpk,jpt) , e3vw  (jpi,jpj,jpk,jpt) ,    STAT=ierr(5) )                        
     273         ! 
     274      ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj)    , hv_0(jpi,jpj)     ,                                             & 
     275         &      ht  (jpi,jpj) , hu(  jpi,jpj,jpt), hv(  jpi,jpj,jpt) , r1_hu(jpi,jpj,jpt) , r1_hv(jpi,jpj,jpt) ,   & 
     276         &                      STAT=ierr(6)  ) 
    276277         ! 
    277278         ! 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/domain.F90

    r10978 r11480  
    163163      ! 
    164164         !       before        !          now          !       after         ! 
    165             gdept(:,:,:,Kbb) = gdept_0  ;   gdept(:,:,:,Kmm) = gdept_0   !        ---          ! depth of grid-points 
    166             gdepw(:,:,:,Kbb) = gdepw_0  ;   gdepw(:,:,:,Kmm) = gdepw_0   !        ---          ! 
     165            gdept(:,:,:,Kbb) = gdept_0  ;   gdept(:,:,:,Kmm) = gdept_0   ;   gdept(:,:,:,Kaa) = gdept_0   ! depth of grid-points 
     166            gdepw(:,:,:,Kbb) = gdepw_0  ;   gdepw(:,:,:,Kmm) = gdepw_0   ;   gdepw(:,:,:,Kaa) = gdepw_0   ! 
    167167                                   gde3w = gde3w_0   !        ---          ! 
    168168         !                                                                   
    … …  
    171171              e3v(:,:,:,Kbb) =   e3v_0  ;     e3v(:,:,:,Kmm) =   e3v_0   ;   e3v(:,:,:,Kaa) =  e3v_0    ! 
    172172                                     e3f =   e3f_0   !        ---          ! 
    173               e3w(:,:,:,Kbb) =   e3w_0  ;     e3w(:,:,:,Kmm) =   e3w_0   !        ---          ! 
    174              e3uw(:,:,:,Kbb) =  e3uw_0  ;    e3uw(:,:,:,Kmm) =  e3uw_0   !        ---          ! 
    175              e3vw(:,:,:,Kbb) =  e3vw_0  ;    e3vw(:,:,:,Kmm) =  e3vw_0   !        ---          ! 
     173              e3w(:,:,:,Kbb) =   e3w_0  ;     e3w(:,:,:,Kmm) =   e3w_0   ;    e3w(:,:,:,Kaa) =   e3w_0   !  
     174             e3uw(:,:,:,Kbb) =  e3uw_0  ;    e3uw(:,:,:,Kmm) =  e3uw_0   ;   e3uw(:,:,:,Kaa) =  e3uw_0   !   
     175             e3vw(:,:,:,Kbb) =  e3vw_0  ;    e3vw(:,:,:,Kmm) =  e3vw_0   ;   e3vw(:,:,:,Kaa) =  e3vw_0   ! 
    176176         ! 
    177177         z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF 
    … …  
    179179         ! 
    180180         !        before       !          now          !       after         ! 
    181                                       ht_n =    ht_0   !                     ! water column thickness 
    182                hu_b =    hu_0  ;      hu_n =    hu_0   ;    hu_a =    hu_0   !  
    183                hv_b =    hv_0  ;      hv_n =    hv_0   ;    hv_a =    hv_0   ! 
    184             r1_hu_b = z1_hu_0  ;   r1_hu_n = z1_hu_0   ; r1_hu_a = z1_hu_0   ! inverse of water column thickness 
    185             r1_hv_b = z1_hv_0  ;   r1_hv_n = z1_hv_0   ; r1_hv_a = z1_hv_0   ! 
     181                                      ht =    ht_0   !                     ! water column thickness 
     182               hu(:,:,Kbb) =    hu_0  ;      hu(:,:,Kmm) =    hu_0   ;    hu(:,:,Kaa) =    hu_0   !  
     183               hv(:,:,Kbb) =    hv_0  ;      hv(:,:,Kmm) =    hv_0   ;    hv(:,:,Kaa) =    hv_0   ! 
     184            r1_hu(:,:,Kbb) = z1_hu_0  ;   r1_hu(:,:,Kmm) = z1_hu_0   ; r1_hu(:,:,Kaa) = z1_hu_0   ! inverse of water column thickness 
     185            r1_hv(:,:,Kbb) = z1_hv_0  ;   r1_hv(:,:,Kmm) = z1_hv_0   ; r1_hv(:,:,Kaa) = z1_hv_0   ! 
    186186         ! 
    187187         ! 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/domvvl.F90

    r10978 r11480  
    1313   !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
    1414   !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
    15    !!   dom_vvl_sf_swp   : Swap vertical scale factors and update the vertical grid 
     15   !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
    1616   !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
    1717   !!   dom_vvl_rst      : read/write restart file 
    … …  
    3737   PUBLIC  dom_vvl_init       ! called by domain.F90 
    3838   PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
    39    PUBLIC  dom_vvl_sf_swp     ! called by step.F90 
     39   PUBLIC  dom_vvl_sf_update  ! called by step.F90 
    4040   PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
    4141 
    … …  
    181181      ! 
    182182      !                    !==  thickness of the water column  !!   (ocean portion only) 
    183       ht_n(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
    184       hu_b(:,:) = e3u(:,:,1,Kbb) * umask(:,:,1) 
    185       hu_n(:,:) = e3u(:,:,1,Kmm) * umask(:,:,1) 
    186       hv_b(:,:) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
    187       hv_n(:,:) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
     183      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
     184      hu(:,:,Kbb) = e3u(:,:,1,Kbb) * umask(:,:,1) 
     185      hu(:,:,Kmm) = e3u(:,:,1,Kmm) * umask(:,:,1) 
     186      hv(:,:,Kbb) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
     187      hv(:,:,Kmm) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
    188188      DO jk = 2, jpkm1 
    189          ht_n(:,:) = ht_n(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    190          hu_b(:,:) = hu_b(:,:) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
    191          hu_n(:,:) = hu_n(:,:) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
    192          hv_b(:,:) = hv_b(:,:) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
    193          hv_n(:,:) = hv_n(:,:) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
     189         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     190         hu(:,:,Kbb) = hu(:,:,Kbb) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
     191         hu(:,:,Kmm) = hu(:,:,Kmm) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
     192         hv(:,:,Kbb) = hv(:,:,Kbb) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
     193         hv(:,:,Kmm) = hv(:,:,Kmm) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    194194      END DO 
    195195      ! 
    196196      !                    !==  inverse of water column thickness   ==!   (u- and v- points) 
    197       r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
    198       r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
    199       r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
    200       r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     197      r1_hu(:,:,Kbb) = ssumask(:,:) / ( hu(:,:,Kbb) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
     198      r1_hu(:,:,Kmm) = ssumask(:,:) / ( hu(:,:,Kmm) + 1._wp - ssumask(:,:) ) 
     199      r1_hv(:,:,Kbb) = ssvmask(:,:) / ( hv(:,:,Kbb) + 1._wp - ssvmask(:,:) ) 
     200      r1_hv(:,:,Kmm) = ssvmask(:,:) / ( hv(:,:,Kmm) + 1._wp - ssvmask(:,:) ) 
    201201 
    202202      !                    !==   z_tilde coordinate case  ==!   (Restoring frequencies) 
    … …  
    550550      ! *********************************** ! 
    551551 
    552       hu_a(:,:) = e3u(:,:,1,Kaa) * umask(:,:,1) 
    553       hv_a(:,:) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
     552      hu(:,:,Kaa) = e3u(:,:,1,Kaa) * umask(:,:,1) 
     553      hv(:,:,Kaa) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
    554554      DO jk = 2, jpkm1 
    555          hu_a(:,:) = hu_a(:,:) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
    556          hv_a(:,:) = hv_a(:,:) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
     555         hu(:,:,Kaa) = hu(:,:,Kaa) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     556         hv(:,:,Kaa) = hv(:,:,Kaa) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    557557      END DO 
    558558      !                                        ! Inverse of the local depth 
    559559!!gm BUG ?  don't understand the use of umask_i here ..... 
    560       r1_hu_a(:,:) = ssumask(:,:) / ( hu_a(:,:) + 1._wp - ssumask(:,:) ) 
    561       r1_hv_a(:,:) = ssvmask(:,:) / ( hv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     560      r1_hu(:,:,Kaa) = ssumask(:,:) / ( hu(:,:,Kaa) + 1._wp - ssumask(:,:) ) 
     561      r1_hv(:,:,Kaa) = ssvmask(:,:) / ( hv(:,:,Kaa) + 1._wp - ssvmask(:,:) ) 
    562562      ! 
    563563      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_nxt') 
    … …  
    566566 
    567567 
    568    SUBROUTINE dom_vvl_sf_swp( kt, Kbb, Kmm, Kaa ) 
    569       !!---------------------------------------------------------------------- 
    570       !!                ***  ROUTINE dom_vvl_sf_swp  *** 
     568   SUBROUTINE dom_vvl_sf_update( kt, Kbb, Kmm, Kaa ) 
     569      !!---------------------------------------------------------------------- 
     570      !!                ***  ROUTINE dom_vvl_sf_update  *** 
    571571      !!                    
    572       !! ** Purpose :  compute time filter and swap of scale factors  
     572      !! ** Purpose :  for z tilde case: compute time filter and swap of scale factors  
    573573      !!               compute all depths and related variables for next time step 
    574574      !!               write outputs and restart file 
    575575      !! 
    576       !! ** Method  :  - swap of e3t with trick for volume/tracer conservation 
     576      !! ** Method  :  - swap of e3t with trick for volume/tracer conservation (ONLY FOR Z TILDE CASE) 
    577577      !!               - reconstruct scale factor at other grid points (interpolate) 
    578578      !!               - recompute depths and water height fields 
    579579      !! 
    580       !! ** Action  :  - e3t_(b/n), tilde_e3t_(b/n) and e3(u/v)_n ready for next time step 
     580      !! ** Action  :  - tilde_e3t_(b/n) ready for next time step 
    581581      !!               - Recompute: 
    582582      !!                    e3(u/v)_b        
    … …  
    599599      IF( ln_linssh )   RETURN      ! No calculation in linear free surface 
    600600      ! 
    601       IF( ln_timing )   CALL timing_start('dom_vvl_sf_swp') 
     601      IF( ln_timing )   CALL timing_start('dom_vvl_sf_update') 
    602602      ! 
    603603      IF( kt == nit000 )   THEN 
    604604         IF(lwp) WRITE(numout,*) 
    605          IF(lwp) WRITE(numout,*) 'dom_vvl_sf_swp : - time filter and swap of scale factors' 
    606          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~   - interpolate scale factors and compute depths for next time step' 
     605         IF(lwp) WRITE(numout,*) 'dom_vvl_sf_update : - interpolate scale factors and compute depths for next time step' 
     606         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~' 
    607607      ENDIF 
    608608      ! 
    … …  
    619619         tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) 
    620620      ENDIF 
    621       gdept(:,:,:,Kbb) = gdept(:,:,:,Kmm) 
    622       gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm) 
    623  
    624       e3t(:,:,:,Kmm) = e3t(:,:,:,Kaa) 
    625       e3u(:,:,:,Kmm) = e3u(:,:,:,Kaa) 
    626       e3v(:,:,:,Kmm) = e3v(:,:,:,Kaa) 
    627621 
    628622      ! Compute all missing vertical scale factor and depths 
    … …  
    630624      ! Horizontal scale factor interpolations 
    631625      ! -------------------------------------- 
    632       ! - ML - e3u(:,:,:,Kbb) and e3v(:,:,:,Kbb) are allready computed in dynnxt 
    633       ! - JC - hu_b, hv_b, hur_b, hvr_b also 
     626      ! - ML - e3u(:,:,:,Kbb) and e3v(:,:,:,Kbb) are already computed in dynnxt 
     627      ! - JC - hu(:,:,:,Kbb), hv(:,:,:,:,Kbb), hur_b, hvr_b also 
    634628       
    635629      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F'  ) 
    … …  
    663657      ! Local depth and Inverse of the local depth of the water 
    664658      ! ------------------------------------------------------- 
    665       hu_n(:,:) = hu_a(:,:)   ;   r1_hu_n(:,:) = r1_hu_a(:,:) 
    666       hv_n(:,:) = hv_a(:,:)   ;   r1_hv_n(:,:) = r1_hv_a(:,:) 
    667       ! 
    668       ht_n(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
     659      ! 
     660      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
    669661      DO jk = 2, jpkm1 
    670          ht_n(:,:) = ht_n(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     662         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    671663      END DO 
    672664 
    … …  
    675667      IF( lrst_oce  )   CALL dom_vvl_rst( kt, Kbb, Kmm, 'WRITE' ) 
    676668      ! 
    677       IF( ln_timing )   CALL timing_stop('dom_vvl_sf_swp') 
    678       ! 
    679    END SUBROUTINE dom_vvl_sf_swp 
     669      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_update') 
     670      ! 
     671   END SUBROUTINE dom_vvl_sf_update 
    680672 
    681673 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/iscplrst.F90

    r10978 r11480  
    108108      gdept(:,:,:,Kbb) = gdept(:,:,:,Kmm) 
    109109      gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm) 
    110       hu_b   (:,:)   = hu_n   (:,:) 
    111       hv_b   (:,:)   = hv_n   (:,:) 
    112       r1_hu_b(:,:)   = r1_hu_n(:,:) 
    113       r1_hv_b(:,:)   = r1_hv_n(:,:) 
     110      hu   (:,:,Kbb)   = hu   (:,:,Kmm) 
     111      hv   (:,:,Kbb)   = hv   (:,:,Kmm) 
     112      r1_hu(:,:,Kbb)   = r1_hu(:,:,Kmm) 
     113      r1_hv(:,:,Kbb)   = r1_hv(:,:,Kmm) 
    114114      ! 
    115115   END SUBROUTINE iscpl_stp 
    … …  
    240240      ! t-, u- and v- water column thickness 
    241241      ! ------------------------------------ 
    242          ht_n(:,:) = 0._wp ; hu_n(:,:) = 0._wp ; hv_n(:,:) = 0._wp 
     242         ht(:,:) = 0._wp ; hu(:,:,Kmm) = 0._wp ; hv(:,:,Kmm) = 0._wp 
    243243         DO jk = 1, jpkm1 
    244             hu_n(:,:) = hu_n(:,:) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
    245             hv_n(:,:) = hv_n(:,:) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    246             ht_n(:,:) = ht_n(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     244            hu(:,:,Kmm) = hu(:,:,Kmm) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
     245            hv(:,:,Kmm) = hv(:,:,Kmm) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
     246            ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    247247         END DO 
    248248         !                                        ! Inverse of the local depth 
    249          r1_hu_n(:,:) = 1._wp / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) * ssumask(:,:) 
    250          r1_hv_n(:,:) = 1._wp / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) * ssvmask(:,:) 
     249         r1_hu(:,:,Kmm) = 1._wp / ( hu(:,:,Kmm) + 1._wp - ssumask(:,:) ) * ssumask(:,:) 
     250         r1_hv(:,:,Kmm) = 1._wp / ( hv(:,:,Kmm) + 1._wp - ssvmask(:,:) ) * ssvmask(:,:) 
    251251 
    252252      END IF 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/istate.F90

    r10978 r11480  
    175175      END DO 
    176176      ! 
    177       uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * r1_hu_n(:,:) 
    178       vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * r1_hv_n(:,:) 
     177      uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * r1_hu(:,:,Kmm) 
     178      vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * r1_hv(:,:,Kmm) 
    179179      ! 
    180       uu_b(:,:,Kbb) = uu_b(:,:,Kbb) * r1_hu_b(:,:) 
    181       vv_b(:,:,Kbb) = vv_b(:,:,Kbb) * r1_hv_b(:,:) 
     180      uu_b(:,:,Kbb) = uu_b(:,:,Kbb) * r1_hu(:,:,Kbb) 
     181      vv_b(:,:,Kbb) = vv_b(:,:,Kbb) * r1_hv(:,:,Kbb) 
    182182      ! 
    183183   END SUBROUTINE istate_init 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynspg_ts.F90

    r11027 r11480  
    250250               DO jj = 1, jpjm1 
    251251                  DO ji = 1, jpim1 
    252                      zwz(ji,jj) =   ( ht_n(ji  ,jj+1) + ht_n(ji+1,jj+1) +                    & 
    253                         &             ht_n(ji  ,jj  ) + ht_n(ji+1,jj  )   ) * 0.25_wp   
     252                     zwz(ji,jj) =   ( ht(ji  ,jj+1) + ht(ji+1,jj+1) +                    & 
     253                        &             ht(ji  ,jj  ) + ht(ji+1,jj  )   ) * 0.25_wp   
    254254                     IF( zwz(ji,jj) /= 0._wp )   zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 
    255255                  END DO 
    … …  
    258258               DO jj = 1, jpjm1 
    259259                  DO ji = 1, jpim1 
    260                      zwz(ji,jj) =             (  ht_n  (ji  ,jj+1) + ht_n  (ji+1,jj+1)      & 
    261                         &                      + ht_n  (ji  ,jj  ) + ht_n  (ji+1,jj  )  )   & 
     260                     zwz(ji,jj) =             (  ht  (ji  ,jj+1) + ht  (ji+1,jj+1)      & 
     261                        &                      + ht  (ji  ,jj  ) + ht  (ji+1,jj  )  )   & 
    262262                        &       / ( MAX( 1._wp,  ssmask(ji  ,jj+1) + ssmask(ji+1,jj+1)      & 
    263263                        &                      + ssmask(ji  ,jj  ) + ssmask(ji+1,jj  )  )   ) 
    … …  
    282282            DO jj = 2, jpj 
    283283               DO ji = 2, jpi 
    284                   z1_ht = ssmask(ji,jj) / ( ht_n(ji,jj) + 1._wp - ssmask(ji,jj) ) 
     284                  z1_ht = ssmask(ji,jj) / ( ht(ji,jj) + 1._wp - ssmask(ji,jj) ) 
    285285                  ftne(ji,jj) = ( ff_f(ji-1,jj  ) + ff_f(ji  ,jj  ) + ff_f(ji  ,jj-1) ) * z1_ht 
    286286                  ftnw(ji,jj) = ( ff_f(ji-1,jj-1) + ff_f(ji-1,jj  ) + ff_f(ji  ,jj  ) ) * z1_ht 
    … …  
    367367      END DO 
    368368      ! 
    369       zu_frc(:,:) = zu_frc(:,:) * r1_hu_n(:,:) 
    370       zv_frc(:,:) = zv_frc(:,:) * r1_hv_n(:,:) 
     369      zu_frc(:,:) = zu_frc(:,:) * r1_hu(:,:,Kmm) 
     370      zv_frc(:,:) = zv_frc(:,:) * r1_hv(:,:,Kmm) 
    371371      ! 
    372372      ! 
    … …  
    388388      !                                   ! -------------------------------------------------------- 
    389389      ! 
    390       zwx(:,:) = puu_b(:,:,Kmm) * hu_n(:,:) * e2u(:,:)        ! now fluxes  
    391       zwy(:,:) = pvv_b(:,:,Kmm) * hv_n(:,:) * e1v(:,:) 
     390      zwx(:,:) = puu_b(:,:,Kmm) * hu(:,:,Kmm) * e2u(:,:)        ! now fluxes  
     391      zwy(:,:) = pvv_b(:,:,Kmm) * hv(:,:,Kmm) * e1v(:,:) 
    392392      ! 
    393393      SELECT CASE( nvor_scheme ) 
    … …  
    395395         DO jj = 2, jpjm1 
    396396            DO ji = 2, jpim1   ! vector opt. 
    397                zu_trd(ji,jj) = + r1_4 * r1_e1e2u(ji,jj) * r1_hu_n(ji,jj)                    & 
    398                   &               * (  e1e2t(ji+1,jj)*ht_n(ji+1,jj)*ff_t(ji+1,jj) * ( pvv_b(ji+1,jj,Kmm) + pvv_b(ji+1,jj-1,Kmm) )   & 
    399                   &                  + e1e2t(ji  ,jj)*ht_n(ji  ,jj)*ff_t(ji  ,jj) * ( pvv_b(ji  ,jj,Kmm) + pvv_b(ji  ,jj-1,Kmm) )   ) 
     397               zu_trd(ji,jj) = + r1_4 * r1_e1e2u(ji,jj) * r1_hu(ji,jj,Kmm)                    & 
     398                  &               * (  e1e2t(ji+1,jj)*ht(ji+1,jj)*ff_t(ji+1,jj) * ( pvv_b(ji+1,jj,Kmm) + pvv_b(ji+1,jj-1,Kmm) )   & 
     399                  &                  + e1e2t(ji  ,jj)*ht(ji  ,jj)*ff_t(ji  ,jj) * ( pvv_b(ji  ,jj,Kmm) + pvv_b(ji  ,jj-1,Kmm) )   ) 
    400400                  ! 
    401                zv_trd(ji,jj) = - r1_4 * r1_e1e2v(ji,jj) * r1_hv_n(ji,jj)                    & 
    402                   &               * (  e1e2t(ji,jj+1)*ht_n(ji,jj+1)*ff_t(ji,jj+1) * ( puu_b(ji,jj+1,Kmm) + puu_b(ji-1,jj+1,Kmm) )   &  
    403                   &                  + e1e2t(ji,jj  )*ht_n(ji,jj  )*ff_t(ji,jj  ) * ( puu_b(ji,jj  ,Kmm) + puu_b(ji-1,jj  ,Kmm) )   )  
     401               zv_trd(ji,jj) = - r1_4 * r1_e1e2v(ji,jj) * r1_hv(ji,jj,Kmm)                    & 
     402                  &               * (  e1e2t(ji,jj+1)*ht(ji,jj+1)*ff_t(ji,jj+1) * ( puu_b(ji,jj+1,Kmm) + puu_b(ji-1,jj+1,Kmm) )   &  
     403                  &                  + e1e2t(ji,jj  )*ht(ji,jj  )*ff_t(ji,jj  ) * ( puu_b(ji,jj  ,Kmm) + puu_b(ji-1,jj  ,Kmm) )   )  
    404404            END DO   
    405405         END DO   
    … …  
    546546            DO ji = fs_2, fs_jpim1   ! vector opt. 
    547547               zu_frc(ji,jj) = zu_frc(ji,jj) + &  
    548                & MAX(r1_hu_n(ji,jj) * r1_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ), zztmp ) * zwx(ji,jj) *  wdrampu(ji,jj) 
     548               & MAX(r1_hu(ji,jj,Kmm) * r1_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ), zztmp ) * zwx(ji,jj) *  wdrampu(ji,jj) 
    549549               zv_frc(ji,jj) = zv_frc(ji,jj) + &  
    550                & MAX(r1_hv_n(ji,jj) * r1_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ), zztmp ) * zwy(ji,jj) *  wdrampv(ji,jj) 
     550               & MAX(r1_hv(ji,jj,Kmm) * r1_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ), zztmp ) * zwy(ji,jj) *  wdrampv(ji,jj) 
    551551            END DO 
    552552         END DO 
    … …  
    554554         DO jj = 2, jpjm1 
    555555            DO ji = fs_2, fs_jpim1   ! vector opt. 
    556                zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * r1_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 
    557                zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * r1_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 
     556               zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu(ji,jj,Kmm) * r1_2 * ( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zwx(ji,jj) 
     557               zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv(ji,jj,Kmm) * r1_2 * ( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zwy(ji,jj) 
    558558            END DO 
    559559         END DO 
    … …  
    584584         DO jj = 2, jpjm1               
    585585            DO ji = fs_2, fs_jpim1   ! vector opt. 
    586                zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu_n(ji,jj) * r1_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zwx(ji,jj) 
    587                zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv_n(ji,jj) * r1_2 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * zwy(ji,jj) 
     586               zu_frc(ji,jj) = zu_frc(ji,jj) + r1_hu(ji,jj,Kmm) * r1_2 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zwx(ji,jj) 
     587               zv_frc(ji,jj) = zv_frc(ji,jj) + r1_hv(ji,jj,Kmm) * r1_2 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * zwy(ji,jj) 
    588588            END DO 
    589589         END DO 
    … …  
    593593         DO jj = 2, jpjm1 
    594594            DO ji = fs_2, fs_jpim1   ! vector opt. 
    595                zu_frc(ji,jj) =  zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 
    596                zv_frc(ji,jj) =  zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) 
     595               zu_frc(ji,jj) =  zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) 
     596               zv_frc(ji,jj) =  zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) 
    597597            END DO 
    598598         END DO 
    … …  
    601601         DO jj = 2, jpjm1 
    602602            DO ji = fs_2, fs_jpim1   ! vector opt. 
    603                zu_frc(ji,jj) =  zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 
    604                zv_frc(ji,jj) =  zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 
     603               zu_frc(ji,jj) =  zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) 
     604               zv_frc(ji,jj) =  zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv(ji,jj,Kmm) 
    605605            END DO 
    606606         END DO 
    … …  
    681681         vn_e  (:,:) =    pvv_b(:,:,Kmm) 
    682682         ! 
    683          hu_e  (:,:) =    hu_n(:,:)        
    684          hv_e  (:,:) =    hv_n(:,:)  
    685          hur_e (:,:) = r1_hu_n(:,:)     
    686          hvr_e (:,:) = r1_hv_n(:,:) 
     683         hu_e  (:,:) =    hu(:,:,Kmm)        
     684         hv_e  (:,:) =    hv(:,:,Kmm)  
     685         hur_e (:,:) = r1_hu(:,:,Kmm)     
     686         hvr_e (:,:) = r1_hv(:,:,Kmm) 
    687687      ELSE                                ! CENTRED integration: start from BEFORE fields 
    688688         sshn_e(:,:) =    pssh(:,:,Kbb) 
    … …  
    690690         vn_e  (:,:) =    pvv_b(:,:,Kbb) 
    691691         ! 
    692          hu_e  (:,:) =    hu_b(:,:)        
    693          hv_e  (:,:) =    hv_b(:,:)  
    694          hur_e (:,:) = r1_hu_b(:,:)     
    695          hvr_e (:,:) = r1_hv_b(:,:) 
     692         hu_e  (:,:) =    hu(:,:,Kbb)        
     693         hv_e  (:,:) =    hv(:,:,Kbb)  
     694         hur_e (:,:) = r1_hu(:,:,Kbb)     
     695         hvr_e (:,:) = r1_hv(:,:,Kbb) 
    696696      ENDIF 
    697697      ! 
    … …  
    790790            zhtp2_e(:,:) = ht_0(:,:) + zsshp2_e(:,:) 
    791791         ELSE 
    792             zhup2_e(:,:) = hu_n(:,:) 
    793             zhvp2_e(:,:) = hv_n(:,:) 
    794             zhtp2_e(:,:) = ht_n(:,:) 
     792            zhup2_e(:,:) = hu(:,:,Kmm) 
     793            zhvp2_e(:,:) = hv(:,:,Kmm) 
     794            zhtp2_e(:,:) = ht(:,:) 
    795795         ENDIF 
    796796         !                                                !* after ssh 
    … …  
    11381138                            &     + rdtbt * ( zhust_e(ji,jj)  *    zwx(ji,jj)   &  
    11391139                            &               + zhup2_e(ji,jj)  * zu_trd(ji,jj)   & 
    1140                             &               +    hu_n(ji,jj)  * zu_frc(ji,jj) ) & 
     1140                            &               +    hu(ji,jj,Kmm)  * zu_frc(ji,jj) ) & 
    11411141                            &   ) * zhura 
    11421142 
    … …  
    11441144                            &     + rdtbt * ( zhvst_e(ji,jj)  *    zwy(ji,jj)   & 
    11451145                            &               + zhvp2_e(ji,jj)  * zv_trd(ji,jj)   & 
    1146                             &               +    hv_n(ji,jj)  * zv_frc(ji,jj) ) & 
     1146                            &               +    hv(ji,jj,Kmm)  * zv_frc(ji,jj) ) & 
    11471147                            &   ) * zhvra 
    11481148               END DO 
    … …  
    12571257         ! 
    12581258         DO jk=1,jpkm1 
    1259             puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu_n(:,:) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu_b(:,:) ) * r1_2dt_b 
    1260             pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv_n(:,:) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv_b(:,:) ) * r1_2dt_b 
     1259            puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + r1_hu(:,:,Kmm) * ( puu_b(:,:,Kaa) - puu_b(:,:,Kbb) * hu(:,:,Kbb) ) * r1_2dt_b 
     1260            pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + r1_hv(:,:,Kmm) * ( pvv_b(:,:,Kaa) - pvv_b(:,:,Kbb) * hv(:,:,Kbb) ) * r1_2dt_b 
    12611261         END DO 
    12621262         ! Save barotropic velocities not transport: 
    … …  
    12681268      ! Correct velocities so that the barotropic velocity equals (un_adv, vn_adv) (in all cases)   
    12691269      DO jk = 1, jpkm1 
    1270          puu(:,:,jk,Kmm) = ( puu(:,:,jk,Kmm) + un_adv(:,:)*r1_hu_n(:,:) - puu_b(:,:,Kmm) ) * umask(:,:,jk) 
    1271          pvv(:,:,jk,Kmm) = ( pvv(:,:,jk,Kmm) + vn_adv(:,:)*r1_hv_n(:,:) - pvv_b(:,:,Kmm) ) * vmask(:,:,jk) 
     1270         puu(:,:,jk,Kmm) = ( puu(:,:,jk,Kmm) + un_adv(:,:)*r1_hu(:,:,Kmm) - puu_b(:,:,Kmm) ) * umask(:,:,jk) 
     1271         pvv(:,:,jk,Kmm) = ( pvv(:,:,jk,Kmm) + vn_adv(:,:)*r1_hv(:,:,Kmm) - pvv_b(:,:,Kmm) ) * vmask(:,:,jk) 
    12721272      END DO 
    12731273 
    12741274      IF ( ln_wd_dl .and. ln_wd_dl_bc) THEN  
    12751275         DO jk = 1, jpkm1 
    1276             puu(:,:,jk,Kmm) = ( un_adv(:,:)*r1_hu_n(:,:) & 
    1277                        & + zuwdav2(:,:)*(puu(:,:,jk,Kmm) - un_adv(:,:)*r1_hu_n(:,:)) ) * umask(:,:,jk)  
    1278             pvv(:,:,jk,Kmm) = ( vn_adv(:,:)*r1_hv_n(:,:) &  
    1279                        & + zvwdav2(:,:)*(pvv(:,:,jk,Kmm) - vn_adv(:,:)*r1_hv_n(:,:)) ) * vmask(:,:,jk)   
     1276            puu(:,:,jk,Kmm) = ( un_adv(:,:)*r1_hu(:,:,Kmm) & 
     1277                       & + zuwdav2(:,:)*(puu(:,:,jk,Kmm) - un_adv(:,:)*r1_hu(:,:,Kmm)) ) * umask(:,:,jk)  
     1278            pvv(:,:,jk,Kmm) = ( vn_adv(:,:)*r1_hv(:,:,Kmm) &  
     1279                       & + zvwdav2(:,:)*(pvv(:,:,jk,Kmm) - vn_adv(:,:)*r1_hv(:,:,Kmm)) ) * vmask(:,:,jk)   
    12801280         END DO 
    12811281      END IF  
    12821282 
    12831283       
    1284       CALL iom_put(  "ubar", un_adv(:,:)*r1_hu_n(:,:) )    ! barotropic i-current 
    1285       CALL iom_put(  "vbar", vn_adv(:,:)*r1_hv_n(:,:) )    ! barotropic i-current 
     1284      CALL iom_put(  "ubar", un_adv(:,:)*r1_hu(:,:,Kmm) )    ! barotropic i-current 
     1285      CALL iom_put(  "vbar", vn_adv(:,:)*r1_hv(:,:,Kmm) )    ! barotropic i-current 
    12861286      ! 
    12871287#if defined key_agrif 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/sshwzv.F90

    r11421 r11480  
    99   !!             -   !  2010-09  (D.Storkey and E.O'Dea) bug fixes for BDY module 
    1010   !!            3.3  !  2011-10  (M. Leclair) split former ssh_wzv routine and remove all vvl related work 
     11   !!            4.1  !  2019-08  (A. Coward, D. Storkey) Rename ssh_nxt -> ssh_atf. Now only does time filtering. 
    1112   !!---------------------------------------------------------------------- 
    1213 
    1314   !!---------------------------------------------------------------------- 
    1415   !!   ssh_nxt       : after ssh 
    15    !!   ssh_swp       : filter ans swap the ssh arrays 
     16   !!   ssh_atf       : time filter the ssh arrays 
    1617   !!   wzv           : compute now vertical velocity 
    1718   !!---------------------------------------------------------------------- 
    … …  
    4344   PUBLIC   wzv        ! called by step.F90 
    4445   PUBLIC   wAimp      ! called by step.F90 
    45    PUBLIC   ssh_swp    ! called by step.F90 
     46   PUBLIC   ssh_atf    ! called by step.F90 
    4647 
    4748   !! * Substitutions 
    … …  
    218219 
    219220 
    220    SUBROUTINE ssh_swp( kt, Kbb, Kmm, Kaa ) 
    221       !!---------------------------------------------------------------------- 
    222       !!                    ***  ROUTINE ssh_nxt  *** 
    223       !! 
    224       !! ** Purpose :   achieve the sea surface  height time stepping by  
    225       !!              applying Asselin time filter and swapping the arrays 
    226       !!              ssh(:,:,Kaa)  already computed in ssh_nxt   
     221   SUBROUTINE ssh_atf( kt, Kbb, Kmm, Kaa, pssh ) 
     222      !!---------------------------------------------------------------------- 
     223      !!                    ***  ROUTINE ssh_atf  *** 
     224      !! 
     225      !! ** Purpose :   Apply Asselin time filter to now SSH. 
    227226      !! 
    228227      !! ** Method  : - apply Asselin time fiter to now ssh (excluding the forcing 
    229228      !!              from the filter, see Leclair and Madec 2010) and swap : 
    230       !!                ssh(:,:,Kmm) = ssh(:,:,Kaa) + atfp * ( ssh(:,:,Kbb) -2 ssh(:,:,Kmm) + ssh(:,:,Kaa) ) 
     229      !!                pssh(:,:,Kmm) = pssh(:,:,Kaa) + atfp * ( pssh(:,:,Kbb) -2 pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
    231230      !!                            - atfp * rdt * ( emp_b - emp ) / rau0 
    232       !!                ssh(:,:,Kmm) = ssh(:,:,Kaa) 
    233       !! 
    234       !! ** action  : - ssh(:,:,Kbb), ssh(:,:,Kmm)   : before & now sea surface height 
    235       !!                               ready for the next time step 
     231      !! 
     232      !! ** action  : - pssh(:,:,Kmm) time filtered 
    236233      !! 
    237234      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling. 
    238235      !!---------------------------------------------------------------------- 
    239       INTEGER, INTENT(in) ::   kt              ! ocean time-step index 
    240       INTEGER, INTENT(in) ::   Kbb, Kmm, Kaa   ! ocean time-step index 
     236      INTEGER                         , INTENT(in   ) ::   kt             ! ocean time-step index 
     237      INTEGER                         , INTENT(in   ) ::   Kbb, Kmm, Kaa  ! ocean time level indices 
     238      REAL(wp), DIMENSION(jpi,jpj,jpt), INTENT(inout) ::   pssh           ! SSH field 
    241239      ! 
    242240      REAL(wp) ::   zcoef   ! local scalar 
    243241      !!---------------------------------------------------------------------- 
    244242      ! 
    245       IF( ln_timing )   CALL timing_start('ssh_swp') 
     243      IF( ln_timing )   CALL timing_start('ssh_atf') 
    246244      ! 
    247245      IF( kt == nit000 ) THEN 
    248246         IF(lwp) WRITE(numout,*) 
    249          IF(lwp) WRITE(numout,*) 'ssh_swp : Asselin time filter and swap of sea surface height' 
     247         IF(lwp) WRITE(numout,*) 'ssh_atf : Asselin time filter of sea surface height' 
    250248         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    251249      ENDIF 
    252250      !              !==  Euler time-stepping: no filter, just swap  ==! 
    253       IF ( neuler == 0 .AND. kt == nit000 ) THEN 
    254          ssh(:,:,Kmm) = ssh(:,:,Kaa)                              ! now    <-- after  (before already = now) 
    255          ! 
    256       ELSE           !==  Leap-Frog time-stepping: Asselin filter + swap  ==! 
    257          !                                                  ! before <-- now filtered 
    258          ssh(:,:,Kbb) = ssh(:,:,Kmm) + atfp * ( ssh(:,:,Kbb) - 2 * ssh(:,:,Kmm) + ssh(:,:,Kaa) ) 
    259          IF( .NOT.ln_linssh ) THEN                          ! before <-- with forcing removed 
     251      IF ( .NOT.( neuler == 0 .AND. kt == nit000 ) ) THEN   ! Only do time filtering for leapfrog timesteps 
     252         !                                                  ! filtered "now" field 
     253         pssh(:,:,Kmm) = pssh(:,:,Kmm) + atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
     254         IF( .NOT.ln_linssh ) THEN                          ! "now" <-- with forcing removed 
    260255            zcoef = atfp * rdt * r1_rau0 
    261             ssh(:,:,Kbb) = ssh(:,:,Kbb) - zcoef * (     emp_b(:,:) - emp   (:,:)   & 
     256            pssh(:,:,Kmm) = pssh(:,:,Kmm) - zcoef * (     emp_b(:,:) - emp   (:,:)   & 
    262257               &                             -    rnf_b(:,:) + rnf   (:,:)   & 
    263258               &                             + fwfisf_b(:,:) - fwfisf(:,:)   ) * ssmask(:,:) 
    264259         ENDIF 
    265          ssh(:,:,Kmm) = ssh(:,:,Kaa)                              ! now <-- after 
    266       ENDIF 
    267       ! 
    268       IF(ln_ctl)   CALL prt_ctl( tab2d_1=ssh(:,:,Kbb), clinfo1=' ssh(:,:,Kbb)  - : ', mask1=tmask ) 
    269       ! 
    270       IF( ln_timing )   CALL timing_stop('ssh_swp') 
    271       ! 
    272    END SUBROUTINE ssh_swp 
     260      ENDIF 
     261      ! 
     262      IF(ln_ctl)   CALL prt_ctl( tab2d_1=pssh(:,:,Kmm), clinfo1=' pssh(:,:,Kmm)  - : ', mask1=tmask ) 
     263      ! 
     264      IF( ln_timing )   CALL timing_stop('ssh_atf') 
     265      ! 
     266   END SUBROUTINE ssh_atf 
    273267 
    274268   SUBROUTINE wAimp( kt, Kmm ) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/SBC/fldread.F90

    r10922 r11480  
    130130CONTAINS 
    131131 
    132    SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset, jpk_bdy, fvl ) 
     132   SUBROUTINE fld_read( kt, kn_fsbc, sd, map, kit, kt_offset, jpk_bdy, fvl, Kmm ) 
    133133      !!--------------------------------------------------------------------- 
    134134      !!                    ***  ROUTINE fld_read  *** 
    … …  
    153153      INTEGER  , INTENT(in   ), OPTIONAL     ::   jpk_bdy   ! number of vertical levels in the BDY data 
    154154      LOGICAL  , INTENT(in   ), OPTIONAL     ::   fvl   ! number of vertical levels in the BDY data 
     155      INTEGER  , INTENT(in   ), OPTIONAL     ::   Kmm   ! ocean time level index 
    155156      !! 
    156157      INTEGER  ::   itmp         ! local variable 
    … …  
    287288               ! read after data 
    288289               IF( PRESENT(jpk_bdy) ) THEN 
    289                   CALL fld_get( sd(jf), imap, jpk_bdy, fvl) 
     290                  CALL fld_get( sd(jf), imap, jpk_bdy, fvl, Kmm ) 
    290291               ELSE 
    291292                  CALL fld_get( sd(jf), imap ) 
    … …  
    614615 
    615616 
    616    SUBROUTINE fld_get( sdjf, map, jpk_bdy, fvl ) 
     617   SUBROUTINE fld_get( sdjf, map, jpk_bdy, fvl, Kmm ) 
    617618      !!--------------------------------------------------------------------- 
    618619      !!                    ***  ROUTINE fld_get  *** 
    … …  
    624625      INTEGER  , INTENT(in), OPTIONAL  ::   jpk_bdy ! number of vertical levels in the bdy data 
    625626      LOGICAL  , INTENT(in), OPTIONAL  ::   fvl     ! number of vertical levels in the bdy data 
     627      INTEGER  , INTENT(in), OPTIONAL  ::   Kmm     ! ocean time level index 
    626628      ! 
    627629      INTEGER ::   ipk      ! number of vertical levels of sdjf%fdta ( 2D: ipk=1 ; 3D: ipk=jpk ) 
    … …  
    638640         IF( PRESENT(jpk_bdy) ) THEN 
    639641            IF( sdjf%ln_tint ) THEN   ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fdta(:,:,:,2),                & 
    640                                                         sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 
     642                                                        sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl, Kmm ) 
    641643            ELSE                      ;   CALL fld_map( sdjf%num, sdjf%clvar, sdjf%fnow(:,:,:  ),                & 
    642                                                         sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl ) 
     644                                                        sdjf%nrec_a(1), map, sdjf%igrd, sdjf%ibdy, jpk_bdy, fvl, Kmm ) 
    643645            ENDIF 
    644646         ELSE 
    … …  
    701703   END SUBROUTINE fld_get 
    702704 
    703    SUBROUTINE fld_map( num, clvar, dta, nrec, map, igrd, ibdy, jpk_bdy, fvl ) 
     705   SUBROUTINE fld_map( num, clvar, dta, nrec, map, igrd, ibdy, jpk_bdy, fvl, Kmm ) 
    704706      !!--------------------------------------------------------------------- 
    705707      !!                    ***  ROUTINE fld_map  *** 
    … …  
    718720      INTEGER  , INTENT(in), OPTIONAL         ::   igrd, ibdy, jpk_bdy  ! grid type, set number and number of vertical levels in the bdy data 
    719721      LOGICAL  , INTENT(in), OPTIONAL         ::   fvl     ! grid type, set number and number of vertical levels in the bdy data 
     722      INTEGER  , INTENT(in), OPTIONAL         ::   Kmm     ! ocean time level index  
    720723      INTEGER                                 ::   jpkm1_bdy! number of vertical levels in the bdy data minus 1 
    721724      !! 
    … …  
    813816 
    814817      IF ( ln_bdy ) &  
    815          CALL fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta) 
     818         CALL fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta, Kmm) 
    816819 
    817820      ELSE ! boundary data assumed to be on model grid 
    … …  
    838841   END SUBROUTINE fld_map 
    839842    
    840    SUBROUTINE fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta) 
     843   SUBROUTINE fld_bdy_interp(dta_read, dta_read_z, dta_read_dz, map, jpk_bdy, igrd, ibdy, fv, dta, fvl, ilendta, Kmm) 
    841844 
    842845      !!--------------------------------------------------------------------- 
    … …  
    857860      INTEGER  , INTENT(in)                   ::   igrd, ibdy, jpk_bdy        ! number of levels in bdy data 
    858861      INTEGER  , INTENT(in)                   ::   ilendta                    ! length of data in file 
     862      INTEGER  , INTENT(in), OPTIONAL         ::   Kmm                        ! ocean time level index 
    859863      !! 
    860864      INTEGER                                 ::   ipi                        ! length of boundary data on local process 
    … …  
    900904            SELECT CASE( igrd )                          
    901905               CASE(1) 
    902                   IF( ABS( (zh - ht_n(zij,zjj)) / ht_n(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 
     906                  IF( ABS( (zh - ht(zij,zjj)) / ht(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 
    903907                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    904908                     CALL ctl_warn('fld_bdy_interp: T depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    905                  !   IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t(zij,zjj,:,nfld_Nnn), mask=tmask(zij,zjj,:)==1),  ht_n(zij,zjj), map%ptr(ib), ib, zij, zjj 
     909                 !   IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t(zij,zjj,:,Kmm), mask=tmask(zij,zjj,:)==1),  ht(zij,zjj), map%ptr(ib), ib, zij, zjj 
    906910                  ENDIF 
    907911               CASE(2) 
    908                   IF( ABS( (zh - hu_n(zij,zjj)) * r1_hu_n(zij,zjj)) * umask(zij,zjj,1) > 0.01_wp ) THEN 
     912                  IF( ABS( (zh - hu(zij,zjj,Kmm)) * r1_hu(zij,zjj,Kmm)) * umask(zij,zjj,1) > 0.01_wp ) THEN 
    909913                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    910914                     CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    911                      IF(lwp) WRITE(*,*) 'DEPTHU', zh, sum(e3u(zij,zjj,:,nfld_Nnn), mask=umask(zij,zjj,:)==1),  sum(umask(zij,zjj,:)), & 
    912                        &                hu_n(zij,zjj), map%ptr(ib), ib, zij, zjj, narea-1  , & 
     915                     IF(lwp) WRITE(*,*) 'DEPTHU', zh, sum(e3u(zij,zjj,:,Kmm), mask=umask(zij,zjj,:)==1),  sum(umask(zij,zjj,:)), & 
     916                       &                hu(zij,zjj,Kmm), map%ptr(ib), ib, zij, zjj, narea-1  , & 
    913917                        &                dta_read(map%ptr(ib),1,:) 
    914918                  ENDIF 
    915919               CASE(3) 
    916                   IF( ABS( (zh - hv_n(zij,zjj)) * r1_hv_n(zij,zjj)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 
     920                  IF( ABS( (zh - hv(zij,zjj,Kmm)) * r1_hv(zij,zjj,Kmm)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 
    917921                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    918922                     CALL ctl_warn('fld_bdy_interp: V depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    … …  
    922926               SELECT CASE( igrd )                        
    923927                  CASE(1) 
    924                      zl =  gdept(zij,zjj,ik,nfld_Nnn)                                          ! if using in step could use fsdept instead of gdept_n? 
     928                     zl =  gdept(zij,zjj,ik,Kmm)                                          ! if using in step could use fsdept instead of gdept_n? 
    925929                  CASE(2) 
    926930                     IF(ln_sco) THEN 
    927                         zl =  ( gdept(zij,zjj,ik,nfld_Nnn) + gdept(zij+1,zjj,ik,nfld_Nnn) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
     931                        zl =  ( gdept(zij,zjj,ik,Kmm) + gdept(zij+1,zjj,ik,Kmm) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
    928932                     ELSE 
    929                         zl =  MIN( gdept(zij,zjj,ik,nfld_Nnn), gdept(zij+1,zjj,ik,nfld_Nnn) )  
     933                        zl =  MIN( gdept(zij,zjj,ik,Kmm), gdept(zij+1,zjj,ik,Kmm) )  
    930934                     ENDIF 
    931935                  CASE(3) 
    932936                     IF(ln_sco) THEN 
    933                         zl =  ( gdept(zij,zjj,ik,nfld_Nnn) + gdept(zij,zjj+1,ik,nfld_Nnn) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
     937                        zl =  ( gdept(zij,zjj,ik,Kmm) + gdept(zij,zjj+1,ik,Kmm) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
    934938                     ELSE 
    935                         zl =  MIN( gdept(zij,zjj,ik,nfld_Nnn), gdept(zij,zjj+1,ik,nfld_Nnn) ) 
     939                        zl =  MIN( gdept(zij,zjj,ik,Kmm), gdept(zij,zjj+1,ik,Kmm) ) 
    936940                     ENDIF 
    937941               END SELECT 
    … …  
    941945                  dta(ib,1,ik) =  dta_read(map%ptr(ib),1,MAXLOC(dta_read_z(map%ptr(ib),1,:),1)) 
    942946               ELSE                                                                                ! inbetween : vertical interpolation between ikk & ikk+1 
    943                   DO ikk = 1, jpkm1_bdy                                                            ! when  gdept(ikk,nfld_Nnn) < zl < gdept(ikk+1,nfld_Nnn) 
     947                  DO ikk = 1, jpkm1_bdy                                                            ! when  gdept(ikk,Kmm) < zl < gdept(ikk+1,Kmm) 
    944948                     IF( ( (zl-dta_read_z(map%ptr(ib),1,ikk)) * (zl-dta_read_z(map%ptr(ib),1,ikk+1)) <= 0._wp) & 
    945949                    &    .AND. (dta_read_z(map%ptr(ib),1,ikk+1) /= fv_alt)) THEN 
    … …  
    965969              ENDDO 
    966970              DO ik = 1, ipk                                ! calculate transport on model grid 
    967                   ztrans_new = ztrans_new + dta(ib,1,ik) * e3u(zij,zjj,ik,nfld_Nnn) * umask(zij,zjj,ik) 
     971                  ztrans_new = ztrans_new + dta(ib,1,ik) * e3u(zij,zjj,ik,Kmm) * umask(zij,zjj,ik) 
    968972              ENDDO 
    969973              DO ik = 1, ipk                                ! make transport correction 
    970974                 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 
    971                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 
     975                    dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu(zij,zjj,Kmm) ) * umask(zij,zjj,ik) 
    972976                 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 
    973                     IF( ABS(ztrans * r1_hu_n(zij,zjj)) > 0.01_wp ) & 
     977                    IF( ABS(ztrans * r1_hu(zij,zjj,Kmm)) > 0.01_wp ) & 
    974978                   &   CALL ctl_warn('fld_bdy_interp: barotropic component of > 0.01 ms-1 found in baroclinic velocities at') 
    975                     dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hu_n(zij,zjj) * umask(zij,zjj,ik) 
     979                    dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hu(zij,zjj,Kmm) * umask(zij,zjj,ik) 
    976980                 ENDIF 
    977981              ENDDO 
    … …  
    990994              ENDDO 
    991995              DO ik = 1, ipk                                ! calculate transport on model grid 
    992                   ztrans_new = ztrans_new + dta(ib,1,ik) * e3v(zij,zjj,ik,nfld_Nnn) * vmask(zij,zjj,ik) 
     996                  ztrans_new = ztrans_new + dta(ib,1,ik) * e3v(zij,zjj,ik,Kmm) * vmask(zij,zjj,ik) 
    993997              ENDDO 
    994998              DO ik = 1, ipk                                ! make transport correction 
    995999                 IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 
    996                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 
     1000                    dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv(zij,zjj,Kmm) ) * vmask(zij,zjj,ik) 
    9971001                 ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 
    998                     dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hv_n(zij,zjj) * vmask(zij,zjj,ik) 
     1002                    dta(ib,1,ik) = dta(ib,1,ik) + ( 0._wp - ztrans_new ) * r1_hv(zij,zjj,Kmm) * vmask(zij,zjj,ik) 
    9991003                 ENDIF 
    10001004              ENDDO 
    … …  
    10251029            SELECT CASE( igrd )                          
    10261030               CASE(1) 
    1027                   IF( ABS( (zh - ht_n(zij,zjj)) / ht_n(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 
     1031                  IF( ABS( (zh - ht(zij,zjj)) / ht(zij,zjj)) * tmask(zij,zjj,1) > 0.01_wp ) THEN 
    10281032                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    10291033                     CALL ctl_warn('fld_bdy_interp: T depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    1030                  !   IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t(zij,zjj,:,nfld_Nnn), mask=tmask(zij,zjj,:)==1),  ht_n(zij,zjj), map%ptr(ib), ib, zij, zjj 
     1034                 !   IF(lwp) WRITE(*,*) 'DEPTHT', zh, sum(e3t(zij,zjj,:,Kmm), mask=tmask(zij,zjj,:)==1),  ht(zij,zjj), map%ptr(ib), ib, zij, zjj 
    10311035                  ENDIF 
    10321036               CASE(2) 
    1033                   IF( ABS( (zh - hu_n(zij,zjj)) * r1_hu_n(zij,zjj)) * umask(zij,zjj,1) > 0.01_wp ) THEN 
     1037                  IF( ABS( (zh - hu(zij,zjj,Kmm)) * r1_hu(zij,zjj,Kmm)) * umask(zij,zjj,1) > 0.01_wp ) THEN 
    10341038                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    10351039                     CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    10361040                  ENDIF 
    10371041               CASE(3) 
    1038                   IF( ABS( (zh - hv_n(zij,zjj)) * r1_hv_n(zij,zjj)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 
     1042                  IF( ABS( (zh - hv(zij,zjj,Kmm)) * r1_hv(zij,zjj,Kmm)) * vmask(zij,zjj,1) > 0.01_wp ) THEN 
    10391043                     WRITE(ibstr,"(I10.10)") map%ptr(ib)  
    10401044                     CALL ctl_warn('fld_bdy_interp: V depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 
    … …  
    10441048               SELECT CASE( igrd )                          ! coded for sco - need zco and zps option using min 
    10451049                  CASE(1) 
    1046                      zl =  gdept(zij,zjj,ik,nfld_Nnn)              ! if using in step could use fsdept instead of gdept_n? 
     1050                     zl =  gdept(zij,zjj,ik,Kmm)              ! if using in step could use fsdept instead of gdept_n? 
    10471051                  CASE(2) 
    10481052                     IF(ln_sco) THEN 
    1049                         zl =  ( gdept(zij,zjj,ik,nfld_Nnn) + gdept(zij+1,zjj,ik,nfld_Nnn) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
     1053                        zl =  ( gdept(zij,zjj,ik,Kmm) + gdept(zij+1,zjj,ik,Kmm) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
    10501054                     ELSE 
    1051                         zl =  MIN( gdept(zij,zjj,ik,nfld_Nnn), gdept(zij+1,zjj,ik,nfld_Nnn) ) 
     1055                        zl =  MIN( gdept(zij,zjj,ik,Kmm), gdept(zij+1,zjj,ik,Kmm) ) 
    10521056                     ENDIF 
    10531057                  CASE(3) 
    10541058                     IF(ln_sco) THEN 
    1055                         zl =  ( gdept(zij,zjj,ik,nfld_Nnn) + gdept(zij,zjj+1,ik,nfld_Nnn) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
     1059                        zl =  ( gdept(zij,zjj,ik,Kmm) + gdept(zij,zjj+1,ik,Kmm) ) * 0.5_wp  ! if using in step could use fsdept instead of gdept_n? 
    10561060                     ELSE 
    1057                         zl =  MIN( gdept(zij,zjj,ik,nfld_Nnn), gdept(zij,zjj+1,ik,nfld_Nnn) ) 
     1061                        zl =  MIN( gdept(zij,zjj,ik,Kmm), gdept(zij,zjj+1,ik,Kmm) ) 
    10581062                     ENDIF 
    10591063               END SELECT 
    … …  
    10631067                  dta(ib,1,ik) =  dta_read(ji,jj,MAXLOC(dta_read_z(ji,jj,:),1)) 
    10641068               ELSE                                                                     ! inbetween : vertical interpolation between ikk & ikk+1 
    1065                   DO ikk = 1, jpkm1_bdy                                                 ! when  gdept(ikk,nfld_Nnn) < zl < gdept(ikk+1,nfld_Nnn) 
     1069                  DO ikk = 1, jpkm1_bdy                                                 ! when  gdept(ikk,Kmm) < zl < gdept(ikk+1,Kmm) 
    10661070                     IF( ( (zl-dta_read_z(ji,jj,ikk)) * (zl-dta_read_z(ji,jj,ikk+1)) <= 0._wp) & 
    10671071                    &    .AND. (dta_read_z(ji,jj,ikk+1) /= fv_alt)) THEN 
    … …  
    10891093               ENDDO 
    10901094               DO ik = 1, ipk                                ! calculate transport on model grid 
    1091                   ztrans_new = ztrans_new + dta(ib,1,ik) * e3u(zij,zjj,ik,nfld_Nnn) * umask(zij,zjj,ik) 
     1095                  ztrans_new = ztrans_new + dta(ib,1,ik) * e3u(zij,zjj,ik,Kmm) * umask(zij,zjj,ik) 
    10921096               ENDDO 
    10931097               DO ik = 1, ipk                                ! make transport correction 
    10941098                  IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 
    1095                      dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 
     1099                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hu(zij,zjj,Kmm) ) * umask(zij,zjj,ik) 
    10961100                  ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 
    1097                      dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp  - ztrans_new ) * r1_hu_n(zij,zjj) ) * umask(zij,zjj,ik) 
     1101                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp  - ztrans_new ) * r1_hu(zij,zjj,Kmm) ) * umask(zij,zjj,ik) 
    10981102                  ENDIF 
    10991103               ENDDO 
    … …  
    11141118               ENDDO 
    11151119               DO ik = 1, ipk                                ! calculate transport on model grid 
    1116                   ztrans_new = ztrans_new + dta(ib,1,ik) * e3v(zij,zjj,ik,nfld_Nnn) * vmask(zij,zjj,ik) 
     1120                  ztrans_new = ztrans_new + dta(ib,1,ik) * e3v(zij,zjj,ik,Kmm) * vmask(zij,zjj,ik) 
    11171121               ENDDO 
    11181122               DO ik = 1, ipk                                ! make transport correction 
    11191123                  IF(fvl) THEN ! bdy data are total velocity so adjust bt transport term to match input data 
    1120                      dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 
     1124                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( ztrans - ztrans_new ) * r1_hv(zij,zjj,Kmm) ) * vmask(zij,zjj,ik) 
    11211125                  ELSE ! we're just dealing with bc velocity so bt transport term should sum to zero 
    1122                      dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp  - ztrans_new ) * r1_hv_n(zij,zjj) ) * vmask(zij,zjj,ik) 
     1126                     dta(ib,1,ik) = ( dta(ib,1,ik) + ( 0._wp  - ztrans_new ) * r1_hv(zij,zjj,Kmm) ) * vmask(zij,zjj,ik) 
    11231127                  ENDIF 
    11241128               ENDDO 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/SBC/sbcmod.F90

    r11027 r11480  
    442442      CASE(  1 )   ;         CALL sbc_ice_if   ( kt, Kbb, Kmm )   ! Ice-cover climatology ("Ice-if" model) 
    443443#if defined key_si3 
    444       CASE(  2 )   ;         CALL ice_stp  ( kt, Kbb, nsbc )  ! SI3 ice model 
     444      CASE(  2 )   ;         CALL ice_stp  ( kt, Kbb, Kmm, nsbc ) ! SI3 ice model 
    445445#endif 
    446446      CASE(  3 )   ;         CALL sbc_ice_cice ( kt, nsbc )       ! CICE ice model 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trasbc.F90

    r10985 r11480  
    233233            DO jj = 2, jpj  
    234234               DO ji = fs_2, fs_jpim1 
    235                   ztim = ssh_iau(ji,jj) / ( ht_n(ji,jj) + 1. - ssmask(ji, jj) ) 
     235                  ztim = ssh_iau(ji,jj) / ( ht(ji,jj) + 1. - ssmask(ji, jj) ) 
    236236                  pts(ji,jj,:,jp_tem,Krhs) = pts(ji,jj,:,jp_tem,Krhs) + pts(ji,jj,:,jp_tem,Kmm) * ztim 
    237237                  pts(ji,jj,:,jp_sal,Krhs) = pts(ji,jj,:,jp_sal,Krhs) + pts(ji,jj,:,jp_sal,Kmm) * ztim 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trddyn.F90

    r10946 r11480  
    123123                              z3dx(:,:,:) = 0._wp                  ! U.dxU & V.dyV (approximation) 
    124124                              z3dy(:,:,:) = 0._wp 
    125                               DO jk = 1, jpkm1   ! no mask as un,vn are masked 
     125                              DO jk = 1, jpkm1   ! no mask as uu, vv are masked 
    126126                                 DO jj = 2, jpjm1 
    127127                                    DO ji = 2, jpim1 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRD/trdvor.F90

    r10946 r11480  
    189189 
    190190      ! Average except for Beta.V 
    191       zudpvor(:,:) = zudpvor(:,:) * r1_hu_n(:,:) 
    192       zvdpvor(:,:) = zvdpvor(:,:) * r1_hv_n(:,:) 
     191      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm) 
     192      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm) 
    193193    
    194194      ! Curl 
    … …  
    276276         END DO 
    277277         ! Average of the Curl and Surface mask 
    278          vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * r1_hu_n(:,:) * fmask(:,:,1) 
     278         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * r1_hu(:,:,Kmm) * fmask(:,:,1) 
    279279      ENDIF 
    280280      ! 
    281281      ! Average  
    282       zudpvor(:,:) = zudpvor(:,:) * r1_hu_n(:,:) 
    283       zvdpvor(:,:) = zvdpvor(:,:) * r1_hv_n(:,:) 
     282      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm) 
     283      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm) 
    284284      ! 
    285285      ! Curl 
    … …  
    342342      END DO 
    343343  
    344       zuu(:,:) = zuu(:,:) * r1_hu_n(:,:) 
    345       zvv(:,:) = zvv(:,:) * r1_hv_n(:,:) 
     344      zuu(:,:) = zuu(:,:) * r1_hu(:,:,Kmm) 
     345      zvv(:,:) = zvv(:,:) * r1_hv(:,:,Kmm) 
    346346 
    347347      ! Curl 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/ZDF/zdfosm.F90

    r10955 r11480  
    489489 
    490490      zhbl_t(:,:) = hbl(:,:) + (zdhdt(:,:) - ww(ji,jj,ibld(ji,jj)))* rn_rdt ! certainly need wb here, so subtract it 
    491       zhbl_t(:,:) = MIN(zhbl_t(:,:), ht_n(:,:)) 
     491      zhbl_t(:,:) = MIN(zhbl_t(:,:), ht(:,:)) 
    492492      zdhdt(:,:) = MIN(zdhdt(:,:), (zhbl_t(:,:) - hbl(:,:))/rn_rdt + ww(ji,jj,ibld(ji,jj))) ! adjustment to represent limiting by ocean bottom 
    493493 
    … …  
    525525 
    526526                     zhbl_s = zhbl_s + MIN( - zwb_ent(ji,jj) / zdb * rn_rdt / FLOAT(ibld(ji,jj)-imld(ji,jj) ), e3w(ji,jj,jk,Kmm) ) 
    527                      zhbl_s = MIN(zhbl_s, ht_n(ji,jj)) 
     527                     zhbl_s = MIN(zhbl_s, ht(ji,jj)) 
    528528 
    529529                     IF ( zhbl_s >= gdepw(ji,jj,jm+1,Kmm) ) jm = jm + 1 
    … …  
    546546                          &          * zwstrl(ji,jj)**3 / hbli(ji,jj) ) / zdb * e3w(ji,jj,jk,Kmm) / zdhdt(ji,jj)  ! ALMG to investigate whether need to include ww here 
    547547 
    548                      zhbl_s = MIN(zhbl_s, ht_n(ji,jj)) 
     548                     zhbl_s = MIN(zhbl_s, ht(ji,jj)) 
    549549                     IF ( zhbl_s >= gdepw(ji,jj,jm,Kmm) ) jm = jm + 1 
    550550                  END DO 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/nemogcm.F90

    r11421 r11480  
    473473#if defined key_top 
    474474      !                                      ! Passive tracers 
    475                            CALL     trc_init( Nbb, Nnn, Naa ) 
     475                           CALL     trc_init 
    476476#endif 
    477477      IF( l_ldfslp     )   CALL ldf_slp_init    ! slope of lateral mixing 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/oce.F90

    r11047 r11480  
    88   !!            3.3  !  2010-09  (C. Ethe) TRA-TRC merge: add ts, gtsu, gtsv 4D arrays 
    99   !!            3.7  !  2014-01  (G. Madec) suppression of curl and before hdiv from in-core memory 
     10   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rename prognostic variables in preparation for new time scheme 
    1011   !!---------------------------------------------------------------------- 
    1112   USE par_oce        ! ocean parameters 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90

    r11421 r11480  
    3131   !!             -   !  2015-11  (J. Chanut) free surface simplification (remove filtered free surface) 
    3232   !!            4.0  !  2017-05  (G. Madec)  introduction of the vertical physics manager (zdfphy) 
     33   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rewrite in preparation for new timestepping scheme 
    3334   !!---------------------------------------------------------------------- 
    3435 
    … …  
    265266 
    266267      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    267       ! Set boundary conditions and Swap 
     268      ! Set boundary conditions, time filter and swap time levels 
    268269      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    269270!!jc1: For agrif, it would be much better to finalize tracers/momentum here (e.g. bdy conditions) and move the swap  
    … …  
    281282!!  
    282283!!jc2: dynnxt must be the latest call. e3t(:,:,:,Nbb) are indeed updated in that routine 
    283                          CALL tra_nxt       ( kstp, Nbb, Nnn, Nrhs, Naa )  ! finalize (bcs) tracer fields at next time step and swap 
    284                          CALL dyn_nxt       ( kstp, Nbb, Nnn, Nrhs, uu, vv, Naa  )  ! finalize (bcs) velocities at next time step and swap (always called after tra_nxt) 
    285                          CALL ssh_swp       ( kstp, Nbb, Nnn, Naa )  ! swap of sea surface height 
    286       IF(.NOT.ln_linssh) CALL dom_vvl_sf_swp( kstp, Nbb, Nnn, Naa )  ! swap of vertical scale factors 
     284                         CALL tra_atf       ( kstp, Nbb, Nnn, Naa, ts )                      ! time filtering of "now" tracer arrays 
     285                         CALL dyn_atf       ( kstp, Nbb, Nnn, Naa, uu, vv, e3t, e3u, e3v  )  ! time filtering of "now" velocities and scale factors 
     286                         CALL ssh_atf       ( kstp, Nbb, Nnn, Naa, ssh )                     ! time filtering of "now" sea surface height 
     287      ! 
     288      ! Swap time levels 
     289      Nrhs = Nbb 
     290      Nbb = Nnn 
     291      Nnn = Naa 
     292      Naa = Nrhs 
     293      ! 
     294      IF(.NOT.ln_linssh) CALL dom_vvl_sf_update( kstp, Nbb, Nnn, Naa )  ! recompute vertical scale factors 
    287295      ! 
    288296      IF( ln_diahsb  )   CALL dia_hsb       ( kstp, Nbb, Nnn )  ! - ML - global conservation diagnostics 
    … …  
    300308      ! AGRIF 
    301309      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<       
     310                         Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs      ! agrif_oce module copies of time level indices 
    302311                         CALL Agrif_Integrate_ChildGrids( stp )       ! allows to finish all the Child Grids before updating 
    303312 
    304313                         IF( Agrif_NbStepint() == 0 ) THEN 
    305                             Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices 
    306314                            CALL Agrif_update_all( )                  ! Update all components 
    307315                         ENDIF 
    … …  
    312320      ! Control 
    313321      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    314                          CALL stp_ctl      ( kstp, Nnn, indic ) 
     322                         CALL stp_ctl      ( kstp, Nbb, Nnn, indic ) 
    315323                          
    316324      IF( kstp == nit000 ) THEN                          ! 1st time step only 
    … …  
    337345      ! 
    338346   END SUBROUTINE stp 
    339     
     347   ! 
    340348   !!====================================================================== 
    341349END MODULE step 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step_oce.F90

    r10927 r11480  
    3030   USE traldf          ! lateral mixing                   (tra_ldf routine) 
    3131   USE trazdf          ! vertical mixing                  (tra_zdf routine) 
    32    USE tranxt          ! time-stepping                    (tra_nxt routine) 
     32   USE traatf          ! time filtering                   (tra_atf routine) 
    3333   USE tranpc          ! non-penetrative convection       (tra_npc routine) 
    3434 
    … …  
    4343   USE dynspg          ! surface pressure gradient        (dyn_spg routine) 
    4444 
    45    USE dynnxt          ! time-stepping                    (dyn_nxt routine) 
     45   USE dynatf          ! time-filtering                   (dyn_atf routine) 
    4646 
    4747   USE stopar          ! Stochastic parametrization       (sto_par routine) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/stpctl.F90

    r10965 r11480  
    4242CONTAINS 
    4343 
    44    SUBROUTINE stp_ctl( kt, Kmm, kindic ) 
     44   SUBROUTINE stp_ctl( kt, Kbb, Kmm, kindic ) 
    4545      !!---------------------------------------------------------------------- 
    4646      !!                    ***  ROUTINE stp_ctl  *** 
    … …  
    6060      !!---------------------------------------------------------------------- 
    6161      INTEGER, INTENT(in   ) ::   kt       ! ocean time-step index 
    62       INTEGER, INTENT(in   ) ::   Kmm      ! ocean time level index 
     62      INTEGER, INTENT(in   ) ::   Kbb, Kmm      ! ocean time level index 
    6363      INTEGER, INTENT(inout) ::   kindic   ! error indicator 
    6464      !! 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OFF/dtadyn.F90

    r11027 r11480  
    4646   PRIVATE 
    4747 
    48    PUBLIC   dta_dyn_init       ! called by opa.F90 
    49    PUBLIC   dta_dyn            ! called by step.F90 
    50    PUBLIC   dta_dyn_sed_init   ! called by opa.F90 
    51    PUBLIC   dta_dyn_sed        ! called by step.F90 
    52    PUBLIC   dta_dyn_swp        ! called by step.F90 
     48   PUBLIC   dta_dyn_init       ! called by nemo_init 
     49   PUBLIC   dta_dyn            ! called by nemo_gcm 
     50   PUBLIC   dta_dyn_sed_init   ! called by nemo_init 
     51   PUBLIC   dta_dyn_sed        ! called by nemo_gcm 
     52   PUBLIC   dta_dyn_atf        ! called by nemo_gcm 
     53   PUBLIC   dta_dyn_sf_interp  ! called by nemo_gcm 
    5354 
    5455   CHARACTER(len=100) ::   cn_dir          !: Root directory for location of ssr files 
    … …  
    535536   END SUBROUTINE dta_dyn_sed_init 
    536537 
    537    SUBROUTINE dta_dyn_swp( kt, Kbb, Kmm, Kaa ) 
     538   SUBROUTINE dta_dyn_atf( kt, Kbb, Kmm, Kaa ) 
    538539     !!--------------------------------------------------------------------- 
    539540      !!                    ***  ROUTINE dta_dyn_swp  *** 
    540541      !! 
    541       !! ** Purpose :   Swap and the data and compute the vertical scale factor  
    542       !!              at U/V/W pointand the depht 
     542      !! ** Purpose :   Asselin time filter of now SSH 
    543543      !!--------------------------------------------------------------------- 
    544544      INTEGER, INTENT(in) :: kt             ! time step 
    545545      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa  ! ocean time level indices 
    546546      ! 
     547      !!--------------------------------------------------------------------- 
     548 
     549      IF( kt == nit000 ) THEN 
     550         IF(lwp) WRITE(numout,*) 
     551         IF(lwp) WRITE(numout,*) 'dta_dyn_atf : Asselin time filter of sea surface height' 
     552         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
     553      ENDIF 
     554 
     555      ssh(:,:,Kmm) = ssh(:,:,Kmm) + atfp * ( ssh(:,:,Kbb) - 2 * ssh(:,:,Kmm) + ssh(:,:,Kaa))   
     556 
     557      !! Do we also need to time filter e3t?? 
     558      ! 
     559   END SUBROUTINE dta_dyn_atf 
     560    
     561   SUBROUTINE dta_dyn_sf_interp( kt, Kmm ) 
     562      !!--------------------------------------------------------------------- 
     563      !!                    ***  ROUTINE dta_dyn_sf_interp  *** 
     564      !! 
     565      !! ** Purpose :   Calculate scale factors at U/V/W points and depths 
     566      !!                given the after e3t field 
     567      !!--------------------------------------------------------------------- 
     568      INTEGER, INTENT(in) :: kt   ! time step 
     569      INTEGER, INTENT(in) :: Kmm  ! ocean time level indices 
     570      ! 
    547571      INTEGER             :: ji, jj, jk 
    548572      REAL(wp)            :: zcoef 
    549573      !!--------------------------------------------------------------------- 
    550  
    551       IF( kt == nit000 ) THEN 
    552          IF(lwp) WRITE(numout,*) 
    553          IF(lwp) WRITE(numout,*) 'ssh_swp : Asselin time filter and swap of sea surface height' 
    554          IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    555       ENDIF 
    556  
    557       ssh(:,:,Kbb) = ssh(:,:,Kmm) + atfp * ( ssh(:,:,Kbb) - 2 * ssh(:,:,Kmm) + ssh(:,:,Kaa))  ! before <-- now filtered 
    558       ssh(:,:,Kmm) = ssh(:,:,Kaa) 
    559  
    560       e3t(:,:,:,Kmm) = e3t(:,:,:,Kaa) 
    561  
    562       ! Reconstruction of all vertical scale factors at now and before time steps 
    563       ! ============================================================================= 
    564574 
    565575      ! Horizontal scale factor interpolations 
    … …  
    571581      ! ------------------------------------ 
    572582      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W' ) 
    573  
    574       e3t(:,:,:,Kbb)  = e3t(:,:,:,Kmm) 
    575       e3u(:,:,:,Kbb)  = e3u(:,:,:,Kmm) 
    576       e3v(:,:,:,Kbb)  = e3v(:,:,:,Kmm) 
    577583 
    578584      ! t- and w- points depth 
    … …  
    592598      END DO 
    593599      ! 
    594       gdept(:,:,:,Kbb) = gdept(:,:,:,Kmm) 
    595       gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm) 
    596       ! 
    597    END SUBROUTINE dta_dyn_swp 
    598     
     600   END SUBROUTINE dta_dyn_sf_interp 
    599601 
    600602   SUBROUTINE dta_dyn_ssh( kt, phdivtr, psshb,  pemp, pssha, pe3ta ) 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OFF/nemogcm.F90

    r11027 r11480  
    119119#else 
    120120                                CALL dta_dyn    ( istp, Nbb, Nnn, Naa )       ! Interpolation of the dynamical fields 
    121          IF( .NOT.ln_linssh )   CALL dta_dyn_swp( istp, Nbb, Nnn, Naa )       ! swap of sea  surface height and vertical scale factors 
    122121#endif 
    123122                                CALL trc_stp    ( istp, Nbb, Nnn, Nrhs, Naa ) ! time-stepping 
     123#if ! defined key_sed_off 
     124         IF( .NOT.ln_linssh )   CALL dta_dyn_atf( istp, Nbb, Nnn, Naa )       ! time filter of sea  surface height and vertical scale factors 
     125#endif 
     126         ! Swap time levels 
     127         Nrhs = Nbb 
     128         Nbb = Nnn 
     129         Nnn = Naa 
     130         Naa = Nrhs 
     131         ! 
     132#if ! defined key_sed_off 
     133         IF( .NOT.ln_linssh )   CALL dta_dyn_sf_interp( istp, Nnn )  ! calculate now grid parameters 
     134#endif 
    124135                                CALL stp_ctl    ( istp, indic )  ! Time loop: control and print 
    125136         istp = istp + 1 
    … …  
    333344#endif 
    334345 
    335                            CALL     trc_init( Nbb, Nnn, Naa )   ! Passive tracers initialization 
     346                           CALL     trc_init                         ! Passive tracers initialization 
    336347                           CALL dia_ptr_init   ! Poleward TRansports initialization 
    337348                            

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/PISCES/SED/trcdmp_sed.F90

    r10975 r11480  
    149149   !!---------------------------------------------------------------------- 
    150150CONTAINS 
    151    SUBROUTINE trc_dmp_sed( kt )        ! Empty routine 
     151   SUBROUTINE trc_dmp_sed( kt, Kbb, Kmm, Krhs )   ! Empty routine 
    152152      INTEGER, INTENT(in) :: kt 
     153      INTEGER, INTENT(in) :: Kbb, Kmm, Krhs 
    153154      WRITE(*,*) 'trc_dmp_sed: You should not have seen this print! error?', kt 
    154155   END SUBROUTINE trc_dmp_sed 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcrad.F90

    r10985 r11480  
    3737CONTAINS 
    3838 
    39    SUBROUTINE trc_rad( kt, Kbb, Kmm, Krhs, ptr ) 
     39   SUBROUTINE trc_rad( kt, Kbb, Kmm, ptr ) 
    4040      !!---------------------------------------------------------------------- 
    4141      !!                  ***  ROUTINE trc_rad  *** 
    … …  
    5252      !!                (the total CFC content is not strictly preserved) 
    5353      !!---------------------------------------------------------------------- 
    54       INTEGER,                                    INTENT(in   ) :: kt              ! ocean time-step index 
    55       INTEGER,                                    INTENT(in   ) :: Kbb, Kmm, Krhs  ! time level indices 
    56       REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr             ! passive tracers and RHS of tracer equation 
     54      INTEGER,                                    INTENT(in   ) :: kt         ! ocean time-step index 
     55      INTEGER,                                    INTENT(in   ) :: Kbb, Kmm   ! time level indices 
     56      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr        ! passive tracers and RHS of tracer equation 
    5757      ! 
    5858      CHARACTER (len=22) :: charout 
    … …  
    6161      IF( ln_timing )   CALL timing_start('trc_rad') 
    6262      ! 
    63       IF( ln_age     )   CALL trc_rad_sms( kt, Kmm, Krhs, ptr(:,:,:,:,Kbb), ptr(:,:,:,:,Kmm), jp_age , jp_age                )  !  AGE 
    64       IF( ll_cfc     )   CALL trc_rad_sms( kt, Kmm, Krhs, ptr(:,:,:,:,Kbb), ptr(:,:,:,:,Kmm), jp_cfc0, jp_cfc1               )  !  CFC model 
    65       IF( ln_c14     )   CALL trc_rad_sms( kt, Kmm, Krhs, ptr(:,:,:,:,Kbb), ptr(:,:,:,:,Kmm), jp_c14 , jp_c14                )  !  C14 
    66       IF( ln_pisces  )   CALL trc_rad_sms( kt, Kmm, Krhs, ptr(:,:,:,:,Kbb), ptr(:,:,:,:,Kmm), jp_pcs0, jp_pcs1, cpreserv='Y' )  !  PISCES model 
    67       IF( ln_my_trc  )   CALL trc_rad_sms( kt, Kmm, Krhs, ptr(:,:,:,:,Kbb), ptr(:,:,:,:,Kmm), jp_myt0, jp_myt1               )  !  MY_TRC model 
     63      IF( ln_age     )   CALL trc_rad_sms( kt, Kbb, Kmm, ptr, jp_age , jp_age                )  !  AGE 
     64      IF( ll_cfc     )   CALL trc_rad_sms( kt, Kbb, Kmm, ptr, jp_cfc0, jp_cfc1               )  !  CFC model 
     65      IF( ln_c14     )   CALL trc_rad_sms( kt, Kbb, Kmm, ptr, jp_c14 , jp_c14                )  !  C14 
     66      IF( ln_pisces  )   CALL trc_rad_sms( kt, Kbb, Kmm, ptr, jp_pcs0, jp_pcs1, cpreserv='Y' )  !  PISCES model 
     67      IF( ln_my_trc  )   CALL trc_rad_sms( kt, Kbb, Kmm, ptr, jp_myt0, jp_myt1               )  !  MY_TRC model 
    6868      ! 
    6969      IF(ln_ctl) THEN      ! print mean trends (used for debugging) 
    7070         WRITE(charout, FMT="('rad')") 
    7171         CALL prt_ctl_trc_info( charout ) 
    72          CALL prt_ctl_trc( tab4d=ptr(:,:,:,:,Kmm), mask=tmask, clinfo=ctrcnm ) 
     72         CALL prt_ctl_trc( tab4d=ptr(:,:,:,:,Kbb), mask=tmask, clinfo=ctrcnm ) 
    7373      ENDIF 
    7474      ! 
    … …  
    115115 
    116116 
    117    SUBROUTINE trc_rad_sms( kt, Kmm, Krhs, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv ) 
    118       !!----------------------------------------------------------------------------- 
    119       !!                  ***  ROUTINE trc_rad_sms  *** 
    120       !! 
    121       !! ** Purpose :   "crappy" routine to correct artificial negative 
    122       !!              concentrations due to isopycnal scheme 
    123       !! 
    124       !! ** Method  : 2 cases : 
    125       !!                - Set negative concentrations to zero while computing 
    126       !!                  the corresponding tracer content that is added to the 
    127       !!                  tracers. Then, adjust the tracer concentration using 
    128       !!                  a multiplicative factor so that the total tracer  
    129       !!                  concentration is preserved. 
    130       !!                - simply set to zero the negative CFC concentration 
    131       !!                  (the total content of concentration is not strictly preserved) 
    132       !!-------------------------------------------------------------------------------- 
    133       INTEGER                                , INTENT(in   ) ::   kt                 ! ocean time-step index 
    134       INTEGER                                , INTENT(in   ) ::   Kmm, Krhs          ! time level indices 
    135       INTEGER                                , INTENT(in   ) ::   jp_sms0, jp_sms1   ! First & last index of the passive tracer model 
    136       REAL(wp), DIMENSION (jpi,jpj,jpk,jptra), INTENT(inout) ::   ptrb    , ptrn     ! before and now traceur concentration 
    137       CHARACTER( len = 1), OPTIONAL          , INTENT(in   ) ::   cpreserv           ! flag to preserve content or not 
    138       ! 
    139       INTEGER ::   ji, ji2, jj, jj2, jk, jn     ! dummy loop indices 
    140       INTEGER ::   icnt 
    141       LOGICAL ::   lldebug = .FALSE.            ! local logical 
    142       REAL(wp)::   zcoef, zs2rdt, ztotmass 
    143       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrneg, ztrpos 
    144       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrtrd   ! workspace arrays 
    145       !!---------------------------------------------------------------------- 
    146       ! 
    147       IF( l_trdtrc )   ALLOCATE( ztrtrd(jpi,jpj,jpk) ) 
    148       zs2rdt = 1. / ( 2. * rdt * REAL( nn_dttrc, wp ) ) 
    149       ! 
    150       IF( PRESENT( cpreserv )  ) THEN     !==  total tracer concentration is preserved  ==! 
    151          ! 
    152          ALLOCATE( ztrneg(1:jpi,1:jpj,jp_sms0:jp_sms1), ztrpos(1:jpi,1:jpj,jp_sms0:jp_sms1) ) 
    153  
    154          DO jn = jp_sms0, jp_sms1 
    155             ztrneg(:,:,jn) = SUM( MIN( 0., ptrb(:,:,:,jn) ) * cvol(:,:,:), dim = 3 )   ! sum of the negative values 
    156             ztrpos(:,:,jn) = SUM( MAX( 0., ptrb(:,:,:,jn) ) * cvol(:,:,:), dim = 3 )   ! sum of the positive values 
    157          END DO 
    158          CALL sum3x3( ztrneg ) 
    159          CALL sum3x3( ztrpos ) 
    160           
    161          DO jn = jp_sms0, jp_sms1 
    162             ! 
    163             IF( l_trdtrc )   ztrtrd(:,:,:) = ptrb(:,:,:,jn)                            ! save input tr(:,:,:,:,Kbb) for trend computation            
    164             ! 
    165             DO jk = 1, jpkm1 
    166                DO jj = 1, jpj 
    167                   DO ji = 1, jpi 
    168                      IF( ztrneg(ji,jj,jn) /= 0. ) THEN                                 ! if negative values over the 3x3 box 
    169                         ! 
    170                         ptrb(ji,jj,jk,jn) = ptrb(ji,jj,jk,jn) * tmask(ji,jj,jk)   ! really needed? 
    171                         IF( ptrb(ji,jj,jk,jn) < 0. ) ptrb(ji,jj,jk,jn) = 0.       ! supress negative values 
    172                         IF( ptrb(ji,jj,jk,jn) > 0. ) THEN                         ! use positive values to compensate mass gain 
    173                            zcoef = 1. + ztrneg(ji,jj,jn) / ztrpos(ji,jj,jn)       ! ztrpos > 0 as ptrb > 0 
    174                            ptrb(ji,jj,jk,jn) = ptrb(ji,jj,jk,jn) * zcoef 
    175                            IF( zcoef < 0. ) THEN                                  ! if the compensation exceed the positive value 
    176                               gainmass(jn,1) = gainmass(jn,1) - ptrb(ji,jj,jk,jn) * cvol(ji,jj,jk)   ! we are adding mass... 
    177                               ptrb(ji,jj,jk,jn) = 0.                              ! limit the compensation to keep positive value 
    178                            ENDIF 
    179                         ENDIF 
    180                         ! 
    181                      ENDIF 
    182                   END DO 
    183                END DO 
    184             END DO 
    185             ! 
    186             IF( l_trdtrc ) THEN 
    187                ztrtrd(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt 
    188                CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling 
    189             ENDIF 
    190             ! 
    191          END DO 
    192   
    193          IF( kt == nitend ) THEN 
    194             CALL mpp_sum( 'trcrad', gainmass(:,1) ) 
    195             DO jn = jp_sms0, jp_sms1 
    196                IF( gainmass(jn,1) > 0. ) THEN 
    197                   ztotmass = glob_sum( 'trcrad', ptrb(:,:,:,jn) * cvol(:,:,:) ) 
    198                   IF(lwp) WRITE(numout, '(a, i2, a, D23.16, a, D23.16)') 'trcrad ptrb, traceur ', jn  & 
    199                      &        , ' total mass : ', ztotmass, ', mass gain : ',  gainmass(jn,1) 
    200                END IF 
    201             END DO 
    202          ENDIF 
    203  
    204          DO jn = jp_sms0, jp_sms1 
    205             ztrneg(:,:,jn) = SUM( MIN( 0., ptrn(:,:,:,jn) ) * cvol(:,:,:), dim = 3 )   ! sum of the negative values 
    206             ztrpos(:,:,jn) = SUM( MAX( 0., ptrn(:,:,:,jn) ) * cvol(:,:,:), dim = 3 )   ! sum of the positive values 
    207          END DO 
    208          CALL sum3x3( ztrneg ) 
    209          CALL sum3x3( ztrpos ) 
    210           
    211          DO jn = jp_sms0, jp_sms1 
    212             ! 
    213             IF( l_trdtrc )   ztrtrd(:,:,:) = ptrn(:,:,:,jn)                            ! save input tr for trend computation 
    214             ! 
    215             DO jk = 1, jpkm1 
    216                DO jj = 1, jpj 
    217                   DO ji = 1, jpi 
    218                      IF( ztrneg(ji,jj,jn) /= 0. ) THEN                                 ! if negative values over the 3x3 box 
    219                         ! 
    220                         ptrn(ji,jj,jk,jn) = ptrn(ji,jj,jk,jn) * tmask(ji,jj,jk)   ! really needed? 
    221                         IF( ptrn(ji,jj,jk,jn) < 0. ) ptrn(ji,jj,jk,jn) = 0.       ! supress negative values 
    222                         IF( ptrn(ji,jj,jk,jn) > 0. ) THEN                         ! use positive values to compensate mass gain 
    223                            zcoef = 1. + ztrneg(ji,jj,jn) / ztrpos(ji,jj,jn)       ! ztrpos > 0 as ptrb > 0 
    224                            ptrn(ji,jj,jk,jn) = ptrn(ji,jj,jk,jn) * zcoef 
    225                            IF( zcoef < 0. ) THEN                                  ! if the compensation exceed the positive value 
    226                               gainmass(jn,2) = gainmass(jn,2) - ptrn(ji,jj,jk,jn) * cvol(ji,jj,jk)   ! we are adding mass... 
    227                               ptrn(ji,jj,jk,jn) = 0.                              ! limit the compensation to keep positive value 
    228                            ENDIF 
    229                         ENDIF 
    230                         ! 
    231                      ENDIF 
    232                   END DO 
    233                END DO 
    234             END DO 
    235             ! 
    236             IF( l_trdtrc ) THEN 
    237                ztrtrd(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt 
    238                CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling 
    239             ENDIF 
    240             ! 
    241          END DO 
    242   
    243          IF( kt == nitend ) THEN 
    244             CALL mpp_sum( 'trcrad', gainmass(:,2) ) 
    245             DO jn = jp_sms0, jp_sms1 
    246                IF( gainmass(jn,2) > 0. ) THEN 
    247                   ztotmass = glob_sum( 'trcrad', ptrn(:,:,:,jn) * cvol(:,:,:) ) 
    248                   WRITE(numout, '(a, i2, a, D23.16, a, D23.16)') 'trcrad ptrn, traceur ', jn  & 
    249                      &        , ' total mass : ', ztotmass, ', mass gain : ',  gainmass(jn,1) 
    250                END IF 
    251             END DO 
    252          ENDIF 
    253  
    254          DEALLOCATE( ztrneg, ztrpos ) 
    255          ! 
    256       ELSE                                !==  total CFC content is NOT strictly preserved  ==! 
    257          ! 
    258          DO jn = jp_sms0, jp_sms1   
    259             ! 
    260             IF( l_trdtrc )   ztrtrd(:,:,:) = ptrb(:,:,:,jn)                        ! save input tr for trend computation 
    261             ! 
    262             WHERE( ptrb(:,:,:,jn) < 0. )   ptrb(:,:,:,jn) = 0. 
    263             ! 
    264             IF( l_trdtrc ) THEN 
    265                ztrtrd(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt 
    266                CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling 
    267             ENDIF 
    268             ! 
    269             IF( l_trdtrc )   ztrtrd(:,:,:) = ptrn(:,:,:,jn)                        ! save input tr for trend computation 
    270             ! 
    271             WHERE( ptrn(:,:,:,jn) < 0. )   ptrn(:,:,:,jn) = 0. 
    272             ! 
    273             IF( l_trdtrc ) THEN 
    274                ztrtrd(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrd(:,:,:) ) * zs2rdt 
    275                CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_radn, ztrtrd )       ! standard     trend handling 
    276             ENDIF 
    277             ! 
    278          END DO 
    279          ! 
    280       ENDIF 
     117   SUBROUTINE trc_rad_sms( kt, Kbb, Kmm, ptr, jp_sms0, jp_sms1, cpreserv ) 
     118     !!----------------------------------------------------------------------------- 
     119     !!                  ***  ROUTINE trc_rad_sms  *** 
     120     !! 
     121     !! ** Purpose :   "crappy" routine to correct artificial negative 
     122     !!              concentrations due to isopycnal scheme 
     123     !! 
     124     !! ** Method  : 2 cases : 
     125     !!                - Set negative concentrations to zero while computing 
     126     !!                  the corresponding tracer content that is added to the 
     127     !!                  tracers. Then, adjust the tracer concentration using 
     128     !!                  a multiplicative factor so that the total tracer  
     129     !!                  concentration is preserved. 
     130     !!                - simply set to zero the negative CFC concentration 
     131     !!                  (the total content of concentration is not strictly preserved) 
     132     !!-------------------------------------------------------------------------------- 
     133     INTEGER                                    , INTENT(in   ) ::   kt                 ! ocean time-step index 
     134     INTEGER                                    , INTENT(in   ) ::   Kbb, Kmm           ! time level indices 
     135     INTEGER                                    , INTENT(in   ) ::   jp_sms0, jp_sms1   ! First & last index of the passive tracer model 
     136     REAL(wp), DIMENSION (jpi,jpj,jpk,jptra,jpt), INTENT(inout) ::   ptr                ! before and now traceur concentration 
     137     CHARACTER( len = 1), OPTIONAL              , INTENT(in   ) ::   cpreserv           ! flag to preserve content or not 
     138     ! 
     139     INTEGER ::   ji, ji2, jj, jj2, jk, jn, jt ! dummy loop indices 
     140     INTEGER ::   icnt, itime 
     141     LOGICAL ::   lldebug = .FALSE.            ! local logical 
     142     REAL(wp)::   zcoef, zs2rdt, ztotmass 
     143     REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrneg, ztrpos 
     144     REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrtrd   ! workspace arrays 
     145     !!---------------------------------------------------------------------- 
     146     ! 
     147     IF( l_trdtrc )   ALLOCATE( ztrtrd(jpi,jpj,jpk) ) 
     148     zs2rdt = 1. / ( 2. * rdt * REAL( nn_dttrc, wp ) ) 
     149     ! 
     150     DO jt = 1,2  ! Loop over time indices since exactly the same fix is applied to "now" and "after" fields 
     151        IF( jt == 1 ) itime = Kbb 
     152        IF( jt == 2 ) itime = Kmm 
     153 
     154        IF( PRESENT( cpreserv )  ) THEN     !==  total tracer concentration is preserved  ==! 
     155           ! 
     156           ALLOCATE( ztrneg(1:jpi,1:jpj,jp_sms0:jp_sms1), ztrpos(1:jpi,1:jpj,jp_sms0:jp_sms1) ) 
     157 
     158           DO jn = jp_sms0, jp_sms1 
     159              ztrneg(:,:,jn) = SUM( MIN( 0., ptr(:,:,:,jn,itime) ) * cvol(:,:,:), dim = 3 )   ! sum of the negative values 
     160              ztrpos(:,:,jn) = SUM( MAX( 0., ptr(:,:,:,jn,itime) ) * cvol(:,:,:), dim = 3 )   ! sum of the positive values 
     161           END DO 
     162           CALL sum3x3( ztrneg ) 
     163           CALL sum3x3( ztrpos ) 
     164 
     165           DO jn = jp_sms0, jp_sms1 
     166              ! 
     167              IF( l_trdtrc )   ztrtrd(:,:,:) = ptr(:,:,:,jn,itime)                       ! save input tr(:,:,:,:,Kbb) for trend computation            
     168              ! 
     169              DO jk = 1, jpkm1 
     170                 DO jj = 1, jpj 
     171                    DO ji = 1, jpi 
     172                       IF( ztrneg(ji,jj,jn) /= 0. ) THEN                                 ! if negative values over the 3x3 box 
     173                          ! 
     174                          ptr(ji,jj,jk,jn,itime) = ptr(ji,jj,jk,jn,itime) * tmask(ji,jj,jk)   ! really needed? 
     175                          IF( ptr(ji,jj,jk,jn,itime) < 0. ) ptr(ji,jj,jk,jn,itime) = 0.       ! suppress negative values 
     176                          IF( ptr(ji,jj,jk,jn,itime) > 0. ) THEN                    ! use positive values to compensate mass gain 
     177                             zcoef = 1. + ztrneg(ji,jj,jn) / ztrpos(ji,jj,jn)       ! ztrpos > 0 as ptr > 0 
     178                             ptr(ji,jj,jk,jn,itime) = ptr(ji,jj,jk,jn,itime) * zcoef 
     179                             IF( zcoef < 0. ) THEN                                  ! if the compensation exceed the positive value 
     180                                gainmass(jn,1) = gainmass(jn,1) - ptr(ji,jj,jk,jn,itime) * cvol(ji,jj,jk)   ! we are adding mass... 
     181                                ptr(ji,jj,jk,jn,itime) = 0.                         ! limit the compensation to keep positive value 
     182                             ENDIF 
     183                          ENDIF 
     184                          ! 
     185                       ENDIF 
     186                    END DO 
     187                 END DO 
     188              END DO 
     189              ! 
     190              IF( l_trdtrc ) THEN 
     191                 ztrtrd(:,:,:) = ( ptr(:,:,:,jn,itime) - ztrtrd(:,:,:) ) * zs2rdt 
     192                 CALL trd_tra( kt, Kbb, Kmm, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling 
     193              ENDIF 
     194              ! 
     195           END DO 
     196 
     197           IF( kt == nitend ) THEN 
     198              CALL mpp_sum( 'trcrad', gainmass(:,1) ) 
     199              DO jn = jp_sms0, jp_sms1 
     200                 IF( gainmass(jn,1) > 0. ) THEN 
     201                    ztotmass = glob_sum( 'trcrad', ptr(:,:,:,jn,itime) * cvol(:,:,:) ) 
     202                    IF(lwp) WRITE(numout, '(a, i2, a, D23.16, a, D23.16)') 'trcrad ptrb, traceur ', jn  & 
     203                         &        , ' total mass : ', ztotmass, ', mass gain : ',  gainmass(jn,1) 
     204                 END IF 
     205              END DO 
     206           ENDIF 
     207 
     208           DEALLOCATE( ztrneg, ztrpos ) 
     209           ! 
     210        ELSE                                !==  total CFC content is NOT strictly preserved  ==! 
     211           ! 
     212           DO jn = jp_sms0, jp_sms1   
     213              ! 
     214              IF( l_trdtrc )   ztrtrd(:,:,:) = ptr(:,:,:,jn,itime)                 ! save input tr for trend computation 
     215              ! 
     216              WHERE( ptr(:,:,:,jn,itime) < 0. )   ptr(:,:,:,jn,itime) = 0. 
     217              ! 
     218              IF( l_trdtrc ) THEN 
     219                 ztrtrd(:,:,:) = ( ptr(:,:,:,jn,itime) - ztrtrd(:,:,:) ) * zs2rdt 
     220                 CALL trd_tra( kt, Kbb, Kmm, 'TRC', jn, jptra_radb, ztrtrd )       ! Asselin-like trend handling 
     221              ENDIF 
     222              ! 
     223           END DO 
     224           ! 
     225        ENDIF 
     226        ! 
     227      END DO 
    281228      ! 
    282229      IF( l_trdtrc )  DEALLOCATE( ztrtrd ) 
    … …  
    289236   !!---------------------------------------------------------------------- 
    290237CONTAINS 
    291    SUBROUTINE trc_rad( kt, Kbb, Kmm, Krhs )              ! Empty routine 
     238   SUBROUTINE trc_rad( kt, Kbb, Kmm )              ! Empty routine 
    292239      INTEGER, INTENT(in) ::   kt 
    293       INTEGER, INTENT(in) ::   Kbb, Kmm, Krhs  ! time level indices 
     240      INTEGER, INTENT(in) ::   Kbb, Kmm  ! time level indices 
    294241      WRITE(*,*) 'trc_rad: You should not have seen this print! error?', kt 
    295242   END SUBROUTINE trc_rad 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcsbc.F90

    r10985 r11480  
    197197   !!   Dummy module :                      NO passive tracer 
    198198   !!---------------------------------------------------------------------- 
     199   USE par_oce 
     200   USE par_trc 
    199201CONTAINS 
    200202   SUBROUTINE trc_sbc ( kt, Kmm, ptr, Krhs )      ! Empty routine 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90

    r10985 r11480  
    2020   USE trcadv          ! advection                           (trc_adv routine) 
    2121   USE trczdf          ! vertical diffusion                  (trc_zdf routine) 
    22    USE trcnxt          ! time-stepping                       (trc_nxt routine) 
     22   USE trcatf          ! time filtering                      (trc_atf routine) 
    2323   USE trcrad          ! positivity                          (trc_rad routine) 
    2424   USE trcsbc          ! surface boundary condition          (trc_sbc routine) 
    … …  
    5353      !!              - Update the passive tracers 
    5454      !!---------------------------------------------------------------------- 
    55       INTEGER, INTENT( in ) :: kt                  ! ocean time-step index 
    56       INTEGER, INTENT( in ) :: Kbb, Kmm, Krhs, Kaa ! time level indices 
     55      INTEGER, INTENT( in    ) :: kt                  ! ocean time-step index 
     56      INTEGER, INTENT( inout ) :: Kbb, Kmm, Krhs, Kaa ! TOP time level indices (swapped in this routine) 
    5757      !! --------------------------------------------------------------------- 
    5858      ! 
    … …  
    7878#endif 
    7979                                CALL trc_zdf    ( kt, Kbb, Kmm, Krhs, tr, Kaa  )  ! vert. mixing & after tracer   ==> after 
    80                                 CALL trc_nxt    ( kt, Kbb, Kmm, Krhs )            ! tracer fields at next time step      
    81          IF( ln_trcrad )        CALL trc_rad    ( kt, Kbb, Kmm, Krhs, tr       )  ! Correct artificial negative concentrations 
    82          IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kt, Kbb, Kmm )                  ! internal damping trends on closed seas only 
     80                                CALL trc_atf    ( kt, Kbb, Kmm, Kaa , tr )        ! time filtering of "now" tracer fields     
     81         ! 
     82         ! Swap TOP time levels (= Nrhs_trc, Nbb_trc etc) 
     83         Krhs = Kbb 
     84         Kbb = Kmm 
     85         Kmm = Kaa 
     86         Kaa = Krhs 
     87         ! 
     88         IF( ln_trcrad )        CALL trc_rad    ( kt, Kbb, Kmm, tr       )    ! Correct artificial negative concentrations 
     89         IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kt, Kbb, Kmm )              ! internal damping trends on closed seas only 
    8390 
    8491         ! 
    … …  
    8794         IF( ln_trcdmp )        CALL trc_dmp( kt, Kbb, Kmm,       tr, Krhs )  ! internal damping trends 
    8895                                CALL trc_zdf( kt, Kbb, Kmm, Krhs, tr, Kaa  )  ! vert. mixing & after tracer ==> after 
    89                                 CALL trc_nxt( kt, Kbb, Kmm, Krhs )            ! tracer fields at next time step      
    90           IF( ln_trcrad )       CALL trc_rad( kt, Kbb, Kmm, Krhs, tr       )  ! Correct artificial negative concentrations 
     96                                CALL trc_atf( kt, Kbb, Kmm, Kaa , tr )        ! time filtering of "now" tracer fields 
     97         ! 
     98         ! Swap TOP time levels (= Nrhs_trc, Nbb_trc etc) 
     99         Krhs = Kbb 
     100         Kbb = Kmm 
     101         Kmm = Kaa 
     102         Kaa = Krhs 
     103         ! 
     104         IF( ln_trcrad )       CALL trc_rad( kt, Kbb, Kmm, tr       )  ! Correct artificial negative concentrations 
    91105         ! 
    92106      END IF 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trcini.F90

    r11027 r11480  
    4040CONTAINS 
    4141    
    42    SUBROUTINE trc_init( Kbb, Kmm, Kaa ) 
     42   SUBROUTINE trc_init 
    4343      !!--------------------------------------------------------------------- 
    4444      !!                     ***  ROUTINE trc_init  *** 
    … …  
    5252      !!                or read data or analytical formulation 
    5353      !!--------------------------------------------------------------------- 
    54       INTEGER, INTENT( in ) :: Kbb, Kmm, Kaa ! time level indices 
    5554      ! 
    5655      IF( ln_timing )   CALL timing_start('trc_init') 
    … …  
    5958      IF(lwp) WRITE(numout,*) 'trc_init : initial set up of the passive tracers' 
    6059      IF(lwp) WRITE(numout,*) '~~~~~~~~' 
     60      ! 
     61      ! Initialise time level indices 
     62      Nbb_trc = 1; Nnn_trc = 2; Naa_trc = 3; Nrhs_trc = Naa_trc 
    6163      ! 
    6264      CALL trc_ini_ctl   ! control  
    … …  
    7173      IF(lwp) WRITE(numout,*) 
    7274      ! 
    73       CALL trc_ini_sms( Kmm )   ! SMS 
     75      CALL trc_ini_sms( Nnn_trc )   ! SMS 
    7476      CALL trc_ini_trp          ! passive tracers transport 
    7577      CALL trc_ice_ini          ! Tracers in sea ice 
    … …  
    7981      ENDIF 
    8082      ! 
    81       CALL trc_ini_state( Kbb, Kmm, Kaa )  !  passive tracers initialisation : from a restart or from clim 
     83      CALL trc_ini_state( Nbb_trc, Nnn_trc, Naa_trc )  !  passive tracers initialisation : from a restart or from clim 
    8284      IF( nn_dttrc /= 1 ) & 
    83       CALL trc_sub_ini( Kmm )    ! Initialize variables for substepping passive tracers 
    84       ! 
    85       CALL trc_ini_inv( Kmm )    ! Inventories 
     85      CALL trc_sub_ini( Nnn_trc )    ! Initialize variables for substepping passive tracers 
     86      ! 
     87      CALL trc_ini_inv( Nnn_trc )    ! Inventories 
    8688      ! 
    8789      IF( ln_timing )   CALL timing_stop('trc_init') 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trcstp.F90

    r11027 r11480  
    3838 
    3939   !!---------------------------------------------------------------------- 
     40   !! time level indices 
     41   !!---------------------------------------------------------------------- 
     42   INTEGER, PUBLIC :: Nbb_trc, Nnn_trc, Naa_trc, Nrhs_trc      !! used by trc_init 
     43 
     44   !!---------------------------------------------------------------------- 
    4045   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
    4146   !! $Id$  
    … …  
    4449CONTAINS 
    4550 
    46    SUBROUTINE trc_stp( kt, Kbb, Kmm, Krhs, Kaa ) 
     51   SUBROUTINE trc_stp( kt, Kbb_oce, Kmm_oce, Krhs_oce, Kaa_oce ) 
    4752      !!------------------------------------------------------------------- 
    4853      !!                     ***  ROUTINE trc_stp  *** 
    … …  
    5459      !!------------------------------------------------------------------- 
    5560      INTEGER, INTENT( in ) :: kt                  ! ocean time-step index 
    56       INTEGER, INTENT( in ) :: Kbb, Kmm, Krhs, Kaa ! time level indices 
     61      INTEGER, INTENT( in ) :: Kbb_oce, Kmm_oce, Krhs_oce, Kaa_oce ! time level indices 
    5762      ! 
    5863      INTEGER ::   jk, jn   ! dummy loop indices 
    … …  
    7681      IF( .NOT.ln_linssh ) THEN                                           ! update ocean volume due to ssh temporal evolution 
    7782         DO jk = 1, jpk 
    78             cvol(:,:,jk) = e1e2t(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     83            cvol(:,:,jk) = e1e2t(:,:) * e3t(:,:,jk,Kmm_oce) * tmask(:,:,jk) 
    7984         END DO 
    8085         IF ( ln_ctl .OR. kt == nitrst .OR. ( ln_check_mass .AND. kt == nitend )              & 
    … …  
    8691      IF( l_trcdm2dc )   CALL trc_mean_qsr( kt ) 
    8792      !     
    88       IF( nn_dttrc /= 1 )   CALL trc_sub_stp( kt, Kbb, Kmm, Krhs )  ! averaging physical variables for sub-stepping 
     93      IF( nn_dttrc == 1 )  THEN 
     94         IF(lwp) WRITE(numout,*) "Kbb_oce, Kmm_oce, Kaa_oce, Krhs_oce : ",Kbb_oce, Kmm_oce, Kaa_oce, Krhs_oce 
     95         IF(lwp) WRITE(numout,*) "Nbb_trc, Nnn_trc, Naa_trc, Nrhs_trc : ",Nbb_trc, Nnn_trc, Naa_trc, Nrhs_trc 
     96         IF(lwp) CALL FLUSH(numout) 
     97         CALL mppsync()       
     98         IF( Kmm_oce /= Nnn_trc .OR. Kaa_oce /= Naa_trc .OR. Krhs_oce /= Nrhs_trc ) THEN 
     99            ! The nn_dttrc == 1 case depends on the OCE and TRC time indices being the same always.  
     100            ! If this is not the case then something has gone wrong. 
     101            CALL ctl_stop( 'trc_stp : nn_dttrc = 1 but OCE and TRC time indices are different! Something has gone wrong.' ) 
     102         ENDIF 
     103      ELSE 
     104         CALL trc_sub_stp( kt, Nbb_trc, Nnn_trc, Nrhs_trc )  ! averaging physical variables for sub-stepping 
     105      ENDIF 
    89106      !     
    90107      IF( MOD( kt , nn_dttrc ) == 0 ) THEN      ! only every nn_dttrc time step 
    … …  
    95112         ENDIF 
    96113         ! 
    97          tr(:,:,:,:,Krhs) = 0.e0 
     114         tr(:,:,:,:,Nrhs_trc) = 0.e0 
    98115         ! 
    99116                                   CALL trc_rst_opn  ( kt )       ! Open tracer restart file  
    100117         IF( lrst_trc )            CALL trc_rst_cal  ( kt, 'WRITE' )   ! calendar 
    101                                    CALL trc_wri      ( kt,      Kmm       )       ! output of passive tracers with iom I/O manager 
    102                                    CALL trc_sms      ( kt, Kbb, Kmm, Krhs )       ! tracers: sinks and sources 
    103                                    CALL trc_trp      ( kt, Kbb, Kmm, Krhs, Kaa )  ! transport of passive tracers 
     118                                   CALL trc_wri      ( kt,          Nnn_trc                    )  ! output of passive tracers with iom I/O manager 
     119                                   CALL trc_sms      ( kt, Nbb_trc, Nnn_trc, Nrhs_trc          )  ! tracers: sinks and sources 
     120                                   CALL trc_trp      ( kt, Nbb_trc, Nnn_trc, Nrhs_trc, Naa_trc )  ! transport of passive tracers 
    104121         IF( kt == nittrc000 ) THEN 
    105122            CALL iom_close( numrtr )       ! close input tracer restart file 
    106123            IF(lwm) CALL FLUSH( numont )   ! flush namelist output 
    107124         ENDIF 
    108          IF( lrst_trc )            CALL trc_rst_wri  ( kt, Kbb, Kmm, Krhs )       ! write tracer restart file 
    109          IF( lk_trdmxl_trc  )      CALL trd_mxl_trc  ( kt,      Kmm       )       ! trends: Mixed-layer 
    110          ! 
    111          IF( nn_dttrc /= 1   )     CALL trc_sub_reset( kt, Kbb, Kmm, Krhs )       ! resetting physical variables when sub-stepping 
     125         IF( lrst_trc )            CALL trc_rst_wri  ( kt, Nbb_trc, Nnn_trc, Nrhs_trc )       ! write tracer restart file 
     126         IF( lk_trdmxl_trc  )      CALL trd_mxl_trc  ( kt,          Nnn_trc           )       ! trends: Mixed-layer 
     127         ! 
     128         IF( nn_dttrc /= 1   )     CALL trc_sub_reset( kt, Nbb_trc, Nnn_trc, Nrhs_trc )       ! resetting physical variables when sub-stepping 
    112129         ! 
    113130      ENDIF 
    … …  
    116133         ztrai = 0._wp                                                   !  content of all tracers 
    117134         DO jn = 1, jptra 
    118             ztrai = ztrai + glob_sum( 'trcstp', tr(:,:,:,jn,Kmm) * cvol(:,:,:)   ) 
     135            ztrai = ztrai + glob_sum( 'trcstp', tr(:,:,:,jn,Nnn_trc) * cvol(:,:,:)   ) 
    119136         END DO 
    120137         IF( lwm ) WRITE(numstr,9300) kt,  ztrai / areatot 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trcsub.F90

    r10969 r11480  
    9292      !!------------------------------------------------------------------- 
    9393      INTEGER, INTENT( in ) ::   kt              ! ocean time-step index 
    94       INTEGER, INTENT( in ) ::   Kbb, Kmm, Krhs  ! ocean time-level index 
     94      INTEGER, INTENT( in ) ::   Kbb, Kmm, Krhs  !       time-level index 
    9595      ! 
    9696      INTEGER ::   ji, jj, jk   ! dummy loop indices 

  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/tests/VORTEX/MY_SRC/domvvl.F90

    r11412 r11480  
    1313   !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
    1414   !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
    15    !!   dom_vvl_sf_swp   : Swap vertical scale factors and update the vertical grid 
     15   !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
    1616   !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
    1717   !!   dom_vvl_rst      : read/write restart file 
    … …  
    3737   PUBLIC  dom_vvl_init       ! called by domain.F90 
    3838   PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
    39    PUBLIC  dom_vvl_sf_swp     ! called by step.F90 
     39   PUBLIC  dom_vvl_sf_update  ! called by step.F90 
    4040   PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
    4141 
    … …  
    181181      ! 
    182182      !                    !==  thickness of the water column  !!   (ocean portion only) 
    183       ht_n(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
    184       hu_b(:,:) = e3u(:,:,1,Kbb) * umask(:,:,1) 
    185       hu_n(:,:) = e3u(:,:,1,Kmm) * umask(:,:,1) 
    186       hv_b(:,:) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
    187       hv_n(:,:) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
     183      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
     184      hu(:,:,Kbb) = e3u(:,:,1,Kbb) * umask(:,:,1) 
     185      hu(:,:,Kmm) = e3u(:,:,1,Kmm) * umask(:,:,1) 
     186      hv(:,:,Kbb) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
     187      hv(:,:,Kmm) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
    188188      DO jk = 2, jpkm1 
    189          ht_n(:,:) = ht_n(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    190          hu_b(:,:) = hu_b(:,:) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
    191          hu_n(:,:) = hu_n(:,:) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
    192          hv_b(:,:) = hv_b(:,:) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
    193          hv_n(:,:) = hv_n(:,:) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
     189         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     190         hu(:,:,Kbb) = hu(:,:,Kbb) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
     191         hu(:,:,Kmm) = hu(:,:,Kmm) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
     192         hv(:,:,Kbb) = hv(:,:,Kbb) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
     193         hv(:,:,Kmm) = hv(:,:,Kmm) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    194194      END DO 
    195195      ! 
    196196      !                    !==  inverse of water column thickness   ==!   (u- and v- points) 
    197       r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
    198       r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
    199       r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
    200       r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     197      r1_hu(:,:,Kbb) = ssumask(:,:) / ( hu(:,:,Kbb) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
     198      r1_hu(:,:,Kmm) = ssumask(:,:) / ( hu(:,:,Kmm) + 1._wp - ssumask(:,:) ) 
     199      r1_hv(:,:,Kbb) = ssvmask(:,:) / ( hv(:,:,Kbb) + 1._wp - ssvmask(:,:) ) 
     200      r1_hv(:,:,Kmm) = ssvmask(:,:) / ( hv(:,:,Kmm) + 1._wp - ssvmask(:,:) ) 
    201201 
    202202      !                    !==   z_tilde coordinate case  ==!   (Restoring frequencies) 
    … …  
    550550      ! *********************************** ! 
    551551 
    552       hu_a(:,:) = e3u(:,:,1,Kaa) * umask(:,:,1) 
    553       hv_a(:,:) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
     552      hu(:,:,Kaa) = e3u(:,:,1,Kaa) * umask(:,:,1) 
     553      hv(:,:,Kaa) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
    554554      DO jk = 2, jpkm1 
    555          hu_a(:,:) = hu_a(:,:) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
    556          hv_a(:,:) = hv_a(:,:) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
     555         hu(:,:,Kaa) = hu(:,:,Kaa) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     556         hv(:,:,Kaa) = hv(:,:,Kaa) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    557557      END DO 
    558558      !                                        ! Inverse of the local depth 
    559559!!gm BUG ?  don't understand the use of umask_i here ..... 
    560       r1_hu_a(:,:) = ssumask(:,:) / ( hu_a(:,:) + 1._wp - ssumask(:,:) ) 
    561       r1_hv_a(:,:) = ssvmask(:,:) / ( hv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     560      r1_hu(:,:,Kaa) = ssumask(:,:) / ( hu(:,:,Kaa) + 1._wp - ssumask(:,:) ) 
     561      r1_hv(:,:,Kaa) = ssvmask(:,:) / ( hv(:,:,Kaa) + 1._wp - ssvmask(:,:) ) 
    562562      ! 
    563563      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_nxt') 
    … …  
    566566 
    567567 
    568    SUBROUTINE dom_vvl_sf_swp( kt, Kbb, Kmm, Kaa ) 
    569       !!---------------------------------------------------------------------- 
    570       !!                ***  ROUTINE dom_vvl_sf_swp  *** 
     568   SUBROUTINE dom_vvl_sf_update( kt, Kbb, Kmm, Kaa ) 
     569      !!---------------------------------------------------------------------- 
     570      !!                ***  ROUTINE dom_vvl_sf_update  *** 
    571571      !!                    
    572       !! ** Purpose :  compute time filter and swap of scale factors  
     572      !! ** Purpose :  for z tilde case: compute time filter and swap of scale factors  
    573573      !!               compute all depths and related variables for next time step 
    574574      !!               write outputs and restart file 
    575575      !! 
    576       !! ** Method  :  - swap of e3t with trick for volume/tracer conservation 
     576      !! ** Method  :  - swap of e3t with trick for volume/tracer conservation (ONLY FOR Z TILDE CASE) 
    577577      !!               - reconstruct scale factor at other grid points (interpolate) 
    578578      !!               - recompute depths and water height fields 
    579579      !! 
    580       !! ** Action  :  - e3t_(b/n), tilde_e3t_(b/n) and e3(u/v)_n ready for next time step 
     580      !! ** Action  :  - tilde_e3t_(b/n) ready for next time step 
    581581      !!               - Recompute: 
    582582      !!                    e3(u/v)_b        
    … …  
    599599      IF( ln_linssh )   RETURN      ! No calculation in linear free surface 
    600600      ! 
    601       IF( ln_timing )   CALL timing_start('dom_vvl_sf_swp') 
     601      IF( ln_timing )   CALL timing_start('dom_vvl_sf_update') 
    602602      ! 
    603603      IF( kt == nit000 )   THEN 
    604604         IF(lwp) WRITE(numout,*) 
    605          IF(lwp) WRITE(numout,*) 'dom_vvl_sf_swp : - time filter and swap of scale factors' 
    606          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~   - interpolate scale factors and compute depths for next time step' 
     605         IF(lwp) WRITE(numout,*) 'dom_vvl_sf_update : - interpolate scale factors and compute depths for next time step' 
     606         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~' 
    607607      ENDIF 
    608608      ! 
    … …  
    619619         tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) 
    620620      ENDIF 
    621       gdept(:,:,:,Kbb) = gdept(:,:,:,Kmm) 
    622       gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm) 
    623  
    624       e3t(:,:,:,Kmm) = e3t(:,:,:,Kaa) 
    625       e3u(:,:,:,Kmm) = e3u(:,:,:,Kaa) 
    626       e3v(:,:,:,Kmm) = e3v(:,:,:,Kaa) 
    627621 
    628622      ! Compute all missing vertical scale factor and depths 
    … …  
    631625      ! -------------------------------------- 
    632626      ! - ML - e3u(:,:,:,Kbb) and e3v(:,:,:,Kbb) are allready computed in dynnxt 
    633       ! - JC - hu_b, hv_b, hur_b, hvr_b also 
     627      ! - JC - hu(:,:,:,Kbb), hv(:,:,:,:,Kbb), hur_b, hvr_b also 
    634628       
    635629      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F'  ) 
    … …  
    663657      ! Local depth and Inverse of the local depth of the water 
    664658      ! ------------------------------------------------------- 
    665       hu_n(:,:) = hu_a(:,:)   ;   r1_hu_n(:,:) = r1_hu_a(:,:) 
    666       hv_n(:,:) = hv_a(:,:)   ;   r1_hv_n(:,:) = r1_hv_a(:,:) 
    667       ! 
    668       ht_n(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
     659      ! 
     660      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
    669661      DO jk = 2, jpkm1 
    670          ht_n(:,:) = ht_n(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     662         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    671663      END DO 
    672664 
    … …  
    675667      IF( lrst_oce  )   CALL dom_vvl_rst( kt, Kbb, Kmm, 'WRITE' ) 
    676668      ! 
    677       IF( ln_timing )   CALL timing_stop('dom_vvl_sf_swp') 
    678       ! 
    679    END SUBROUTINE dom_vvl_sf_swp 
     669      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_update') 
     670      ! 
     671   END SUBROUTINE dom_vvl_sf_update 
    680672 
    681673 
    … …  
    931923            ELSE 
    932924               ! 
    933 !               ! 
    934                ! usr_def_istate called here only to get ssh, that is needed to initialize e3t_b and e3t_n 
    935                CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  ) 
    936                ! usr_def_istate will be called again in istate_init to initialize ts(bn), ssh(bn), u(bn) and v(bn) 
     925               ! usr_def_istate called here only to get ssh(Kbb) needed to initialize e3t(Kbb) and e3t(Kmm) 
     926               ! 
     927               CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )   
     928               ! 
     929               ! usr_def_istate will be called again in istate_init to initialize ts, ssh, u and v 
    937930               ! 
    938931               DO jk=1,jpk 
    939                   e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) )           & 
    940                     &               / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    941                     &               + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )             ! make sure e3t(Kbb) != 0 on land points 
     932                  e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) ) & 
     933                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
     934                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )   ! make sure e3t(:,:,:,Kbb) != 0 on land points 
    942935               END DO 
    943936               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
    944                ssh(:,:  ,Kmm) = ssh(:,:  ,Kbb) ! needed later for gde3w 
     937               ssh(:,:,Kmm) = ssh(:,:,Kbb)                                    ! needed later for gde3w 
    945938               ! 
    946939            END IF           ! end of ll_wd edits 
Note: See TracChangeset for help on using the changeset viewer.