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 11099 for NEMO/branches/2019 – NEMO

Changeset 11099 for NEMO/branches/2019


Ignore:
Timestamp:
2019-06-11T15:59:58+02:00 (5 years ago)
Author:
davestorkey
Message:

dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap : Updates and bug fixes. This version still doesn't pass all SETTE tests.

Location:
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/cfgs/SHARED/field_def_nemo-oce.xml

    r11057 r11099  
    6262        <field id="mldr10_1max"  long_name="Max of Mixed Layer Depth (dsigma = 0.01 wrt 10m)"   field_ref="mldr10_1"   operation="maximum"                                                                          /> 
    6363        <field id="mldr10_1min"  long_name="Min of Mixed Layer Depth (dsigma = 0.01 wrt 10m)"   field_ref="mldr10_1"   operation="minimum"                                                                          /> 
     64         <field id="mldzint_1"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
     65         <field id="mldzint_2"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
     66         <field id="mldzint_3"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
     67         <field id="mldzint_4"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
     68         <field id="mldzint_5"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
     69         <field id="mldhtc_1"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
     70         <field id="mldhtc_2"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
     71         <field id="mldhtc_3"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
     72         <field id="mldhtc_4"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
     73         <field id="mldhtc_5"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    6474        <field id="heatc"        long_name="Heat content vertically integrated"                 standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    6575        <field id="saltc"        long_name="Salt content vertically integrated"                                                                                                                   unit="1e-3*kg/m2" /> 
     
    7787        <!-- variables available with MLE --> 
    7888        <field id="Lf_NHpf"      long_name="MLE: Lf = N H / f"   unit="m" /> 
    79  
    80         <!-- next variables available with ln_zad_Aimp=.true. --> 
    81         <field id="Courant"      long_name="Courant number"                                                                 unit="#"   grid_ref="grid_T_3D" /> 
    82         <field id="wimp"         long_name="Implicit vertical velocity"                                                     unit="m/s" grid_ref="grid_T_3D" /> 
    83         <field id="wexp"         long_name="Explicit vertical velocity"                                                     unit="m/s" grid_ref="grid_T_3D" /> 
    84         <field id="wi_cff"       long_name="Fraction of implicit vertical velocity"                                         unit="#"   grid_ref="grid_T_3D" /> 
    8589 
    8690        <!-- next variables available with key_diahth --> 
     
    351355        <field id="uoce"         long_name="ocean current along i-axis"                             standard_name="sea_water_x_velocity"        unit="m/s"        grid_ref="grid_U_3D" /> 
    352356        <field id="uoce_e3u"     long_name="ocean current along i-axis  (thickness weighted)"                                                   unit="m/s"        grid_ref="grid_U_3D"  > uoce * e3u </field> 
     357        <field id="uoce2_e3u"    long_name="ocean current along i-axis squared (thickness weighted)"                                            unit="m3/s2"      grid_ref="grid_U_3D"  > uoce * uoce * e3u </field> 
    353358        <field id="ssu"          long_name="ocean surface current along i-axis"                                                                 unit="m/s"                             /> 
    354359        <field id="sbu"          long_name="ocean bottom current along i-axis"                                                                  unit="m/s"                             /> 
     
    405410        <field id="voce"         long_name="ocean current along j-axis"                             standard_name="sea_water_y_velocity"        unit="m/s"        grid_ref="grid_V_3D" /> 
    406411        <field id="voce_e3v"     long_name="ocean current along j-axis  (thickness weighted)"                                                   unit="m/s"        grid_ref="grid_V_3D"  > voce * e3v </field> 
     412        <field id="voce2_e3v"    long_name="ocean current along j-axis squared (thickness weighted)"                                            unit="m3/s2"      grid_ref="grid_V_3D"  > voce * voce * e3v </field> 
    407413        <field id="ssv"          long_name="ocean surface current along j-axis"                                                                 unit="m/s"                             /> 
    408414        <field id="sbv"          long_name="ocean bottom current along j-axis"                                                                  unit="m/s"                             /> 
     
    975981 
    976982    <field_group id="25h_grid_W" grid_ref="grid_W_3D" operation="instant"> 
    977       <field id="vovecrtz25h"         name="k current 25h mean"                 unit="m/s"      /> 
     983      <field id="vomecrtz25h"         name="k current 25h mean"                 unit="m/s"      /> 
    978984      <field id="avt25h"              name="vertical diffusivity25h mean"       unit="m2/s" /> 
    979985      <field id="avm25h"              name="vertical viscosity 25h mean"        unit="m2/s" /> 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/NST/agrif_oce_interp.F90

    r11053 r11099  
    115115            END DO 
    116116            DO jj = 1, jpj 
    117                uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu_a(ibdy1:ibdy2,jj) 
     117               uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    118118            END DO 
    119119         ENDIF 
     
    135135         END DO 
    136136         DO jj=1,jpj 
    137             zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 
     137            zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    138138         END DO 
    139139             
     
    154154            END DO 
    155155            DO jj = 1, jpj 
    156                zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj) 
     156               zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv(ibdy1:ibdy2,jj,Krhs_a) 
    157157            END DO 
    158158            DO jk = 1, jpkm1 
     
    186186            END DO 
    187187            DO jj = 1, jpj 
    188                uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu_a(ibdy1:ibdy2,jj) 
     188               uu_b(ibdy1:ibdy2,jj,Krhs_a) = uu_b(ibdy1:ibdy2,jj,Krhs_a) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    189189            END DO 
    190190         ENDIF 
     
    206206         END DO 
    207207         DO jj=1,jpj 
    208             zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 
     208            zub(ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu(ibdy1:ibdy2,jj,Krhs_a) 
    209209         END DO 
    210210             
     
    227227            END DO 
    228228            DO jj = 1, jpj 
    229                zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj) 
     229               zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv(ibdy1:ibdy2,jj,Krhs_a) 
    230230            END DO 
    231231            DO jk = 1, jpkm1 
     
    259259            END DO 
    260260            DO ji=1,jpi 
    261                vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv_a(ji,jbdy1:jbdy2) 
     261               vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    262262            END DO 
    263263         ENDIF 
     
    279279         END DO 
    280280         DO ji = 1, jpi 
    281             zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 
     281            zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    282282         END DO 
    283283 
     
    298298            END DO 
    299299            DO ji = 1, jpi 
    300                zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2) 
     300               zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu(ji,jbdy1:jbdy2,Krhs_a) 
    301301            END DO 
    302302                
     
    331331            END DO 
    332332            DO ji=1,jpi 
    333                vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv_a(ji,jbdy1:jbdy2) 
     333               vv_b(ji,jbdy1:jbdy2,Krhs_a) = vv_b(ji,jbdy1:jbdy2,Krhs_a) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    334334            END DO 
    335335         ENDIF 
     
    351351         END DO 
    352352         DO ji = 1, jpi 
    353             zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 
     353            zvb(ji,jbdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv(ji,jbdy1:jbdy2,Krhs_a) 
    354354         END DO 
    355355 
     
    372372            END DO 
    373373            DO ji = 1, jpi 
    374                zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2) 
     374               zub(ji,jbdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu(ji,jbdy1:jbdy2,Krhs_a) 
    375375            END DO 
    376376                
     
    10601060      ! 
    10611061      IF( before ) THEN  
    1062          ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * uu_b(i1:i2,j1:j2,Kmm_a) 
     1062         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) 
    10631063      ELSE 
    10641064         western_side  = (nb == 1).AND.(ndir == 1) 
     
    11131113      !  
    11141114      IF( before ) THEN  
    1115          ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vv_b(i1:i2,j1:j2,Kmm_a) 
     1115         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) 
    11161116      ELSE 
    11171117         western_side  = (nb == 1).AND.(ndir == 1) 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/NST/agrif_oce_update.F90

    r11053 r11099  
    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      ! 
     
    632632         IF (western_side) THEN 
    633633            DO jj=j1,j2 
    634                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) 
     634               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) 
    635635               uu_b(i1-1,jj,Kmm_a) = uu_b(i1-1,jj,Kmm_a) + zcor 
    636636               DO jk=1,jpkm1 
     
    642642         IF (eastern_side) THEN 
    643643            DO jj=j1,j2 
    644                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) 
     644               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) 
    645645               uu_b(i2+1,jj,Kmm_a) = uu_b(i2+1,jj,Kmm_a) + zcor 
    646646               DO jk=1,jpkm1 
     
    822822         IF (southern_side) THEN 
    823823            DO ji=i1,i2 
    824                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) 
     824               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) 
    825825               vv_b(ji,j1-1,Kmm_a) = vv_b(ji,j1-1,Kmm_a) + zcor 
    826826               DO jk=1,jpkm1 
     
    832832         IF (northern_side) THEN 
    833833            DO ji=i1,i2 
    834                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) 
     834               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) 
    835835               vv_b(ji,j2+1,Kmm_a) = vv_b(ji,j2+1,Kmm_a) + zcor 
    836836               DO jk=1,jpkm1 
     
    862862         DO jj=j1,j2 
    863863            DO ji=i1,i2 
    864                tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu_n(ji,jj) * e2u(ji,jj) 
     864               tabres(ji,jj) = zrhoy * uu_b(ji,jj,Kmm_a) * hu(ji,jj,Kmm_a) * e2u(ji,jj) 
    865865            END DO 
    866866         END DO 
     
    876876               END DO 
    877877               ! 
    878                zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu_n(ji,jj) 
     878               zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a) 
    879879               DO jk=1,jpkm1               
    880880                  uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk)            
     
    884884               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 
    885885                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    886                      zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu_a(ji,jj)) * r1_hu_b(ji,jj) 
     886                     zcorr = (tabres(ji,jj) - uu_b(ji,jj,Kmm_a) * hu(ji,jj,Krhs_a)) * r1_hu(ji,jj,Kbb_a) 
    887887                     uu_b(ji,jj,Kbb_a) = uu_b(ji,jj,Kbb_a) + atfp * zcorr * umask(ji,jj,1) 
    888888                  END IF 
    889889               ENDIF     
    890                uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu_n(ji,jj) * umask(ji,jj,1) 
     890               uu_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hu(ji,jj,Kmm_a) * umask(ji,jj,1) 
    891891               !        
    892892               ! Correct "before" velocities to hold correct bt component: 
     
    896896               END DO 
    897897               ! 
    898                zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu_b(ji,jj) 
     898               zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a) 
    899899               DO jk=1,jpkm1               
    900900                  uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk)            
     
    928928         DO jj=j1,j2 
    929929            DO ji=i1,i2 
    930                tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv_n(ji,jj) * e1v(ji,jj)  
     930               tabres(ji,jj) = zrhox * vv_b(ji,jj,Kmm_a) * hv(ji,jj,Kmm_a) * e1v(ji,jj)  
    931931            END DO 
    932932         END DO 
     
    942942               END DO 
    943943               ! 
    944                zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv_n(ji,jj) 
     944               zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a) 
    945945               DO jk=1,jpkm1               
    946946                  vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk)            
     
    950950               IF ( .NOT.ln_dynspg_ts .OR. (ln_dynspg_ts.AND.(.NOT.ln_bt_fw)) ) THEN 
    951951                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    952                      zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv_a(ji,jj)) * r1_hv_b(ji,jj) 
     952                     zcorr = (tabres(ji,jj) - vv_b(ji,jj,Kmm_a) * hv(ji,jj,Krhs_a)) * r1_hv(ji,jj,Kbb_a) 
    953953                     vv_b(ji,jj,Kbb_a) = vv_b(ji,jj,Kbb_a) + atfp * zcorr * vmask(ji,jj,1) 
    954954                  END IF 
    955955               ENDIF               
    956                vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv_n(ji,jj) * vmask(ji,jj,1) 
     956               vv_b(ji,jj,Kmm_a) = tabres(ji,jj) * r1_hv(ji,jj,Kmm_a) * vmask(ji,jj,1) 
    957957               !        
    958958               ! Correct "before" velocities to hold correct bt component: 
     
    962962               END DO 
    963963               ! 
    964                zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv_b(ji,jj) 
     964               zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a) 
    965965               DO jk=1,jpkm1               
    966966                  vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk)            
     
    13731373         ! 
    13741374         ! Update total depth: 
    1375          ht_n(i1:i2,j1:j2) = 0._wp 
     1375         ht(i1:i2,j1:j2) = 0._wp 
    13761376         DO jk = 1, jpkm1 
    1377             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) 
     1377            ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + e3t(i1:i2,j1:j2,jk,Kmm_a) * tmask(i1:i2,j1:j2,jk) 
    13781378         END DO 
    13791379         ! 
     
    13821382         gdept(i1:i2,j1:j2,1,Kmm_a) = 0.5_wp * e3w(i1:i2,j1:j2,1,Kmm_a) 
    13831383         gdepw(i1:i2,j1:j2,1,Kmm_a) = 0.0_wp 
    1384          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 
     1384         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 
    13851385         ! 
    13861386         DO jk = 2, jpk 
     
    13931393               gdept(ji,jj,jk,Kmm_a) =      zcoef  * ( gdepw(ji,jj,jk  ,Kmm_a) + 0.5 * e3w(ji,jj,jk,Kmm_a))  & 
    13941394                   &               + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm_a) +       e3w(ji,jj,jk,Kmm_a))  
    1395                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 
     1395               gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm_a) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh 
    13961396               END DO 
    13971397            END DO 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/OCE/ASM/asminc.F90

    r10969 r11099  
    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_rewrite_time_filterswap/src/OCE/DOM/domain.F90

    r11053 r11099  
    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 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/OCE/TRA/traatf.F90

    r11057 r11099  
    146146      ENDIF 
    147147 
    148       IF( neuler == 0 .AND. kt == nit000 ) THEN       ! Euler time-stepping at first time-step (only swap) 
    149          DO jn = 1, jpts 
    150             DO jk = 1, jpkm1 
    151                pts(:,:,jk,jn,Kmm) = pts(:,:,jk,jn,Kaa)     
    152             END DO 
    153          END DO 
     148      IF( neuler == 0 .AND. kt == nit000 ) THEN       ! Euler time-stepping  
     149         ! 
    154150         IF (l_trdtra .AND. .NOT. ln_linssh ) THEN   ! Zero Asselin filter contribution must be explicitly written out since for vvl 
    155151            !                                        ! Asselin filter is output by tra_atf_vvl that is not called on this time step 
     
    175171         zfact = 1._wp / r2dt              
    176172         DO jk = 1, jpkm1 
    177             ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kbb) - ztrdt(:,:,jk) ) * zfact 
    178             ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kbb) - ztrds(:,:,jk) ) * zfact 
     173            ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kmm) - ztrdt(:,:,jk) ) * zfact 
     174            ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kmm) - ztrds(:,:,jk) ) * zfact 
    179175         END DO 
    180176         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_atf, ztrdt ) 
     
    192188 
    193189 
    194    SUBROUTINE tra_atf_fix( kt, kit000, Kbb, Kmm, Kaa, cdtype, pt, kjpt ) 
     190   SUBROUTINE tra_atf_fix( kt, Kbb, Kmm, Kaa, kit000, cdtype, pt, kjpt ) 
    195191      !!---------------------------------------------------------------------- 
    196192      !!                   ***  ROUTINE tra_atf_fix  *** 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/OCE/TRA/trasbc.F90

    r10985 r11099  
    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_rewrite_time_filterswap/src/TOP/PISCES/SED/trcdmp_sed.F90

    r10975 r11099  
    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_rewrite_time_filterswap/src/TOP/TRP/trcatf.F90

    r11057 r11099  
    3030   USE trd_oce 
    3131   USE trdtra 
    32    USE tranxt 
     32   USE traatf 
    3333   USE bdy_oce   , ONLY: ln_bdy 
    3434   USE trcbdy          ! BDY open boundaries 
     
    143143      ENDIF 
    144144      !                                ! Leap-Frog + Asselin filter time stepping 
    145       IF( .NOT.( (neuler == 0 .AND. kt == nittrc000) .OR. ln_top_euler ) ) THEN    ! Only time filter if not an Euler timestep 
     145      IF( (neuler == 0 .AND. kt == nittrc000) .OR. ln_top_euler ) THEN    ! Euler time-stepping  
     146         ! 
     147         IF (l_trdtrc .AND. .NOT. ln_linssh ) THEN   ! Zero Asselin filter contribution must be explicitly written out since for vvl 
     148            !                                        ! Asselin filter is output by tra_nxt_vvl that is not called on this time step 
     149            ztrdt(:,:,:,:) = 0._wp             
     150            DO jn = 1, jptra 
     151               CALL trd_tra( kt, Kmm, Kaa, 'TRC', jn, jptra_atf, ztrdt(:,:,:,jn) ) 
     152            ENDDO 
     153         END IF 
     154         ! 
     155      ELSE      
    146156         IF( .NOT. l_offline ) THEN ! Leap-Frog + Asselin filter time stepping 
    147             IF( ln_linssh ) THEN   ;   CALL tra_nxt_fix( kt, Kbb, Kmm, Kaa, nittrc000,         'TRC', ptr, jptra )                     !     linear ssh 
    148             ELSE                   ;   CALL tra_nxt_vvl( kt, Kbb, Kmm, Kaa, nittrc000, rdttrc, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh 
     157            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix( kt, Kbb, Kmm, Kaa, nittrc000,         'TRC', ptr, jptra )                     !     linear ssh 
     158            ELSE                   ;   CALL tra_atf_vvl( kt, Kbb, Kmm, Kaa, nittrc000, rdttrc, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh 
    149159            ENDIF 
    150160         ELSE 
    151                                        CALL trc_atf_off( kt, Kbb, Kmm, Kaa )       ! offline  
     161                                       CALL trc_atf_off( kt, Kbb, Kmm, Kaa, ptr )       ! offline  
    152162         ENDIF 
    153163         ! 
     
    260270   !!   Default option                                         Empty module 
    261271   !!---------------------------------------------------------------------- 
     272   USE par_oce 
     273   USE par_trc 
    262274CONTAINS 
    263275   SUBROUTINE trc_atf( kt, Kbb, Kmm, Kaa, ptr )   
    264       INTEGER                                   , INTENT(in) :: kt 
    265       INTEGER                                   , INTENT(in) :: Kbb, Kmm, Kaa ! time level indices 
    266       REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr            ! passive tracers 
     276      INTEGER                                   , INTENT(in)    :: kt 
     277      INTEGER,                                    INTENT(in   ) :: Kbb, Kmm, Kaa ! time level indices 
     278      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr           ! passive tracers and RHS of tracer equation 
    267279      WRITE(*,*) 'trc_atf: You should not have seen this print! error?', kt 
    268280   END SUBROUTINE trc_atf 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/TOP/TRP/trcsbc.F90

    r10985 r11099  
    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_rewrite_time_filterswap/tests/VORTEX/MY_SRC/domvvl.F90

    r10572 r11099  
    106106      !!              - Regrid: e3(u/v)_n 
    107107      !!                        e3(u/v)_b        
    108       !!                        e3w_n            
     108      !!                        e3w(:,:,:,Kmm)            
    109109      !!                        e3(u/v)w_b       
    110110      !!                        e3(u/v)w_n       
    111       !!                        gdept_n, gdepw_n and gde3w_n 
     111      !!                        gdept(:,:,:,Kmm), gdepw(:,:,:,Kmm) and gde3w 
    112112      !!              - h(t/u/v)_0 
    113113      !!              - frq_rst_e3t and frq_rst_hdv 
     
    131131      !                    ! Read or initialize e3t_(b/n), tilde_e3t_(b/n) and hdiv_lf 
    132132      CALL dom_vvl_rst( nit000, 'READ' ) 
    133       e3t_a(:,:,jpk) = e3t_0(:,:,jpk)  ! last level always inside the sea floor set one for all 
     133      e3t(:,:,jpk,Kaa) = e3t_0(:,:,jpk)  ! last level always inside the sea floor set one for all 
    134134      ! 
    135135      !                    !== Set of all other vertical scale factors  ==!  (now and before) 
    136136      !                                ! Horizontal interpolation of e3t 
    137       CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:), 'U' )    ! from T to U 
    138       CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:), 'U' ) 
    139       CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:), 'V' )    ! from T to V  
    140       CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:), 'V' ) 
    141       CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:), 'F' )    ! from U to F 
     137      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' )    ! from T to U 
     138      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
     139      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' )    ! from T to V  
     140      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
     141      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' )    ! from U to F 
    142142      !                                ! Vertical interpolation of e3t,u,v  
    143       CALL dom_vvl_interpol( e3t_n(:,:,:), e3w_n (:,:,:), 'W'  )  ! from T to W 
    144       CALL dom_vvl_interpol( e3t_b(:,:,:), e3w_b (:,:,:), 'W'  ) 
    145       CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' )  ! from U to UW 
    146       CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' ) 
    147       CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' )  ! from V to UW 
    148       CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' ) 
     143      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  )  ! from T to W 
     144      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3w (:,:,:,Kbb), 'W'  ) 
     145      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' )  ! from U to UW 
     146      CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     147      CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' )  ! from V to UW 
     148      CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    149149 
    150150      ! We need to define e3[tuv]_a for AGRIF initialisation (should not be a problem for the restartability...) 
    151       e3t_a(:,:,:) = e3t_n(:,:,:) 
    152       e3u_a(:,:,:) = e3u_n(:,:,:) 
    153       e3v_a(:,:,:) = e3v_n(:,:,:) 
     151      e3t(:,:,:,Kaa) = e3t(:,:,:,Kmm) 
     152      e3u(:,:,:,Kaa) = e3u(:,:,:,Kmm) 
     153      e3v(:,:,:,Kaa) = e3v(:,:,:,Kmm) 
    154154      ! 
    155155      !                    !==  depth of t and w-point  ==!   (set the isf depth as it is in the initial timestep) 
    156       gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1)       ! reference to the ocean surface (used for MLD and light penetration) 
    157       gdepw_n(:,:,1) = 0.0_wp 
    158       gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:)  ! reference to a common level z=0 for hpg 
    159       gdept_b(:,:,1) = 0.5_wp * e3w_b(:,:,1) 
    160       gdepw_b(:,:,1) = 0.0_wp 
     156      gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm)       ! reference to the ocean surface (used for MLD and light penetration) 
     157      gdepw(:,:,1,Kmm) = 0.0_wp 
     158      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm)  ! reference to a common level z=0 for hpg 
     159      gdept(:,:,1,Kbb) = 0.5_wp * e3w(:,:,1,Kbb) 
     160      gdepw(:,:,1,Kbb) = 0.0_wp 
    161161      DO jk = 2, jpk                               ! vertical sum 
    162162         DO jj = 1,jpj 
     
    165165               !                             ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 
    166166               !                             ! 0.5 where jk = mikt      
    167 !!gm ???????   BUG ?  gdept_n as well as gde3w_n  does not include the thickness of ISF ?? 
     167!!gm ???????   BUG ?  gdept(:,:,:,Kmm) as well as gde3w  does not include the thickness of ISF ?? 
    168168               zcoef = ( tmask(ji,jj,jk) - wmask(ji,jj,jk) ) 
    169                gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
    170                gdept_n(ji,jj,jk) =      zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk))  & 
    171                   &                + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk))  
    172                gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - sshn(ji,jj) 
    173                gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1) 
    174                gdept_b(ji,jj,jk) =      zcoef  * ( gdepw_b(ji,jj,jk  ) + 0.5 * e3w_b(ji,jj,jk))  & 
    175                   &                + (1-zcoef) * ( gdept_b(ji,jj,jk-1) +       e3w_b(ji,jj,jk))  
     169               gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     170               gdept(ji,jj,jk,Kmm) =      zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm))  & 
     171                  &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm))  
     172               gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
     173               gdepw(ji,jj,jk,Kbb) = gdepw(ji,jj,jk-1,Kbb) + e3t(ji,jj,jk-1,Kbb) 
     174               gdept(ji,jj,jk,Kbb) =      zcoef  * ( gdepw(ji,jj,jk  ,Kbb) + 0.5 * e3w(ji,jj,jk,Kbb))  & 
     175                  &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb) +       e3w(ji,jj,jk,Kbb))  
    176176            END DO 
    177177         END DO 
     
    179179      ! 
    180180      !                    !==  thickness of the water column  !!   (ocean portion only) 
    181       ht_n(:,:) = e3t_n(:,:,1) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
    182       hu_b(:,:) = e3u_b(:,:,1) * umask(:,:,1) 
    183       hu_n(:,:) = e3u_n(:,:,1) * umask(:,:,1) 
    184       hv_b(:,:) = e3v_b(:,:,1) * vmask(:,:,1) 
    185       hv_n(:,:) = e3v_n(:,:,1) * vmask(:,:,1) 
     181      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
     182      hu(:,:,Kbb) = e3u(:,:,1,Kbb) * umask(:,:,1) 
     183      hu(:,:,Kmm) = e3u(:,:,1,Kmm) * umask(:,:,1) 
     184      hv(:,:,Kbb) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
     185      hv(:,:,Kmm) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
    186186      DO jk = 2, jpkm1 
    187          ht_n(:,:) = ht_n(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
    188          hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk) 
    189          hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk) 
    190          hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk) 
    191          hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk) 
     187         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     188         hu(:,:,Kbb) = hu(:,:,Kbb) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
     189         hu(:,:,Kmm) = hu(:,:,Kmm) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
     190         hv(:,:,Kbb) = hv(:,:,Kbb) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
     191         hv(:,:,Kmm) = hv(:,:,Kmm) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    192192      END DO 
    193193      ! 
    194194      !                    !==  inverse of water column thickness   ==!   (u- and v- points) 
    195       r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
    196       r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
    197       r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
    198       r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     195      r1_hu(:,:,Kbb) = ssumask(:,:) / ( hu(:,:,Kbb) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
     196      r1_hu(:,:,Kmm) = ssumask(:,:) / ( hu(:,:,Kmm) + 1._wp - ssumask(:,:) ) 
     197      r1_hv(:,:,Kbb) = ssvmask(:,:) / ( hv(:,:,Kbb) + 1._wp - ssvmask(:,:) ) 
     198      r1_hv(:,:,Kmm) = ssvmask(:,:) / ( hv(:,:,Kmm) + 1._wp - ssvmask(:,:) ) 
    199199 
    200200      !                    !==   z_tilde coordinate case  ==!   (Restoring frequencies) 
     
    321321      !                                                ! --------------------------------------------- ! 
    322322      ! 
    323       z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * ssmask(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - ssmask(:,:) ) 
     323      z_scale(:,:) = ( ssh(:,:,Kaa) - ssh(:,:,Kbb) ) * ssmask(:,:) / ( ht_0(:,:) + ssh(:,:,Kmm) + 1. - ssmask(:,:) ) 
    324324      DO jk = 1, jpkm1 
    325          ! formally this is the same as e3t_a = e3t_0*(1+ssha/ht_0) 
    326          e3t_a(:,:,jk) = e3t_b(:,:,jk) + e3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
     325         ! formally this is the same as e3t(:,:,:,Kaa) = e3t_0*(1+ssha/ht_0) 
     326         e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kbb) + e3t(:,:,jk,Kmm) * z_scale(:,:) * tmask(:,:,jk) 
    327327      END DO 
    328328      ! 
     
    337337         zht(:,:)   = 0._wp 
    338338         DO jk = 1, jpkm1 
    339             zhdiv(:,:) = zhdiv(:,:) + e3t_n(:,:,jk) * hdivn(:,:,jk) 
    340             zht  (:,:) = zht  (:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
     339            zhdiv(:,:) = zhdiv(:,:) + e3t(:,:,jk,Kmm) * hdiv(:,:,jk) 
     340            zht  (:,:) = zht  (:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    341341         END DO 
    342342         zhdiv(:,:) = zhdiv(:,:) / ( zht(:,:) + 1. - tmask_i(:,:) ) 
     
    348348               DO jk = 1, jpkm1 
    349349                  hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - rdt * frq_rst_hdv(:,:)   & 
    350                      &          * ( hdiv_lf(:,:,jk) - e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) ) 
     350                     &          * ( hdiv_lf(:,:,jk) - e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) ) 
    351351               END DO 
    352352            ENDIF 
     
    361361         IF( ln_vvl_ztilde ) THEN     ! z_tilde case 
    362362            DO jk = 1, jpkm1 
    363                tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - ( e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 
     363               tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - ( e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 
    364364            END DO 
    365365         ELSE                         ! layer case 
    366366            DO jk = 1, jpkm1 
    367                tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) -   e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) * tmask(:,:,jk) 
     367               tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) -   e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) * tmask(:,:,jk) 
    368368            END DO 
    369369         ENDIF 
     
    476476            zht(:,:)  = zht(:,:) + tilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
    477477         END DO 
    478          z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - ssmask(:,:) ) 
     478         z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + ssh(:,:,Kmm) + 1. - ssmask(:,:) ) 
    479479         DO jk = 1, jpkm1 
    480             dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + e3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
     480            dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + e3t(:,:,jk,Kmm) * z_scale(:,:) * tmask(:,:,jk) 
    481481         END DO 
    482482 
     
    486486      !                                           ! ---baroclinic part--------- ! 
    487487         DO jk = 1, jpkm1 
    488             e3t_a(:,:,jk) = e3t_a(:,:,jk) + dtilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
     488            e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kaa) + dtilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
    489489         END DO 
    490490      ENDIF 
     
    501501         zht(:,:) = 0.0_wp 
    502502         DO jk = 1, jpkm1 
    503             zht(:,:) = zht(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
    504          END DO 
    505          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + sshn(:,:) - zht(:,:) ) ) 
     503            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     504         END DO 
     505         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kmm) - zht(:,:) ) ) 
    506506         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    507          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshn-SUM(e3t_n))) =', z_tmax 
     507         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshn-SUM(e3t(:,:,:,Kmm)))) =', z_tmax 
    508508         ! 
    509509         zht(:,:) = 0.0_wp 
    510510         DO jk = 1, jpkm1 
    511             zht(:,:) = zht(:,:) + e3t_a(:,:,jk) * tmask(:,:,jk) 
    512          END DO 
    513          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssha(:,:) - zht(:,:) ) ) 
     511            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kaa) * tmask(:,:,jk) 
     512         END DO 
     513         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kaa) - zht(:,:) ) ) 
    514514         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    515          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+ssha-SUM(e3t_a))) =', z_tmax 
     515         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+ssha-SUM(e3t(:,:,:,Kaa)))) =', z_tmax 
    516516         ! 
    517517         zht(:,:) = 0.0_wp 
    518518         DO jk = 1, jpkm1 
    519             zht(:,:) = zht(:,:) + e3t_b(:,:,jk) * tmask(:,:,jk) 
    520          END DO 
    521          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + sshb(:,:) - zht(:,:) ) ) 
     519            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kbb) * tmask(:,:,jk) 
     520         END DO 
     521         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kbb) - zht(:,:) ) ) 
    522522         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    523          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshb-SUM(e3t_b))) =', z_tmax 
    524          ! 
    525          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( sshb(:,:) ) ) 
     523         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshb-SUM(e3t(:,:,:,Kbb)))) =', z_tmax 
     524         ! 
     525         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kbb) ) ) 
    526526         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    527          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(sshb))) =', z_tmax 
    528          ! 
    529          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( sshn(:,:) ) ) 
     527         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kbb)))) =', z_tmax 
     528         ! 
     529         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kmm) ) ) 
    530530         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    531          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(sshn))) =', z_tmax 
    532          ! 
    533          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssha(:,:) ) ) 
     531         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kmm)))) =', z_tmax 
     532         ! 
     533         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kaa) ) ) 
    534534         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    535          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssha))) =', z_tmax 
     535         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kaa)))) =', z_tmax 
    536536      END IF 
    537537 
     
    540540      ! *********************************** ! 
    541541 
    542       CALL dom_vvl_interpol( e3t_a(:,:,:), e3u_a(:,:,:), 'U' ) 
    543       CALL dom_vvl_interpol( e3t_a(:,:,:), e3v_a(:,:,:), 'V' ) 
     542      CALL dom_vvl_interpol( e3t(:,:,:,Kaa), e3u(:,:,:,Kaa), 'U' ) 
     543      CALL dom_vvl_interpol( e3t(:,:,:,Kaa), e3v(:,:,:,Kaa), 'V' ) 
    544544 
    545545      ! *********************************** ! 
     
    547547      ! *********************************** ! 
    548548 
    549       hu_a(:,:) = e3u_a(:,:,1) * umask(:,:,1) 
    550       hv_a(:,:) = e3v_a(:,:,1) * vmask(:,:,1) 
     549      hu(:,:,Kaa) = e3u(:,:,1,Kaa) * umask(:,:,1) 
     550      hv(:,:,Kaa) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
    551551      DO jk = 2, jpkm1 
    552          hu_a(:,:) = hu_a(:,:) + e3u_a(:,:,jk) * umask(:,:,jk) 
    553          hv_a(:,:) = hv_a(:,:) + e3v_a(:,:,jk) * vmask(:,:,jk) 
     552         hu(:,:,Kaa) = hu(:,:,Kaa) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     553         hv(:,:,Kaa) = hv(:,:,Kaa) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    554554      END DO 
    555555      !                                        ! Inverse of the local depth 
    556556!!gm BUG ?  don't understand the use of umask_i here ..... 
    557       r1_hu_a(:,:) = ssumask(:,:) / ( hu_a(:,:) + 1._wp - ssumask(:,:) ) 
    558       r1_hv_a(:,:) = ssvmask(:,:) / ( hv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     557      r1_hu(:,:,Kaa) = ssumask(:,:) / ( hu(:,:,Kaa) + 1._wp - ssumask(:,:) ) 
     558      r1_hv(:,:,Kaa) = ssvmask(:,:) / ( hv(:,:,Kaa) + 1._wp - ssvmask(:,:) ) 
    559559      ! 
    560560      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_nxt') 
     
    578578      !!               - Recompute: 
    579579      !!                    e3(u/v)_b        
    580       !!                    e3w_n            
     580      !!                    e3w(:,:,:,Kmm)            
    581581      !!                    e3(u/v)w_b       
    582582      !!                    e3(u/v)w_n       
    583       !!                    gdept_n, gdepw_n  and gde3w_n 
     583      !!                    gdept(:,:,:,Kmm), gdepw(:,:,:,Kmm)  and gde3w 
    584584      !!                    h(u/v) and h(u/v)r 
    585585      !! 
     
    615615         tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) 
    616616      ENDIF 
    617       gdept_b(:,:,:) = gdept_n(:,:,:) 
    618       gdepw_b(:,:,:) = gdepw_n(:,:,:) 
    619  
    620       e3t_n(:,:,:) = e3t_a(:,:,:) 
    621       e3u_n(:,:,:) = e3u_a(:,:,:) 
    622       e3v_n(:,:,:) = e3v_a(:,:,:) 
     617      gdept(:,:,:,Kbb) = gdept(:,:,:,Kmm) 
     618      gdepw(:,:,:,Kbb) = gdepw(:,:,:,Kmm) 
     619 
     620      e3t(:,:,:,Kmm) = e3t(:,:,:,Kaa) 
     621      e3u(:,:,:,Kmm) = e3u(:,:,:,Kaa) 
     622      e3v(:,:,:,Kmm) = e3v(:,:,:,Kaa) 
    623623 
    624624      ! Compute all missing vertical scale factor and depths 
     
    626626      ! Horizontal scale factor interpolations 
    627627      ! -------------------------------------- 
    628       ! - ML - e3u_b and e3v_b are allready computed in dynnxt 
    629       ! - JC - hu_b, hv_b, hur_b, hvr_b also 
     628      ! - ML - e3u(:,:,:,Kbb) and e3v(:,:,:,Kbb) are allready computed in dynnxt 
     629      ! - JC - hu(:,:,Kbb), hv(:,:,Kbb), hur_b, hvr_b also 
    630630       
    631       CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:), 'F'  ) 
     631      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F'  ) 
    632632       
    633633      ! Vertical scale factor interpolations 
    634       CALL dom_vvl_interpol( e3t_n(:,:,:),  e3w_n(:,:,:), 'W'  ) 
    635       CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' ) 
    636       CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' ) 
    637       CALL dom_vvl_interpol( e3t_b(:,:,:),  e3w_b(:,:,:), 'W'  ) 
    638       CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' ) 
    639       CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' ) 
     634      CALL dom_vvl_interpol( e3t(:,:,:,Kmm),  e3w(:,:,:,Kmm), 'W'  ) 
     635      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
     636      CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
     637      CALL dom_vvl_interpol( e3t(:,:,:,Kbb),  e3w(:,:,:,Kbb), 'W'  ) 
     638      CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     639      CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    640640 
    641641      ! t- and w- points depth (set the isf depth as it is in the initial step) 
    642       gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1) 
    643       gdepw_n(:,:,1) = 0.0_wp 
    644       gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:) 
     642      gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
     643      gdepw(:,:,1,Kmm) = 0.0_wp 
     644      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    645645      DO jk = 2, jpk 
    646646         DO jj = 1,jpj 
     
    649649                                                                 ! 1 for jk = mikt 
    650650               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
    651                gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
    652                gdept_n(ji,jj,jk) =    zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk) )  & 
    653                    &             + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk) )  
    654                gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - sshn(ji,jj) 
     651               gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     652               gdept(ji,jj,jk,Kmm) =    zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm) )  & 
     653                   &             + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm) )  
     654               gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
    655655            END DO 
    656656         END DO 
     
    659659      ! Local depth and Inverse of the local depth of the water 
    660660      ! ------------------------------------------------------- 
    661       hu_n(:,:) = hu_a(:,:)   ;   r1_hu_n(:,:) = r1_hu_a(:,:) 
    662       hv_n(:,:) = hv_a(:,:)   ;   r1_hv_n(:,:) = r1_hv_a(:,:) 
    663       ! 
    664       ht_n(:,:) = e3t_n(:,:,1) * tmask(:,:,1) 
     661      hu(:,:,Kmm) = hu(:,:,Kaa)   ;   r1_hu(:,:,Kmm) = r1_hu(:,:,Kaa) 
     662      hv(:,:,Kmm) = hv(:,:,Kaa)   ;   r1_hv(:,:,Kmm) = r1_hv(:,:,Kaa) 
     663      ! 
     664      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
    665665      DO jk = 2, jpkm1 
    666          ht_n(:,:) = ht_n(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
     666         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    667667      END DO 
    668668 
     
    806806         IF( ln_rstart ) THEN                   !* Read the restart file 
    807807            CALL rst_read_open                  !  open the restart file if necessary 
    808             CALL iom_get( numror, jpdom_autoglo, 'sshn'   , sshn, ldxios = lrxios    ) 
     808            CALL iom_get( numror, jpdom_autoglo, 'sshn'   , ssh(:,:,Kmm), ldxios = lrxios    ) 
    809809            ! 
    810810            id1 = iom_varid( numror, 'e3t_b', ldstop = .FALSE. ) 
     
    817817            !                             ! --------- ! 
    818818            IF( MIN( id1, id2 ) > 0 ) THEN       ! all required arrays exist 
    819                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t_b(:,:,:), ldxios = lrxios ) 
    820                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t_n(:,:,:), ldxios = lrxios ) 
     819               CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     820               CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
    821821               ! needed to restart if land processor not computed  
    822                IF(lwp) write(numout,*) 'dom_vvl_rst : e3t_b and e3t_n found in restart files' 
     822               IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files' 
    823823               WHERE ( tmask(:,:,:) == 0.0_wp )  
    824                   e3t_n(:,:,:) = e3t_0(:,:,:) 
    825                   e3t_b(:,:,:) = e3t_0(:,:,:) 
     824                  e3t(:,:,:,Kmm) = e3t_0(:,:,:) 
     825                  e3t(:,:,:,Kbb) = e3t_0(:,:,:) 
    826826               END WHERE 
    827827               IF( neuler == 0 ) THEN 
    828                   e3t_b(:,:,:) = e3t_n(:,:,:) 
     828                  e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    829829               ENDIF 
    830830            ELSE IF( id1 > 0 ) THEN 
    831                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_n not found in restart files' 
     831               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart files' 
    832832               IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.' 
    833833               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    834                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t_b(:,:,:), ldxios = lrxios ) 
    835                e3t_n(:,:,:) = e3t_b(:,:,:) 
     834               CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     835               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
    836836               neuler = 0 
    837837            ELSE IF( id2 > 0 ) THEN 
    838                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_b not found in restart files' 
     838               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kbb) not found in restart files' 
    839839               IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.' 
    840840               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    841                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t_n(:,:,:), ldxios = lrxios ) 
    842                e3t_b(:,:,:) = e3t_n(:,:,:) 
     841               CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
     842               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    843843               neuler = 0 
    844844            ELSE 
    845                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_n not found in restart file' 
     845               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart file' 
    846846               IF(lwp) write(numout,*) 'Compute scale factor from sshn' 
    847847               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    848848               DO jk = 1, jpk 
    849                   e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:) ) & 
     849                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & 
    850850                      &                          / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    851851                      &          + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk)) 
    852852               END DO 
    853                e3t_b(:,:,:) = e3t_n(:,:,:) 
     853               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    854854               neuler = 0 
    855855            ENDIF 
     
    888888               IF( cn_cfg == 'wad' ) THEN 
    889889                  ! Wetting and drying test case 
    890                   CALL usr_def_istate( gdept_b, tmask, tsb, ub, vb, sshb  ) 
    891                   tsn  (:,:,:,:) = tsb (:,:,:,:)       ! set now values from to before ones 
    892                   sshn (:,:)     = sshb(:,:) 
    893                   un   (:,:,:)   = ub  (:,:,:) 
    894                   vn   (:,:,:)   = vb  (:,:,:) 
     890                  CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  ) 
     891                  ts  (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb)       ! set now values from to before ones 
     892                  ssh (:,:,Kmm)     = ssh(:,:,Kbb) 
     893                  uu   (:,:,:,Kmm)   = uu  (:,:,:,Kbb) 
     894                  vv   (:,:,:,Kmm)   = vv  (:,:,:,Kbb) 
    895895               ELSE 
    896896                  ! if not test case 
    897                   sshn(:,:) = -ssh_ref 
    898                   sshb(:,:) = -ssh_ref 
     897                  ssh(:,:,Kmm) = -ssh_ref 
     898                  ssh(:,:,Kbb) = -ssh_ref 
    899899 
    900900                  DO jj = 1, jpj 
     
    902902                        IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth 
    903903 
    904                            sshb(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
    905                            sshn(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
    906                            ssha(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
     904                           ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 
     905                           ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 
     906                           ssh(ji,jj,Kaa) = rn_wdmin1 - (ht_0(ji,jj) ) 
    907907                        ENDIF 
    908908                     ENDDO 
     
    912912               ! Adjust vertical metrics for all wad 
    913913               DO jk = 1, jpk 
    914                   e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:)  ) & 
     914                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm)  ) & 
    915915                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    916916                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) ) 
    917917               END DO 
    918                e3t_b(:,:,:) = e3t_n(:,:,:) 
     918               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    919919 
    920920               DO ji = 1, jpi 
     
    928928            ELSE 
    929929               ! 
    930                ! usr_def_istate called here only to get sshb, that is needed to initialize e3t_b and e3t_n 
    931                CALL usr_def_istate( gdept_0, tmask, tsb, ub, vb, sshb  )   
     930               ! usr_def_istate called here only to get ssh(:,:,Kbb), that is needed to initialize e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) 
     931               CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )   
    932932               ! usr_def_istate will be called again in istate_init to initialize ts(bn), ssh(bn), u(bn) and v(bn) 
    933933               ! 
    934934               DO jk=1,jpk 
    935                   e3t_b(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshb(:,:) ) & 
     935                  e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) ) & 
    936936                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    937                     &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )   ! make sure e3t_b != 0 on land points 
     937                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )   ! make sure e3t(:,:,:,Kbb) != 0 on land points 
    938938               END DO 
    939                e3t_n(:,:,:) = e3t_b(:,:,:) 
    940                sshn(:,:) = sshb(:,:)   ! needed later for gde3w 
    941 !!$                e3t_n(:,:,:)=e3t_0(:,:,:) 
    942 !!$                e3t_b(:,:,:)=e3t_0(:,:,:) 
     939               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
     940               ssh(:,:,Kmm) = ssh(:,:,Kbb)   ! needed later for gde3w 
     941!!$                e3t(:,:,:,Kmm)=e3t_0(:,:,:) 
     942!!$                e3t(:,:,:,Kbb)=e3t_0(:,:,:) 
    943943               ! 
    944944            END IF           ! end of ll_wd edits 
     
    958958         !                                           ! all cases ! 
    959959         !                                           ! --------- ! 
    960          CALL iom_rstput( kt, nitrst, numrow, 'e3t_b', e3t_b(:,:,:), ldxios = lwxios ) 
    961          CALL iom_rstput( kt, nitrst, numrow, 'e3t_n', e3t_n(:,:,:), ldxios = lwxios ) 
     960         CALL iom_rstput( kt, nitrst, numrow, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lwxios ) 
     961         CALL iom_rstput( kt, nitrst, numrow, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lwxios ) 
    962962         !                                           ! ----------------------- ! 
    963963         IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN  ! z_tilde and layer cases ! 
Note: See TracChangeset for help on using the changeset viewer.