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 13237 for NEMO – NEMO

Changeset 13237 for NEMO


Ignore:
Timestamp:
2020-07-03T11:12:53+02:00 (4 years ago)
Author:
smasson
Message:

trunk: Mid-year merge, merge back KERNEL-06_techene_e3

Location:
NEMO/trunk/src
Files:
128 edited
6 copied

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/ICE/icedyn_rhg_evp.F90

    r13226 r13237  
    4949   !! * Substitutions 
    5050#  include "do_loop_substitute.h90" 
     51#  include "domzgr_substitute.h90" 
    5152   !!---------------------------------------------------------------------- 
    5253   !! NEMO/ICE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/ICE/iceistate.F90

    r13216 r13237  
    6666   INTEGER , PARAMETER ::   jp_hpd = 9           ! index of pnd depth        (m) 
    6767   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   si  ! structure of input fields (file informations, fields read) 
    68    !    
     68 
    6969   !! * Substitutions 
    7070#  include "do_loop_substitute.h90" 
     
    108108      REAL(wp), DIMENSION(jpi,jpj)     ::   zt_su_ini, zht_s_ini, zsm_i_ini, ztm_i_ini !data from namelist or nc file 
    109109      REAL(wp), DIMENSION(jpi,jpj)     ::   zapnd_ini, zhpnd_ini                       !data from namelist or nc file 
    110       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zti_3d , zts_3d                            !temporary arrays 
     110      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zti_3d , zts_3d                            !locak arrays 
    111111      !! 
    112112      REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   zhi_2d, zhs_2d, zai_2d, zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d 
     
    427427         ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rho0 
    428428         ! 
    429          IF( .NOT.ln_linssh ) THEN 
    430             ! 
    431             WHERE( ht_0(:,:) > 0 )   ;   z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 
    432             ELSEWHERE                ;   z2d(:,:) = 1._wp   ;   END WHERE 
    433             ! 
    434             DO jk = 1,jpkm1                     ! adjust initial vertical scale factors                 
    435                e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 
    436                e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 
    437                e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 
    438             END DO 
    439             ! 
    440             ! Reconstruction of all vertical scale factors at now and before time-steps 
    441             ! ========================================================================= 
    442             ! Horizontal scale factor interpolations 
    443             ! -------------------------------------- 
    444             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
    445             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
    446             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
    447             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
    448             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
    449             ! Vertical scale factor interpolations 
    450             ! ------------------------------------ 
    451             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
    452             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
    453             CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
    454             CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
    455             CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    456             ! t- and w- points depth 
    457             ! ---------------------- 
    458             !!gm not sure of that.... 
    459             gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
    460             gdepw(:,:,1,Kmm) = 0.0_wp 
    461             gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    462             DO jk = 2, jpk 
    463                gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk  ,Kmm) 
    464                gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
    465                gde3w(:,:,jk) = gdept(:,:,jk  ,Kmm) - ssh (:,:,Kmm) 
    466             END DO 
    467          ENDIF 
     429         IF( .NOT.ln_linssh )   CALL dom_vvl_zgr( Kbb, Kmm, Kaa )   ! interpolation scale factor, depth and water column 
     430! !!st 
     431!          IF( .NOT.ln_linssh ) THEN 
     432!             ! 
     433!             WHERE( ht_0(:,:) > 0 )   ;   z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 
     434!             ELSEWHERE                ;   z2d(:,:) = 1._wp   ;   END WHERE 
     435!             ! 
     436!             DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
     437!                e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 
     438!                e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 
     439!                e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 
     440!             END DO 
     441!             ! 
     442!             ! Reconstruction of all vertical scale factors at now and before time-steps 
     443!             ! ========================================================================= 
     444!             ! Horizontal scale factor interpolations 
     445!             ! -------------------------------------- 
     446!             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
     447!             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
     448!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
     449!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
     450!             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
     451!             ! Vertical scale factor interpolations 
     452!             ! ------------------------------------ 
     453!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
     454!             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
     455!             CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
     456!             CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     457!             CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
     458!             ! t- and w- points depth 
     459!             ! ---------------------- 
     460!             !!gm not sure of that.... 
     461!             gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
     462!             gdepw(:,:,1,Kmm) = 0.0_wp 
     463!             gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
     464!             DO jk = 2, jpk 
     465!                gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk  ,Kmm) 
     466!                gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
     467!                gde3w(:,:,jk) = gdept(:,:,jk  ,Kmm) - ssh (:,:,Kmm) 
     468!             END DO 
     469!          ENDIF 
    468470      ENDIF 
    469471       
     
    503505      !! 
    504506      !!----------------------------------------------------------------------------- 
    505       INTEGER ::   ios   ! Local integer output status for namelist read 
    506       INTEGER ::   ifpr, ierror 
     507      INTEGER ::   ios, ifpr, ierror   ! Local integers 
     508 
    507509      ! 
    508510      CHARACTER(len=256) ::  cn_dir          ! Root directory for location of ice files 
  • NEMO/trunk/src/OCE/ASM/asminc.F90

    r13226 r13237  
    9595   !! * Substitutions 
    9696#  include "do_loop_substitute.h90" 
     97#  include "domzgr_substitute.h90" 
    9798   !!---------------------------------------------------------------------- 
    9899   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    417418                     &            - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * u_bkginc(ji-1,jj,jk)    & 
    418419                     &            + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * v_bkginc(ji,jj  ,jk)    & 
    419                      &            - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * v_bkginc(ji,jj-1,jk)  ) / e3t(ji,jj,jk,Kmm) 
     420                     &            - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * v_bkginc(ji,jj-1,jk)  ) & 
     421                     &            / e3t(ji,jj,jk,Kmm) 
    420422               END_2D 
    421423               CALL lbc_lnk( 'asminc', zhdiv, 'T', 1.0_wp )   ! lateral boundary cond. (no sign change) 
     
    758760            ! 
    759761            ssh(:,:,Kbb) = ssh(:,:,Kmm)                        ! Update before fields 
     762#if ! defined key_qco 
    760763            e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
     764#endif 
    761765!!gm why not e3u(:,:,:,Kbb), e3v(:,:,:,Kbb), gdept(:,:,:,Kbb) ???? 
    762766            ! 
     
    970974!           ! set to bottom of a level  
    971975!                 DO jk = jpk-1, 2, -1 
    972 !                   IF ((mld > gdepw(ji,jj,jk)) .and. (mld < gdepw(ji,jj,jk+1))) THEN  
    973 !                     mld=gdepw(ji,jj,jk+1) 
     976!                   IF ((mld > gdepw(ji,jj,jk,Kmm)) .and. (mld < gdepw(ji,jj,jk+1,Kmm))) THEN 
     977!                     mld=gdepw(ji,jj,jk+1,Kmm) 
    974978!                     jkmax=jk 
    975979!                   ENDIF 
  • NEMO/trunk/src/OCE/BDY/bdydta.F90

    r12921 r13237  
    7070   !! * Substitutions 
    7171#  include "do_loop_substitute.h90" 
     72#  include "domzgr_substitute.h90" 
    7273   !!---------------------------------------------------------------------- 
    7374   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    247248               dta_alias%u2d(ib) = 0._wp   ! compute barotrope zonal velocity and put it in u2d 
    248249               DO ik = 1, jpkm1 
    249                   dta_alias%u2d(ib) = dta_alias%u2d(ib) + e3u(ii,ij,ik,Kmm) * umask(ii,ij,ik) * dta_alias%u3d(ib,ik) 
     250                  dta_alias%u2d(ib) = dta_alias%u2d(ib)   & 
     251                     &              + e3u(ii,ij,ik,Kmm) * umask(ii,ij,ik) * dta_alias%u3d(ib,ik) 
    250252               END DO 
    251253               dta_alias%u2d(ib) =  dta_alias%u2d(ib) * r1_hu(ii,ij,Kmm) 
     
    260262               dta_alias%v2d(ib) = 0._wp   ! compute barotrope meridional velocity and put it in v2d 
    261263               DO ik = 1, jpkm1 
    262                   dta_alias%v2d(ib) = dta_alias%v2d(ib) + e3v(ii,ij,ik,Kmm) * vmask(ii,ij,ik) * dta_alias%v3d(ib,ik) 
     264                  dta_alias%v2d(ib) = dta_alias%v2d(ib)   & 
     265                     &              + e3v(ii,ij,ik,Kmm) * vmask(ii,ij,ik) * dta_alias%v3d(ib,ik) 
    263266               END DO 
    264267               dta_alias%v2d(ib) =  dta_alias%v2d(ib) * r1_hv(ii,ij,Kmm) 
  • NEMO/trunk/src/OCE/BDY/bdydyn.F90

    r12377 r13237  
    3030   PUBLIC   bdy_dyn    ! routine called in dyn_nxt 
    3131 
     32   !! * Substitutions 
     33#  include "domzgr_substitute.h90" 
    3234   !!---------------------------------------------------------------------- 
    3335   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/C1D/step_c1d.F90

    r12933 r13237  
    8080      IF(.NOT.ln_linssh )   CALL dom_vvl_sf_nxt( kstp, Nbb, Nnn,      Naa )  ! after vertical scale factors  
    8181 
    82       IF(.NOT.ln_linssh )   CALL wzv           ( kstp, Nbb, Nnn, ww,  Naa )  ! now cross-level velocity  
     82      IF(.NOT.ln_linssh )   CALL wzv           ( kstp, Nbb, Nnn, Naa, ww )  ! now cross-level velocity  
    8383      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    8484      ! diagnostics and outputs        
  • NEMO/trunk/src/OCE/CRS/crsfld.F90

    r13226 r13237  
    3333   !! * Substitutions 
    3434#  include "do_loop_substitute.h90" 
     35#  include "domzgr_substitute.h90" 
    3536   !!---------------------------------------------------------------------- 
    3637   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6869 
    6970      ! Depth work arrrays 
    70       ze3t(:,:,:) = e3t(:,:,:,Kmm) 
    71       ze3u(:,:,:) = e3u(:,:,:,Kmm) 
    72       ze3v(:,:,:) = e3v(:,:,:,Kmm) 
    73       ze3w(:,:,:) = e3w(:,:,:,Kmm) 
     71      DO jk = 1 , jpk  
     72         ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     73         ze3u(:,:,jk) = e3u(:,:,jk,Kmm) 
     74         ze3v(:,:,jk) = e3v(:,:,jk,Kmm) 
     75         ze3w(:,:,jk) = e3w(:,:,jk,Kmm) 
     76      END DO 
    7477 
    7578      IF( kt == nit000  ) THEN 
  • NEMO/trunk/src/OCE/CRS/crsini.F90

    r13226 r13237  
    2828   PUBLIC   crs_init   ! called by nemogcm.F90 module 
    2929 
     30   !! * Substitutions 
     31#  include "domzgr_substitute.h90" 
    3032   !!---------------------------------------------------------------------- 
    3133   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    174176      
    175177     ! 
    176      ze3t(:,:,:) = e3t(:,:,:,Kmm) 
    177      ze3u(:,:,:) = e3u(:,:,:,Kmm) 
    178      ze3v(:,:,:) = e3v(:,:,:,Kmm) 
    179      ze3w(:,:,:) = e3w(:,:,:,Kmm) 
     178     DO jk = 1, jpk 
     179        ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     180        ze3u(:,:,jk) = e3u(:,:,jk,Kmm) 
     181        ze3v(:,:,jk) = e3v(:,:,jk,Kmm) 
     182        ze3w(:,:,jk) = e3w(:,:,jk,Kmm) 
     183     END DO   
    180184 
    181185     !    3.d.2   Surfaces  
  • NEMO/trunk/src/OCE/DIA/diaar5.F90

    r13226 r13237  
    3939   !! * Substitutions 
    4040#  include "do_loop_substitute.h90" 
     41#  include "domzgr_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    4243   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    7677      ! 
    7778      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: zarea_ssh , zbotpres       ! 2D workspace  
    78       REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: zpe, z2d                   ! 2D workspace  
    79       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   :: zrhd , ztpot               ! 3D workspace 
     79      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     :: z2d, zpe                   ! 2D workspace  
     80      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   :: z3d, zrhd, ztpot, zgdept   ! 3D workspace (zgdept: needed to use the substitute) 
    8081      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn                       ! 4D workspace 
    8182 
     
    101102            zrhd(:,:,jk) = e1e2t(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    102103         END DO 
     104         DO jk = 1, jpk 
     105            z3d(:,:,jk) =  rho0 * e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     106         END DO  
    103107         CALL iom_put( 'volcello'  , zrhd(:,:,:)  )  ! WARNING not consistent with CMIP DR where volcello is at ca. 2000 
    104          CALL iom_put( 'masscello' , rho0 * e3t(:,:,:,Kmm) * tmask(:,:,:) )  ! ocean mass 
     108         CALL iom_put( 'masscello' , z3d (:,:,:) )   ! ocean mass 
    105109      ENDIF  
    106110      ! 
     
    128132         ztsn(:,:,:,jp_tem) = ts(:,:,:,jp_tem,Kmm)                    ! thermosteric ssh 
    129133         ztsn(:,:,:,jp_sal) = sn0(:,:,:) 
    130          CALL eos( ztsn, zrhd, gdept(:,:,:,Kmm) )                       ! now in situ density using initial salinity 
     134         ALLOCATE( zgdept(jpi,jpj,jpk) ) 
     135         DO jk = 1, jpk 
     136            zgdept(:,:,jk) = gdept(:,:,jk,Kmm) 
     137         END DO 
     138         CALL eos( ztsn, zrhd, zgdept)                       ! now in situ density using initial salinity 
    131139         ! 
    132140         zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice 
     
    180188         CALL iom_put( 'botpres', zbotpres ) 
    181189         ! 
     190         DEALLOCATE( zgdept ) 
     191         ! 
    182192      ENDIF 
    183193 
  • NEMO/trunk/src/OCE/DIA/diacfl.F90

    r12489 r13237  
    3434   !! * Substitutions 
    3535#  include "do_loop_substitute.h90" 
     36#  include "domzgr_substitute.h90" 
    3637   !!---------------------------------------------------------------------- 
    3738   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/DIA/diadct.F90

    r12489 r13237  
    1111   !!            3.4  ! 09/2011 (C Bricaud) 
    1212   !!---------------------------------------------------------------------- 
    13    !! does not work with agrif 
    1413#if ! defined key_agrif 
     14   !!                        ==>>  CAUTION: does not work with agrif 
    1515   !!---------------------------------------------------------------------- 
    1616   !!   dia_dct      :  Compute the transport through a sec. 
     
    6666   TYPE SECTION 
    6767      CHARACTER(len=60)                            :: name              ! name of the sec 
    68       LOGICAL                                      :: llstrpond         ! true if you want the computation of salt and 
    69                                                                        ! heat transports 
     68      LOGICAL                                      :: llstrpond         ! true if you want the computation of salt and heat transports 
    7069      LOGICAL                                      :: ll_ice_section    ! ice surface and ice volume computation 
    7170      LOGICAL                                      :: ll_date_line      ! = T if the section crosses the date-line 
     
    7473      INTEGER, DIMENSION(nb_point_max)             :: direction         ! vector direction of the point in the section 
    7574      CHARACTER(len=40),DIMENSION(nb_class_max)    :: classname         ! characteristics of the class 
    76       REAL(wp), DIMENSION(nb_class_max)            :: zsigi           ,&! in-situ   density classes    (99 if you don't want) 
    77                                                       zsigp           ,&! potential density classes    (99 if you don't want) 
    78                                                       zsal            ,&! salinity classes   (99 if you don't want) 
    79                                                       ztem            ,&! temperature classes(99 if you don't want) 
    80                                                       zlay              ! level classes      (99 if you don't want) 
     75      REAL(wp), DIMENSION(nb_class_max)            :: zsigi             ! in-situ   density classes    (99 if you don't want) 
     76      REAL(wp), DIMENSION(nb_class_max)            :: zsigp             ! potential density classes    (99 if you don't want) 
     77      REAL(wp), DIMENSION(nb_class_max)            :: zsal              ! salinity classes   (99 if you don't want) 
     78      REAL(wp), DIMENSION(nb_class_max)            :: ztem              ! temperature classes(99 if you don't want) 
     79      REAL(wp), DIMENSION(nb_class_max)            :: zlay              ! level classes      (99 if you don't want) 
    8180      REAL(wp), DIMENSION(nb_type_class,nb_class_max)  :: transport     ! transport output 
    8281      REAL(wp)                                         :: slopeSection  ! slope of the section 
     
    9089   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::  transports_2d   
    9190 
     91 
     92   !! * Substitutions 
     93#  include "domzgr_substitute.h90" 
    9294   !!---------------------------------------------------------------------- 
    9395   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9597   !! Software governed by the CeCILL license (see ./LICENSE) 
    9698   !!---------------------------------------------------------------------- 
     99 
    97100CONTAINS 
    98101  
     
    11191122  !!    |               |                  |       interpolation between ptab(I,J,K) and ptab(I,J,K+1) 
    11201123  !!    |               |                  |       zbis =  
    1121   !!    |               |                  |      [ e3w(I+1,J,K)*ptab(I,J,K) + ( e3w(I,J,K) - e3w(I+1,J,K) ) * ptab(I,J,K-1) ] 
    1122   !!    |               |                  |      /[ e3w(I+1,J,K) + e3w(I,J,K) - e3w(I+1,J,K) ]  
     1124  !!    |               |                  |      [ e3w_n(I+1,J,K,NOW)*ptab(I,J,K) + ( e3w_n(I,J,K,NOW) - e3w_n(I+1,J,K,NOW) ) * ptab(I,J,K-1) ] 
     1125  !!    |               |                  |     /[ e3w_n(I+1,J,K,NOW)             +   e3w_n(I,J,K,NOW) - e3w_n(I+1,J,K,NOW) ]  
    11231126  !!    |               |                  |  
    11241127  !!    |               |                  |    2. Horizontal interpolation: compute value at U/V point 
     
    12131216 
    12141217     ze3t  = e3t(ii2,ij2,kk,Kmm) - e3t(ii1,ij1,kk,Kmm)  
    1215      zwgt1 = ( e3w(ii2,ij2,kk,Kmm) - e3w(ii1,ij1,kk,Kmm) ) / e3w(ii2,ij2,kk,Kmm) 
    1216      zwgt2 = ( e3w(ii1,ij1,kk,Kmm) - e3w(ii2,ij2,kk,Kmm) ) / e3w(ii1,ij1,kk,Kmm) 
     1218     zwgt1 = ( e3w(ii2,ij2,kk,Kmm) - e3w(ii1,ij1,kk,Kmm) )   & 
     1219        &    / e3w(ii2,ij2,kk,Kmm) 
     1220     zwgt2 = ( e3w(ii1,ij1,kk,Kmm) - e3w(ii2,ij2,kk,Kmm) )   & 
     1221        &    / e3w(ii1,ij1,kk,Kmm) 
    12171222 
    12181223     IF(kk .NE. 1)THEN 
  • NEMO/trunk/src/OCE/DIA/diahsb.F90

    r12489 r13237  
    5050   REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   tmask_ini 
    5151 
     52   !! * Substitutions 
     53#  include "domzgr_substitute.h90" 
    5254   !!---------------------------------------------------------------------- 
    5355   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    156158      ! 
    157159      DO jk = 1, jpkm1           ! volume variation (calculated with scale factors) 
    158          zwrk(:,:,jk) = surf(:,:)*e3t(:,:,jk,Kmm)*tmask(:,:,jk) - surf_ini(:,:)*e3t_ini(:,:,jk)*tmask_ini(:,:,jk) 
     160         zwrk(:,:,jk) =   surf    (:,:) * e3t    (:,:,jk,Kmm)*tmask    (:,:,jk)   & 
     161            &           - surf_ini(:,:) * e3t_ini(:,:,jk    )*tmask_ini(:,:,jk) 
    159162      END DO 
    160163      zdiff_v2 = glob_sum_full( 'diahsb', zwrk(:,:,:) )     ! glob_sum_full needed as tmask and tmask_ini could be different 
    161164      DO jk = 1, jpkm1           ! heat content variation 
    162          zwrk(:,:,jk) = ( surf(:,:)*e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_tem,Kmm) - surf_ini(:,:)*hc_loc_ini(:,:,jk) ) 
     165         zwrk(:,:,jk) = ( surf    (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_tem,Kmm)   & 
     166            &           - surf_ini(:,:) *         hc_loc_ini(:,:,jk) ) 
    163167      END DO 
    164168      zdiff_hc = glob_sum_full( 'diahsb', zwrk(:,:,:) ) 
    165169      DO jk = 1, jpkm1           ! salt content variation 
    166          zwrk(:,:,jk) = ( surf(:,:)*e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_sal,Kmm) - surf_ini(:,:)*sc_loc_ini(:,:,jk) ) 
     170         zwrk(:,:,jk) = ( surf    (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_sal,Kmm)   & 
     171            &           - surf_ini(:,:) *         sc_loc_ini(:,:,jk) ) 
    167172      END DO 
    168173      zdiff_sc = glob_sum_full( 'diahsb', zwrk(:,:,:) ) 
  • NEMO/trunk/src/OCE/DIA/diahth.F90

    r12489 r13237  
    4242   !! * Substitutions 
    4343#  include "do_loop_substitute.h90" 
     44#  include "domzgr_substitute.h90" 
    4445   !!---------------------------------------------------------------------- 
    4546   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    361362         ik = ilevel(ji,jj) 
    362363         zthick(ji,jj) = pdep - zthick(ji,jj)   !   remaining thickness to reach depht pdep 
    363          phtc(ji,jj)   = phtc(ji,jj) + pt(ji,jj,ik+1) * MIN( e3t(ji,jj,ik+1,Kmm), zthick(ji,jj) ) & 
     364         phtc(ji,jj)   = phtc(ji,jj)    & 
     365            &           + pt (ji,jj,ik+1) * MIN( e3t(ji,jj,ik+1,Kmm), zthick(ji,jj) ) & 
    364366                                                       * tmask(ji,jj,ik+1) 
    365367      END_2D 
  • NEMO/trunk/src/OCE/DIA/diamlr.F90

    r12642 r13237  
    44   !! Management of the IOM context for multiple-linear-regression analysis 
    55   !!====================================================================== 
    6    !! History :       !  2019  (S. Mueller) 
     6   !! History :  4.0  !  2019  (S. Mueller)   Original code 
    77   !!---------------------------------------------------------------------- 
    88 
    99   USE par_oce        , ONLY :   wp, jpi, jpj 
    1010   USE phycst         , ONLY :   rpi 
     11   USE dom_oce        , ONLY :   adatrj 
     12   USE tide_mod 
     13   ! 
    1114   USE in_out_manager , ONLY :   lwp, numout, ln_timing 
    1215   USE iom            , ONLY :   iom_put, iom_use, iom_update_file_name 
    13    USE dom_oce        , ONLY :   adatrj 
    1416   USE timing         , ONLY :   timing_start, timing_stop 
    1517#if defined key_iomput 
    1618   USE xios 
    1719#endif 
    18    USE tide_mod 
    1920 
    2021   IMPLICIT NONE 
    2122   PRIVATE 
    2223 
    23    LOGICAL, PUBLIC ::   lk_diamlr = .FALSE. 
     24   LOGICAL, PUBLIC ::   lk_diamlr = .FALSE.   !:         ===>>>   NOT a DOCTOR norm name :  use l_diamlr 
     25   !                                                              lk_  is used only for logical controlled by a CPP key 
    2426 
    2527   PUBLIC ::   dia_mlr_init, dia_mlr_iom_init, dia_mlr 
     
    4244      !! 
    4345      !!---------------------------------------------------------------------- 
    44  
     46      ! 
    4547      lk_diamlr = .TRUE. 
    46  
     48      ! 
    4749      IF(lwp) THEN 
    4850         WRITE(numout, *) 
     
    5052         WRITE(numout, *) '~~~~~~~~~~~~   multiple-linear-regression analysis' 
    5153      END IF 
    52  
     54      ! 
    5355   END SUBROUTINE dia_mlr_init 
     56 
    5457 
    5558   SUBROUTINE dia_mlr_iom_init 
     
    397400   END SUBROUTINE dia_mlr_iom_init 
    398401 
     402 
    399403   SUBROUTINE dia_mlr 
    400404      !!---------------------------------------------------------------------- 
     
    404408      !! 
    405409      !!---------------------------------------------------------------------- 
    406  
    407410      REAL(wp), DIMENSION(jpi,jpj) ::   zadatrj2d 
     411      !!---------------------------------------------------------------------- 
    408412 
    409413      IF( ln_timing )   CALL timing_start('dia_mlr') 
     
    412416      ! (value of adatrj converted to time in units of seconds) 
    413417      ! 
    414       ! A 2-dimensional field of constant value is sent, and subsequently used 
    415       ! directly or transformed to a scalar or a constant 3-dimensional field as 
    416       ! required. 
     418      ! A 2-dimensional field of constant value is sent, and subsequently used directly  
     419      ! or transformed to a scalar or a constant 3-dimensional field as required. 
    417420      zadatrj2d(:,:) = adatrj*86400.0_wp 
    418421      IF ( iom_use('diamlr_time') ) CALL iom_put('diamlr_time', zadatrj2d) 
    419        
     422      ! 
    420423      IF( ln_timing )   CALL timing_stop('dia_mlr') 
    421  
     424      ! 
    422425   END SUBROUTINE dia_mlr 
    423426 
     427   !!====================================================================== 
    424428END MODULE diamlr 
  • NEMO/trunk/src/OCE/DIA/diaptr.F90

    r13226 r13237  
    6060 
    6161   LOGICAL ::   ll_init = .TRUE.        !: tracers  trend flag (set from namelist in trdini) 
     62    
    6263   !! * Substitutions 
    6364#  include "do_loop_substitute.h90" 
     65#  include "domzgr_substitute.h90" 
    6466   !!---------------------------------------------------------------------- 
    6567   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/DIA/diawri.F90

    r13226 r13237  
    8585   !! * Substitutions 
    8686#  include "do_loop_substitute.h90" 
     87#  include "domzgr_substitute.h90" 
    8788   !!---------------------------------------------------------------------- 
    8889   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    136137      CALL iom_put("e3v_0", e3v_0(:,:,:) ) 
    137138      ! 
    138       CALL iom_put( "e3t" , e3t(:,:,:,Kmm) ) 
    139       CALL iom_put( "e3u" , e3u(:,:,:,Kmm) ) 
    140       CALL iom_put( "e3v" , e3v(:,:,:,Kmm) ) 
    141       CALL iom_put( "e3w" , e3w(:,:,:,Kmm) ) 
    142       IF( iom_use("e3tdef") )   & 
    143          CALL iom_put( "e3tdef"  , ( ( e3t(:,:,:,Kmm) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 
    144  
    145       IF( ll_wd ) THEN 
    146          CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)*tmask(:,:,1) )   ! sea surface height (brought back to the reference used for wetting and drying) 
     139      IF ( iom_use("e3t") .OR. iom_use("e3tdef") ) THEN  ! time-varying e3t 
     140         DO jk = 1, jpk 
     141            z3d(:,:,jk) =  e3t(:,:,jk,Kmm) 
     142         END DO 
     143         CALL iom_put( "e3t"     ,     z3d(:,:,:) ) 
     144         CALL iom_put( "e3tdef"  , ( ( z3d(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 )  
     145      ENDIF  
     146      IF ( iom_use("e3u") ) THEN                         ! time-varying e3u 
     147         DO jk = 1, jpk 
     148            z3d(:,:,jk) =  e3u(:,:,jk,Kmm) 
     149         END DO  
     150         CALL iom_put( "e3u" , z3d(:,:,:) ) 
     151      ENDIF 
     152      IF ( iom_use("e3v") ) THEN                         ! time-varying e3v 
     153         DO jk = 1, jpk 
     154            z3d(:,:,jk) =  e3v(:,:,jk,Kmm) 
     155         END DO  
     156         CALL iom_put( "e3v" , z3d(:,:,:) ) 
     157      ENDIF 
     158      IF ( iom_use("e3w") ) THEN                         ! time-varying e3w 
     159         DO jk = 1, jpk 
     160            z3d(:,:,jk) =  e3w(:,:,jk,Kmm) 
     161         END DO  
     162         CALL iom_put( "e3w" , z3d(:,:,:) ) 
     163      ENDIF 
     164 
     165      IF( ll_wd ) THEN                                   ! sea surface height (brought back to the reference used for wetting and drying) 
     166         CALL iom_put( "ssh" , (ssh(:,:,Kmm)+ssh_ref)*tmask(:,:,1) ) 
    147167      ELSE 
    148168         CALL iom_put( "ssh" , ssh(:,:,Kmm) )              ! sea surface height 
     
    172192      ENDIF 
    173193 
     194#if ! defined key_qco 
    174195      CALL iom_put( "rhop", rhop(:,:,:) )          ! 3D potential density (sigma0) 
     196#endif 
    175197 
    176198      IF ( iom_use("taubot") ) THEN                ! bottom stress 
     
    210232 
    211233      IF( ln_zad_Aimp ) ww = ww + wi               ! Recombine explicit and implicit parts of vertical velocity for diagnostic output 
    212       ! 
    213234      CALL iom_put( "woce", ww )                   ! vertical velocity 
     235 
    214236      IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN   ! vertical mass transport & its square value 
    215237         ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 
     
    417439      ! 
    418440      REAL(wp), DIMENSION(jpi,jpj)   :: zw2d       ! 2D workspace 
    419       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d       ! 3D workspace 
     441      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zw3d, ze3t, zgdept       ! 3D workspace 
    420442      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d_abl   ! ABL 3D workspace 
    421443      !!---------------------------------------------------------------------- 
     
    457479      it = kt 
    458480      itmod = kt - nit000 + 1 
     481 
     482      ! store e3t for subsitute 
     483      DO jk = 1, jpk 
     484         ze3t  (:,:,jk) =  e3t  (:,:,jk,Kmm) 
     485         zgdept(:,:,jk) =  gdept(:,:,jk,Kmm) 
     486      END DO 
    459487 
    460488 
     
    571599         DEALLOCATE(zw3d_abl) 
    572600         ENDIF 
     601         ! 
    573602 
    574603         ! Declare all the output fields as NETCDF variables 
     
    580609            &          jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    581610         IF(  .NOT.ln_linssh  ) THEN 
    582             CALL histdef( nid_T, "vovvle3t", "Level thickness"                    , "m"      ,&  ! e3t(:,:,:,Kmm) 
     611            CALL histdef( nid_T, "vovvle3t", "Level thickness"                    , "m"      ,&  ! e3t n 
    583612            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    584             CALL histdef( nid_T, "vovvldep", "T point depth"                      , "m"      ,&  ! e3t(:,:,:,Kmm) 
     613            CALL histdef( nid_T, "vovvldep", "T point depth"                      , "m"      ,&  ! e3t n 
    585614            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    586             CALL histdef( nid_T, "vovvldef", "Squared level deformation"          , "%^2"    ,&  ! e3t(:,:,:,Kmm) 
     615            CALL histdef( nid_T, "vovvldef", "Squared level deformation"          , "%^2"    ,&  ! e3t n 
    587616            &             jpi, jpj, nh_T, ipk, 1, ipk, nz_T, 32, clop, zsto, zout ) 
    588617         ENDIF 
     
    768797 
    769798      IF( .NOT.ln_linssh ) THEN 
    770          CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) * e3t(:,:,:,Kmm) , ndim_T , ndex_T  )   ! heat content 
    771          CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) * e3t(:,:,:,Kmm) , ndim_T , ndex_T  )   ! salt content 
    772          CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) * e3t(:,:,1,Kmm) , ndim_hT, ndex_hT )   ! sea surface heat content 
    773          CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) * e3t(:,:,1,Kmm) , ndim_hT, ndex_hT )   ! sea surface salinity content 
     799         CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) * ze3t(:,:,:) , ndim_T , ndex_T  )   ! heat content 
     800         CALL histwrite( nid_T, "vosaline", it, ts(:,:,:,jp_sal,Kmm) * ze3t(:,:,:) , ndim_T , ndex_T  )   ! salt content 
     801         CALL histwrite( nid_T, "sosstsst", it, ts(:,:,1,jp_tem,Kmm) * ze3t(:,:,1) , ndim_hT, ndex_hT )   ! sea surface heat content 
     802         CALL histwrite( nid_T, "sosaline", it, ts(:,:,1,jp_sal,Kmm) * ze3t(:,:,1) , ndim_hT, ndex_hT )   ! sea surface salinity content 
    774803      ELSE 
    775804         CALL histwrite( nid_T, "votemper", it, ts(:,:,:,jp_tem,Kmm) , ndim_T , ndex_T  )   ! temperature 
     
    779808      ENDIF 
    780809      IF( .NOT.ln_linssh ) THEN 
    781          zw3d(:,:,:) = ( ( e3t(:,:,:,Kmm) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 
    782          CALL histwrite( nid_T, "vovvle3t", it, e3t (:,:,:,Kmm) , ndim_T , ndex_T  )   ! level thickness 
    783          CALL histwrite( nid_T, "vovvldep", it, gdept(:,:,:,Kmm) , ndim_T , ndex_T  )   ! t-point depth 
     810         zw3d(:,:,:) = ( ( ze3t(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 
     811         CALL histwrite( nid_T, "vovvle3t", it, ze3t (:,:,:)    , ndim_T , ndex_T  )   ! level thickness 
     812         CALL histwrite( nid_T, "vovvldep", it, zgdept , ndim_T , ndex_T  )   ! t-point depth  
    784813         CALL histwrite( nid_T, "vovvldef", it, zw3d             , ndim_T , ndex_T  )   ! level thickness deformation 
    785814      ENDIF 
     
    920949      !! 
    921950      INTEGER :: inum, jk 
     951      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, zgdept      ! 3D workspace !!st patch to use substitution 
    922952      !!---------------------------------------------------------------------- 
    923953      !  
    924       IF(lwp) WRITE(numout,*) 
    925       IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state' 
    926       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created ' 
    927       IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '...nc' 
     954      IF(lwp) THEN 
     955         WRITE(numout,*) 
     956         WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state' 
     957         WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created ' 
     958         WRITE(numout,*) '                and named :', cdfile_name, '...nc' 
     959      ENDIF  
     960      ! 
     961      DO jk = 1, jpk 
     962         ze3t(:,:,jk) =  e3t(:,:,jk,Kmm) 
     963         zgdept(:,:,jk) =  gdept(:,:,jk,Kmm) 
     964      END DO 
    928965      ! 
    929966      CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE. ) 
     
    940977      ENDIF 
    941978      CALL iom_rstput( 0, 0, inum, 'risfdep', risfdep            )    ! now k-velocity 
    942       CALL iom_rstput( 0, 0, inum, 'ht'     , ht                 )    ! now water column height 
     979      CALL iom_rstput( 0, 0, inum, 'ht'     , ht(:,:)            )    ! now water column height 
    943980      ! 
    944981      IF ( ln_isf ) THEN 
     
    9771014      CALL iom_rstput( 0, 0, inum, 'sometauy', vtau              )    ! j-wind stress 
    9781015      IF(  .NOT.ln_linssh  ) THEN              
    979          CALL iom_rstput( 0, 0, inum, 'vovvldep', gdept(:,:,:,Kmm)        )    !  T-cell depth  
    980          CALL iom_rstput( 0, 0, inum, 'vovvle3t', e3t(:,:,:,Kmm)          )    !  T-cell thickness   
     1016         CALL iom_rstput( 0, 0, inum, 'vovvldep', zgdept        )    !  T-cell depth  
     1017         CALL iom_rstput( 0, 0, inum, 'vovvle3t', ze3t          )    !  T-cell thickness   
    9811018      END IF 
    9821019      IF( ln_wave .AND. ln_sdw ) THEN 
  • NEMO/trunk/src/OCE/DOM/dom_oce.F90

    r13216 r13237  
    22   !!====================================================================== 
    33   !!                       ***  MODULE dom_oce  *** 
    4    !!        
    54   !! ** Purpose :   Define in memory all the ocean space domain variables 
    65   !!====================================================================== 
     
    1312   !!             -   ! 2015-11  (G. Madec, A. Coward)  time varying zgr by default 
    1413   !!            4.1  ! 2019-08  (A. Coward, D. Storkey) rename prognostic variables in preparation for new time scheme. 
     14   !!            4.x  ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio 
    1515   !!---------------------------------------------------------------------- 
    1616 
     
    7474   LOGICAL, PUBLIC ::   l_Iperio, l_Jperio   !   should we explicitely take care I/J periodicity  
    7575 
    76    !                                 ! domain MPP decomposition parameters 
     76   !                             !: domain MPP decomposition parameters 
    7777   INTEGER             , PUBLIC ::   nimpp, njmpp     !: i- & j-indexes for mpp-subdomain left bottom 
    7878   INTEGER             , PUBLIC ::   nreci, nrecj     !: overlap region in i and j 
     
    137137   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::    e3vw_0   !: vw-vert. scale factor [m] 
    138138   !                                                        !  time-dependent scale factors 
     139#if ! defined key_qco 
    139140   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   e3t, e3u, e3v, e3w, e3uw, e3vw  !: vert. scale factor [m] 
    140141   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   e3f                             !: F-point vert. scale factor [m] 
     142#endif 
     143   !                                                        !  time-dependent ratio ssh / h_0 
     144   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   r3t, r3u, r3v                   !: time-dependent    ratio at t-, u- and v-point [-] 
     145   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   r3f                             !: mid-time-level    ratio at f-point            [-] 
     146   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   r3t_f, r3u_f, r3v_f             !: now time-filtered ratio at t-, u- and v-point [-] 
    141147 
    142148   !                                                        !  reference depths of cells 
     
    148154   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::  gde3w   
    149155    
    150    !                                                      !  reference heights of water column 
    151    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_0  !: t-depth              [m] 
    152    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hu_0  !: u-depth              [m] 
    153    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hv_0  !: v-depth              [m] 
    154                                                           ! time-dependent heights of water column 
    155    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ht                     !: height of water column at T-points [m] 
    156    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hu, hv, r1_hu, r1_hv   !: height of water column [m] and reciprocal [1/m] 
     156   !                                                        !  reference heights of ocean water column and its inverse 
     157   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   ht_0, r1_ht_0   !: t-depth        [m] and [1/m] 
     158   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   hu_0, r1_hu_0   !: u-depth        [m] and [1/m] 
     159   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   hv_0, r1_hv_0   !: v-depth        [m] and [1/m] 
     160   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   hf_0, r1_hf_0   !: f-depth        [m] and [1/m] 
     161   !                                                        ! time-dependent heights of ocean water column 
     162#if ! defined key_qco 
     163   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   ht          !: t-points           [m] 
     164#endif 
     165   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hu, r1_hu   !: u-depth            [m] and [1/m] 
     166   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hv, r1_hv   !: v-depth            [m] and [1/m] 
    157167 
    158168   INTEGER, PUBLIC ::   nla10              !: deepest    W level Above  ~10m (nlb10 - 1) 
     
    170180   !! --------------------------------------------------------------------- 
    171181!!gm Proposition of new name for top/bottom vertical indices 
    172 !   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mtk_t, mtk_u, mtk_v   !: top first wet T-, U-, V-, F-level (ISF) 
    173 !   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mbk_t, mbk_u, mbk_v   !: bottom last wet T-, U- and V-level 
     182!   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mtk_t, mtk_u, mtk_v   !: top    first wet T-, U-, and V-level (ISF) 
     183!   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mbk_t, mbk_u, mbk_v   !: bottom last  wet T-, U-, and V-level 
    174184!!gm 
    175185   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mbkt, mbku, mbkv   !: bottom last wet T-, U- and V-level 
     
    179189   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   mikt, miku, mikv, mikf  !: top first wet T-, U-, V-, F-level           (ISF) 
    180190 
    181    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ssmask, ssumask, ssvmask             !: surface mask at T-,U-, V- and F-pts 
     191   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)           ::   ssmask, ssumask, ssvmask, ssfmask   !: surface mask at T-,U-, V- and F-pts 
    182192   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: tmask, umask, vmask, fmask   !: land/ocean mask at T-, U-, V- and F-pts 
    183193   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: wmask, wumask, wvmask        !: land/ocean mask at WT-, WU- and WV-pts 
     
    246256   INTEGER FUNCTION dom_oce_alloc() 
    247257      !!---------------------------------------------------------------------- 
    248       INTEGER, DIMENSION(12) :: ierr 
     258      INTEGER                ::   ii 
     259      INTEGER, DIMENSION(30) :: ierr 
    249260      !!---------------------------------------------------------------------- 
    250       ierr(:) = 0 
     261      ii = 0   ;   ierr(:) = 0 
    251262      ! 
    252       ALLOCATE( mig(jpi), mjg(jpj), STAT=ierr(1) ) 
    253          ! 
     263      ii = ii+1 
     264      ALLOCATE( mig(jpi), mjg(jpj), STAT=ierr(ii) ) 
     265         ! 
     266      ii = ii+1 
    254267      ALLOCATE( mi0(jpiglo)   , mi1 (jpiglo),  mj0(jpjglo)   , mj1 (jpjglo) ,     & 
    255          &      tpol(jpiglo)  , fpol(jpiglo)                                , STAT=ierr(2) ) 
    256          ! 
     268         &      tpol(jpiglo) , fpol(jpiglo)                              , STAT=ierr(ii) ) 
     269         ! 
     270      ii = ii+1 
    257271      ALLOCATE( glamt(jpi,jpj) ,    glamu(jpi,jpj) ,  glamv(jpi,jpj) ,  glamf(jpi,jpj) ,     & 
    258272         &      gphit(jpi,jpj) ,    gphiu(jpi,jpj) ,  gphiv(jpi,jpj) ,  gphif(jpi,jpj) ,     & 
     
    265279         &      e1e2v(jpi,jpj) , r1_e1e2v(jpi,jpj) , e1_e2v(jpi,jpj)                   ,     & 
    266280         &      e1e2f(jpi,jpj) , r1_e1e2f(jpi,jpj)                                     ,     & 
    267          &      ff_f (jpi,jpj) ,    ff_t (jpi,jpj)                                     , STAT=ierr(3) ) 
    268          ! 
     281         &      ff_f (jpi,jpj) ,    ff_t (jpi,jpj)                                     , STAT=ierr(ii) ) 
     282         ! 
     283      ii = ii+1 
    269284      ALLOCATE( gdept_0(jpi,jpj,jpk)     , gdepw_0(jpi,jpj,jpk)     , gde3w_0(jpi,jpj,jpk) ,      & 
    270          &      gdept  (jpi,jpj,jpk,jpt) , gdepw  (jpi,jpj,jpk,jpt) , gde3w  (jpi,jpj,jpk) , STAT=ierr(4) ) 
    271          ! 
    272       ALLOCATE( e3t_0(jpi,jpj,jpk)     , e3u_0(jpi,jpj,jpk)     , e3v_0(jpi,jpj,jpk)     , e3f_0(jpi,jpj,jpk) , e3w_0(jpi,jpj,jpk)     ,   & 
    273          &      e3t  (jpi,jpj,jpk,jpt) , e3u  (jpi,jpj,jpk,jpt) , e3v  (jpi,jpj,jpk,jpt) , e3f  (jpi,jpj,jpk) , e3w  (jpi,jpj,jpk,jpt) ,   &  
    274          &      e3uw_0(jpi,jpj,jpk)     , e3vw_0(jpi,jpj,jpk)     ,         & 
    275          &      e3uw  (jpi,jpj,jpk,jpt) , e3vw  (jpi,jpj,jpk,jpt) ,    STAT=ierr(5) )                        
    276          ! 
    277       ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj)    , hv_0(jpi,jpj)     ,                                             & 
    278          &      ht  (jpi,jpj) , hu(  jpi,jpj,jpt), hv(  jpi,jpj,jpt) , r1_hu(jpi,jpj,jpt) , r1_hv(jpi,jpj,jpt) ,   & 
    279          &                      STAT=ierr(6)  ) 
    280          ! 
    281       ALLOCATE( risfdep(jpi,jpj) , bathy(jpi,jpj) , STAT=ierr(7)  )  
    282          ! 
    283       ALLOCATE( gdept_1d(jpk) , gdepw_1d(jpk) , e3t_1d(jpk) , e3w_1d(jpk) , STAT=ierr(8) ) 
    284          ! 
     285         &      gdept  (jpi,jpj,jpk,jpt) , gdepw  (jpi,jpj,jpk,jpt) , gde3w  (jpi,jpj,jpk) , STAT=ierr(ii) ) 
     286         ! 
     287      ii = ii+1 
     288      ALLOCATE(  e3t_0(jpi,jpj,jpk) , e3u_0 (jpi,jpj,jpk) , e3v_0 (jpi,jpj,jpk) , e3f_0(jpi,jpj,jpk) ,      & 
     289         &       e3w_0(jpi,jpj,jpk) , e3uw_0(jpi,jpj,jpk) , e3vw_0(jpi,jpj,jpk)                      ,  STAT=ierr(ii) ) 
     290         ! 
     291#if ! defined key_qco 
     292      ii = ii+1 
     293      ALLOCATE( e3t(jpi,jpj,jpk,jpt) , e3u (jpi,jpj,jpk,jpt) , e3v (jpi,jpj,jpk,jpt) , e3f(jpi,jpj,jpk) ,      & 
     294         &      e3w(jpi,jpj,jpk,jpt) , e3uw(jpi,jpj,jpk,jpt) , e3vw(jpi,jpj,jpk,jpt)                    ,  STAT=ierr(ii) ) 
     295#endif   
     296         ! 
     297      ii = ii+1 
     298      ALLOCATE( r3t  (jpi,jpj,jpt)   , r3u  (jpi,jpj,jpt)    , r3v  (jpi,jpj,jpt)    , r3f  (jpi,jpj) ,  & 
     299         &      r3t_f(jpi,jpj)       , r3u_f(jpi,jpj)        , r3v_f(jpi,jpj)                         ,  STAT=ierr(ii) )        
     300         ! 
     301      ii = ii+1 
     302      ALLOCATE( ht_0(jpi,jpj) ,    hu_0(jpi,jpj)    ,    hv_0(jpi,jpj)     , hf_0(jpi,jpj) ,       & 
     303         &   r1_ht_0(jpi,jpj) , r1_hu_0(jpi,jpj) ,    r1_hv_0(jpi,jpj),   r1_hf_0(jpi,jpj) ,   STAT=ierr(ii)  ) 
     304         ! 
     305#if ! defined key_qco 
     306      ii = ii+1 
     307      ALLOCATE( ht  (jpi,jpj) ,    hu  (jpi,jpj,jpt),    hv  (jpi,jpj,jpt)                 ,       & 
     308         &                      r1_hu  (jpi,jpj,jpt), r1_hv  (jpi,jpj,jpt)                 ,   STAT=ierr(ii)  ) 
     309#else 
     310      ii = ii+1 
     311      ALLOCATE(                    hu  (jpi,jpj,jpt),    hv  (jpi,jpj,jpt)                 ,       & 
     312         &                      r1_hu  (jpi,jpj,jpt), r1_hv  (jpi,jpj,jpt)                 ,   STAT=ierr(ii)  ) 
     313#endif 
     314         ! 
     315      ii = ii+1 
     316      ALLOCATE( risfdep(jpi,jpj) , bathy(jpi,jpj) , STAT=ierr(ii)  )  
     317         ! 
     318      ii = ii+1 
     319      ALLOCATE( gdept_1d(jpk) , gdepw_1d(jpk) , e3t_1d(jpk) , e3w_1d(jpk) , STAT=ierr(ii) ) 
     320         ! 
     321      ii = ii+1 
    285322      ALLOCATE( tmask_i(jpi,jpj) , tmask_h(jpi,jpj) ,                        &  
    286          &      ssmask (jpi,jpj) , ssumask(jpi,jpj) , ssvmask(jpi,jpj) ,     & 
    287          &      mbkt   (jpi,jpj) , mbku   (jpi,jpj) , mbkv   (jpi,jpj) , STAT=ierr(9) ) 
    288          ! 
    289       ALLOCATE( mikt(jpi,jpj), miku(jpi,jpj), mikv(jpi,jpj), mikf(jpi,jpj), STAT=ierr(10) ) 
    290          ! 
     323         &      ssmask (jpi,jpj) , ssumask(jpi,jpj) , ssvmask(jpi,jpj) , ssfmask(jpi,jpj) ,     & 
     324         &      mbkt   (jpi,jpj) , mbku   (jpi,jpj) , mbkv   (jpi,jpj) ,                    STAT=ierr(ii) ) 
     325         ! 
     326      ii = ii+1 
     327      ALLOCATE( mikt(jpi,jpj), miku(jpi,jpj), mikv(jpi,jpj), mikf(jpi,jpj), STAT=ierr(ii) ) 
     328         ! 
     329      ii = ii+1 
    291330      ALLOCATE( tmask(jpi,jpj,jpk) , umask(jpi,jpj,jpk) ,     &  
    292          &      vmask(jpi,jpj,jpk) , fmask(jpi,jpj,jpk) , STAT=ierr(11) ) 
    293          ! 
    294       ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(12) ) 
     331         &      vmask(jpi,jpj,jpk) , fmask(jpi,jpj,jpk) , STAT=ierr(ii) ) 
     332         ! 
     333      ii = ii+1 
     334      ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(ii) ) 
    295335      ! 
    296336      dom_oce_alloc = MAXVAL(ierr) 
  • NEMO/trunk/src/OCE/DOM/domain.F90

    r13216 r13237  
    1515   !!            3.7  !  2015-11  (G. Madec, A. Coward)  time varying zgr by default 
    1616   !!            4.0  !  2016-10  (G. Madec, S. Flavoni)  domain configuration / user defined interface 
     17   !!            4.x  ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio 
    1718   !!---------------------------------------------------------------------- 
    1819    
     
    3435   USE dommsk         ! domain: set the mask system 
    3536   USE domwri         ! domain: write the meshmask file 
     37#if ! defined key_qco 
    3638   USE domvvl         ! variable volume 
     39#else 
     40   USE domqco          ! variable volume 
     41#endif 
    3742   USE c1d            ! 1D configuration 
    3843   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine) 
     
    7681      CHARACTER (len=*), INTENT(in) :: cdstr                  ! model: NEMO or SAS. Determines core restart variables 
    7782      ! 
    78       INTEGER ::   ji, jj, jk, ik   ! dummy loop indices 
     83      INTEGER ::   ji, jj, jk, jt   ! dummy loop indices 
    7984      INTEGER ::   iconf = 0    ! local integers 
    8085      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))"  
     
    110115         CASE( 7 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north-south)' 
    111116         CASE DEFAULT 
    112             CALL ctl_stop( 'jperio is out of range' ) 
     117            CALL ctl_stop( 'dom_init:   jperio is out of range' ) 
    113118         END SELECT 
    114119         WRITE(numout,*)     '      Ocean model configuration used:' 
     
    140145      IF( ln_closea ) CALL dom_clo      ! Read in masks to define closed seas and lakes 
    141146 
    142       CALL dom_zgr( ik_top, ik_bot )    ! Vertical mesh and bathymetry 
     147      CALL dom_zgr( ik_top, ik_bot )    ! Vertical mesh and bathymetry (return top and bottom ocean t-level indices) 
    143148 
    144149      CALL dom_msk( ik_top, ik_bot )    ! Masks 
     
    147152      hu_0(:,:) = 0._wp 
    148153      hv_0(:,:) = 0._wp 
     154      hf_0(:,:) = 0._wp 
    149155      DO jk = 1, jpk 
    150156         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk) 
    151157         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk) 
    152158         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk) 
     159         hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,jk) * fmask(:,:,jk) 
    153160      END DO 
    154161      ! 
     162      r1_ht_0(:,:) = ssmask (:,:) / ( ht_0(:,:) + 1._wp -  ssmask (:,:) ) 
     163      r1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp -  ssumask(:,:) ) 
     164      r1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp -  ssvmask(:,:) ) 
     165      r1_hf_0(:,:) = ssfmask(:,:) / ( hf_0(:,:) + 1._wp -  ssfmask(:,:) ) 
     166 
     167      ! 
     168#if defined key_qco 
     169      !           !==  initialisation of time varying coordinate  ==!   Quasi-Euerian coordinate case 
     170      ! 
     171      IF( .NOT.l_offline )   CALL dom_qco_init( Kbb, Kmm, Kaa ) 
     172      ! 
     173      IF( ln_linssh )        CALL ctl_stop('STOP','domain: key_qco and ln_linssh = T are incompatible') 
     174      ! 
     175#else 
    155176      !           !==  time varying part of coordinate system  ==! 
    156177      ! 
    157178      IF( ln_linssh ) THEN       != Fix in time : set to the reference one for all 
    158179      ! 
    159          !       before        !          now          !       after         ! 
    160             gdept(:,:,:,Kbb) = gdept_0  ;   gdept(:,:,:,Kmm) = gdept_0   ;   gdept(:,:,:,Kaa) = gdept_0   ! depth of grid-points 
    161             gdepw(:,:,:,Kbb) = gdepw_0  ;   gdepw(:,:,:,Kmm) = gdepw_0   ;   gdepw(:,:,:,Kaa) = gdepw_0   ! 
    162                                    gde3w = gde3w_0   !        ---          ! 
    163          !                                                                   
    164               e3t(:,:,:,Kbb) =   e3t_0  ;     e3t(:,:,:,Kmm) =   e3t_0   ;   e3t(:,:,:,Kaa) =  e3t_0    ! scale factors 
    165               e3u(:,:,:,Kbb) =   e3u_0  ;     e3u(:,:,:,Kmm) =   e3u_0   ;   e3u(:,:,:,Kaa) =  e3u_0    ! 
    166               e3v(:,:,:,Kbb) =   e3v_0  ;     e3v(:,:,:,Kmm) =   e3v_0   ;   e3v(:,:,:,Kaa) =  e3v_0    ! 
    167                                      e3f =   e3f_0   !        ---          ! 
    168               e3w(:,:,:,Kbb) =   e3w_0  ;     e3w(:,:,:,Kmm) =   e3w_0   ;    e3w(:,:,:,Kaa) =   e3w_0   !  
    169              e3uw(:,:,:,Kbb) =  e3uw_0  ;    e3uw(:,:,:,Kmm) =  e3uw_0   ;   e3uw(:,:,:,Kaa) =  e3uw_0   !   
    170              e3vw(:,:,:,Kbb) =  e3vw_0  ;    e3vw(:,:,:,Kmm) =  e3vw_0   ;   e3vw(:,:,:,Kaa) =  e3vw_0   ! 
    171          ! 
    172          z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF 
    173          z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) ) 
    174          ! 
    175          !        before       !          now          !       after         ! 
    176                                       ht =    ht_0   !                     ! water column thickness 
    177                hu(:,:,Kbb) =    hu_0  ;      hu(:,:,Kmm) =    hu_0   ;    hu(:,:,Kaa) =    hu_0   !  
    178                hv(:,:,Kbb) =    hv_0  ;      hv(:,:,Kmm) =    hv_0   ;    hv(:,:,Kaa) =    hv_0   ! 
    179             r1_hu(:,:,Kbb) = z1_hu_0  ;   r1_hu(:,:,Kmm) = z1_hu_0   ; r1_hu(:,:,Kaa) = z1_hu_0   ! inverse of water column thickness 
    180             r1_hv(:,:,Kbb) = z1_hv_0  ;   r1_hv(:,:,Kmm) = z1_hv_0   ; r1_hv(:,:,Kaa) = z1_hv_0   ! 
    181          ! 
     180         DO jt = 1, jpt                         ! depth of t- and w-grid-points 
     181            gdept(:,:,:,jt) = gdept_0(:,:,:) 
     182            gdepw(:,:,:,jt) = gdepw_0(:,:,:) 
     183         END DO 
     184            gde3w(:,:,:)    = gde3w_0(:,:,:)    ! = gdept as the sum of e3t 
     185         ! 
     186         DO jt = 1, jpt                         ! vertical scale factors 
     187            e3t(:,:,:,jt) =  e3t_0(:,:,:) 
     188            e3u(:,:,:,jt) =  e3u_0(:,:,:) 
     189            e3v(:,:,:,jt) =  e3v_0(:,:,:) 
     190            e3w(:,:,:,jt) =  e3w_0(:,:,:) 
     191            e3uw(:,:,:,jt) = e3uw_0(:,:,:) 
     192            e3vw(:,:,:,jt) = e3vw_0(:,:,:) 
     193         END DO 
     194            e3f(:,:,:)    =  e3f_0(:,:,:) 
     195         ! 
     196         DO jt = 1, jpt                         ! water column thickness and its inverse 
     197            hu(:,:,jt)    =    hu_0(:,:) 
     198            hv(:,:,jt)    =    hv_0(:,:) 
     199            r1_hu(:,:,jt) = r1_hu_0(:,:) 
     200            r1_hv(:,:,jt) = r1_hv_0(:,:) 
     201         END DO 
     202            ht(:,:) =    ht_0(:,:) 
    182203         ! 
    183204      ELSE                       != time varying : initialize before/now/after variables 
     
    186207         ! 
    187208      ENDIF 
     209#endif 
     210 
    188211      ! 
    189212 
  • NEMO/trunk/src/OCE/DOM/dommsk.F90

    r13226 r13237  
    1818   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask 
    1919   !!            4.0  ! 2016-06  (G. Madec, S. Flavoni)  domain configuration / user defined interface 
     20   !!            4.x  ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio 
    2021   !!---------------------------------------------------------------------- 
    2122 
     
    192193      ssumask(:,:) = MAXVAL( umask(:,:,:), DIM=3 ) 
    193194      ssvmask(:,:) = MAXVAL( vmask(:,:,:), DIM=3 ) 
     195      ssfmask(:,:) = MAXVAL( fmask(:,:,:), DIM=3 ) 
    194196 
    195197 
  • NEMO/trunk/src/OCE/DOM/domvvl.F90

    r12740 r13237  
    99   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
    1010   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rename dom_vvl_sf_swp -> dom_vvl_sf_update for new timestepping 
     11   !!            4.x  ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio 
    1112   !!---------------------------------------------------------------------- 
    1213 
    13    !!---------------------------------------------------------------------- 
    14    !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
    15    !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
    16    !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
    17    !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
    18    !!   dom_vvl_rst      : read/write restart file 
    19    !!   dom_vvl_ctl      : Check the vvl options 
    20    !!---------------------------------------------------------------------- 
    2114   USE oce             ! ocean dynamics and tracers 
    2215   USE phycst          ! physical constant 
     
    3629   PRIVATE 
    3730 
    38    PUBLIC  dom_vvl_init       ! called by domain.F90 
    39    PUBLIC  dom_vvl_zgr        ! called by isfcpl.F90 
    40    PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
    41    PUBLIC  dom_vvl_sf_update  ! called by step.F90 
    42    PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
    43  
    4431   !                                                      !!* Namelist nam_vvl 
    4532   LOGICAL , PUBLIC :: ln_vvl_zstar           = .FALSE.    ! zstar  vertical coordinate 
     
    6350   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_hdv                 ! retoring period for low freq. divergence 
    6451 
     52#if defined key_qco 
     53   !!---------------------------------------------------------------------- 
     54   !!   'key_qco'      EMPTY MODULE      Quasi-Eulerian vertical coordonate 
     55   !!---------------------------------------------------------------------- 
     56#else 
     57   !!---------------------------------------------------------------------- 
     58   !!   Default key      Old management of time varying vertical coordinate 
     59   !!---------------------------------------------------------------------- 
     60    
     61   !!---------------------------------------------------------------------- 
     62   !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
     63   !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
     64   !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
     65   !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
     66   !!   dom_vvl_rst      : read/write restart file 
     67   !!   dom_vvl_ctl      : Check the vvl options 
     68   !!---------------------------------------------------------------------- 
     69 
     70   PUBLIC  dom_vvl_init       ! called by domain.F90 
     71   PUBLIC  dom_vvl_zgr        ! called by isfcpl.F90 
     72   PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
     73   PUBLIC  dom_vvl_sf_update  ! called by step.F90 
     74   PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
     75    
    6576   !! * Substitutions 
    6677#  include "do_loop_substitute.h90" 
     
    135146      ! 
    136147   END SUBROUTINE dom_vvl_init 
    137    ! 
     148 
     149 
    138150   SUBROUTINE dom_vvl_zgr(Kbb, Kmm, Kaa) 
    139151      !!---------------------------------------------------------------------- 
     
    10291041   END SUBROUTINE dom_vvl_ctl 
    10301042 
     1043#endif 
     1044 
    10311045   !!====================================================================== 
    10321046END MODULE domvvl 
  • NEMO/trunk/src/OCE/DOM/istate.F90

    r13216 r13237  
    4646   !! * Substitutions 
    4747#  include "do_loop_substitute.h90" 
     48#  include "domzgr_substitute.h90" 
    4849   !!---------------------------------------------------------------------- 
    4950   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6263      ! 
    6364      INTEGER ::   ji, jj, jk   ! dummy loop indices 
     65      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zgdept     ! 3D table  !!st patch to use gdept subtitute 
    6466!!gm see comment further down 
    6567      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   zuvd    ! U & V data workspace 
     
    132134             ! 
    133135         ELSE                                 ! user defined initial T and S 
    134             CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )          
     136            DO jk = 1, jpk 
     137               zgdept(:,:,jk) = gdept(:,:,jk,Kbb) 
     138            END DO 
     139            CALL usr_def_istate( zgdept, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )          
    135140         ENDIF 
    136141         ts  (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb)       ! set now values from to before ones 
     
    141146!!gm POTENTIAL BUG : 
    142147!!gm  ISSUE :  if ssh(:,:,Kbb) /= 0  then, in non linear free surface, the e3._n, e3._b should be recomputed 
    143 !!             as well as gdept and gdepw....   !!!!!  
     148!!             as well as gdept_ and gdepw_....   !!!!!  
    144149!!      ===>>>>   probably a call to domvvl initialisation here.... 
    145150 
  • NEMO/trunk/src/OCE/DYN/divhor.F90

    r13226 r13237  
    4040   !! * Substitutions 
    4141#  include "do_loop_substitute.h90" 
     42#  include "domzgr_substitute.h90" 
    4243   !!---------------------------------------------------------------------- 
    4344   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/DYN/dynadv_cen2.F90

    r12377 r13237  
    2828   !! * Substitutions 
    2929#  include "do_loop_substitute.h90" 
     30#  include "domzgr_substitute.h90" 
    3031   !!---------------------------------------------------------------------- 
    3132   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    7980         DO_2D_00_00 
    8081            puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - (  zfu_t(ji+1,jj,jk) - zfu_t(ji,jj  ,jk)    & 
    81                &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
     82               &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj)   & 
     83               &                           / e3u(ji,jj,jk,Kmm) 
    8284            pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - (  zfu_f(ji,jj  ,jk) - zfu_f(ji-1,jj,jk)    & 
    83                &                           + zfv_t(ji,jj+1,jk) - zfv_t(ji  ,jj,jk)  ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     85               &                           + zfv_t(ji,jj+1,jk) - zfv_t(ji  ,jj,jk)  ) * r1_e1e2v(ji,jj)   & 
     86               &                           / e3v(ji,jj,jk,Kmm) 
    8487         END_2D 
    8588      END DO 
     
    115118      END DO 
    116119      DO_3D_00_00( 1, jpkm1 ) 
    117          puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
    118          pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     120         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
     121            &                                      / e3u(ji,jj,jk,Kmm) 
     122         pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj)   & 
     123            &                                      / e3v(ji,jj,jk,Kmm) 
    119124      END_3D 
    120125      ! 
  • NEMO/trunk/src/OCE/DYN/dynadv_ubs.F90

    r13226 r13237  
    3434   !! * Substitutions 
    3535#  include "do_loop_substitute.h90" 
     36#  include "domzgr_substitute.h90" 
    3637   !!---------------------------------------------------------------------- 
    3738   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    169170         DO_2D_00_00 
    170171            puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - (  zfu_t(ji+1,jj,jk) - zfu_t(ji,jj  ,jk)    & 
    171                &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
     172               &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj)   & 
     173               &                           / e3u(ji,jj,jk,Kmm) 
    172174            pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - (  zfu_f(ji,jj  ,jk) - zfu_f(ji-1,jj,jk)    & 
    173                &                           + zfv_t(ji,jj+1,jk) - zfv_t(ji  ,jj,jk)  ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     175               &                           + zfv_t(ji,jj+1,jk) - zfv_t(ji  ,jj,jk)  ) * r1_e1e2v(ji,jj)   & 
     176               &                           / e3v(ji,jj,jk,Kmm) 
    174177         END_2D 
    175178      END DO 
     
    206209      END DO 
    207210      DO_3D_00_00( 1, jpkm1 ) 
    208          puu(ji,jj,jk,Krhs) =  puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
    209          pvv(ji,jj,jk,Krhs) =  pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     211         puu(ji,jj,jk,Krhs) =  puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
     212            &                                       / e3u(ji,jj,jk,Kmm) 
     213         pvv(ji,jj,jk,Krhs) =  pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj)   & 
     214            &                                       / e3v(ji,jj,jk,Kmm) 
    210215      END_3D 
    211216      ! 
  • NEMO/trunk/src/OCE/DYN/dynatf.F90

    r13226 r13237  
    5959   PUBLIC    dyn_atf   ! routine called by step.F90 
    6060 
     61#if defined key_qco 
     62   !!---------------------------------------------------------------------- 
     63   !!   'key_qco'      EMPTY ROUTINE     Quasi-Eulerian vertical coordonate 
     64   !!---------------------------------------------------------------------- 
     65CONTAINS 
     66 
     67   SUBROUTINE dyn_atf ( kt, Kbb, Kmm, Kaa, puu, pvv, pe3t, pe3u, pe3v ) 
     68      INTEGER                             , INTENT(in   ) :: kt               ! ocean time-step index 
     69      INTEGER                             , INTENT(in   ) :: Kbb, Kmm, Kaa    ! before and after time level indices 
     70      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv         ! velocities to be time filtered 
     71      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: pe3t, pe3u, pe3v ! scale factors to be time filtered 
     72 
     73      WRITE(*,*) 'dyn_atf: You should not have seen this print! error?', kt 
     74   END SUBROUTINE dyn_atf 
     75 
     76#else 
     77 
    6178   !! * Substitutions 
    6279#  include "do_loop_substitute.h90" 
     
    198215            zwfld(:,:) = emp_b(:,:) - emp(:,:) 
    199216            IF ( ln_rnf ) zwfld(:,:) =  zwfld(:,:) - ( rnf_b(:,:) - rnf(:,:) ) 
     217 
    200218            DO jk = 1, jpkm1 
    201219               ze3t_f(:,:,jk) = ze3t_f(:,:,jk) - zcoef * zwfld(:,:) * tmask(:,:,jk) & 
     
    312330   END SUBROUTINE dyn_atf 
    313331 
     332#endif 
     333 
    314334   !!========================================================================= 
    315335END MODULE dynatf 
  • NEMO/trunk/src/OCE/DYN/dynhpg.F90

    r13226 r13237  
    7676   !! * Substitutions 
    7777#  include "do_loop_substitute.h90" 
     78#  include "domzgr_substitute.h90" 
     79 
    7880   !!---------------------------------------------------------------------- 
    7981   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    452454      DO_2D_00_00 
    453455         ! hydrostatic pressure gradient along s-surfaces 
    454          zhpi(ji,jj,1) = zcoef0 * (  e3w(ji+1,jj  ,1,Kmm) * ( znad + rhd(ji+1,jj  ,1) )    & 
    455             &                      - e3w(ji  ,jj  ,1,Kmm) * ( znad + rhd(ji  ,jj  ,1) )  ) * r1_e1u(ji,jj) 
    456          zhpj(ji,jj,1) = zcoef0 * (  e3w(ji  ,jj+1,1,Kmm) * ( znad + rhd(ji  ,jj+1,1) )    & 
    457             &                      - e3w(ji  ,jj  ,1,Kmm) * ( znad + rhd(ji  ,jj  ,1) )  ) * r1_e2v(ji,jj) 
     456         zhpi(ji,jj,1) =   & 
     457            &  zcoef0 * (  e3w(ji+1,jj  ,1,Kmm) * ( znad + rhd(ji+1,jj  ,1) )    & 
     458            &            - e3w(ji  ,jj  ,1,Kmm) * ( znad + rhd(ji  ,jj  ,1) )  ) & 
     459            &           * r1_e1u(ji,jj) 
     460         zhpj(ji,jj,1) =   & 
     461            &  zcoef0 * (  e3w(ji  ,jj+1,1,Kmm) * ( znad + rhd(ji  ,jj+1,1) )    & 
     462            &            - e3w(ji  ,jj  ,1,Kmm) * ( znad + rhd(ji  ,jj  ,1) )  ) & 
     463            &           * r1_e2v(ji,jj) 
    458464         ! s-coordinate pressure gradient correction 
    459465         zuap = -zcoef0 * ( rhd    (ji+1,jj,1) + rhd    (ji,jj,1) + 2._wp * znad )   & 
     
    589595         ! hydrostatic pressure gradient along s-surfaces 
    590596         zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 / e1u(ji,jj)   & 
    591             &           * (  e3w(ji+1,jj,jk,Kmm) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) * wmask(ji+1,jj,jk)   & 
    592             &              - e3w(ji  ,jj,jk,Kmm) * ( rhd(ji  ,jj,jk) + rhd(ji  ,jj,jk-1) + 2*znad ) * wmask(ji  ,jj,jk)   ) 
     597            &           * (  e3w(ji+1,jj,jk,Kmm)                   & 
     598            &                  * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) * wmask(ji+1,jj,jk)   & 
     599            &              - e3w(ji  ,jj,jk,Kmm)                   & 
     600            &                  * ( rhd(ji  ,jj,jk) + rhd(ji  ,jj,jk-1) + 2*znad ) * wmask(ji  ,jj,jk)   ) 
    593601         zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 / e2v(ji,jj)   & 
    594             &           * (  e3w(ji,jj+1,jk,Kmm) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) * wmask(ji,jj+1,jk)   & 
    595             &              - e3w(ji,jj  ,jk,Kmm) * ( rhd(ji,jj,  jk) + rhd(ji,jj  ,jk-1) + 2*znad ) * wmask(ji,jj  ,jk)   ) 
     602            &           * (  e3w(ji,jj+1,jk,Kmm)                   & 
     603            &                  * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) * wmask(ji,jj+1,jk)   & 
     604            &              - e3w(ji,jj  ,jk,Kmm)                   & 
     605            &                  * ( rhd(ji,jj,  jk) + rhd(ji,jj  ,jk-1) + 2*znad ) * wmask(ji,jj  ,jk)   ) 
    596606         ! s-coordinate pressure gradient correction 
    597607         zuap = -zcoef0 * ( rhd   (ji+1,jj  ,jk) + rhd   (ji,jj,jk) + 2._wp * znad )   & 
     
    771781      !------------------------------------------------------------- 
    772782 
    773 !!bug gm   :  e3w-gde3w = 0.5*e3w  ....  and gde3w(2)-gde3w(1)=e3w(2) ....   to be verified 
    774 !          true if gde3w is really defined as the sum of the e3w scale factors as, it seems to me, it should be 
     783!!bug gm   :  e3w-gde3w(:,:,:) = 0.5*e3w  ....  and gde3w(:,:,2)-gde3w(:,:,1)=e3w(:,:,2,Kmm) ....   to be verified 
     784!          true if gde3w(:,:,:) is really defined as the sum of the e3w scale factors as, it seems to me, it should be 
    775785 
    776786      DO_2D_00_00 
     
    13591369   !!====================================================================== 
    13601370END MODULE dynhpg 
    1361  
  • NEMO/trunk/src/OCE/DYN/dynldf_iso.F90

    r13226 r13237  
    4242   !! * Substitutions 
    4343#  include "do_loop_substitute.h90" 
     44#  include "domzgr_substitute.h90" 
    4445   !!---------------------------------------------------------------------- 
    4546   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    168169         IF( ln_zps ) THEN      ! z-coordinate - partial steps : min(e3u) 
    169170            DO_2D_00_01 
    170                zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * MIN( e3u(ji,jj,jk,Kmm), e3u(ji-1,jj,jk,Kmm) ) * r1_e1t(ji,jj) 
     171               zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj)   & 
     172                  &    * MIN( e3u(ji  ,jj,jk,Kmm),                & 
     173                  &           e3u(ji-1,jj,jk,Kmm) ) * r1_e1t(ji,jj) 
    171174 
    172175               zmskt = 1._wp / MAX(   umask(ji-1,jj,jk  )+umask(ji,jj,jk+1)     & 
     
    181184         ELSE                   ! other coordinate system (zco or sco) : e3t 
    182185            DO_2D_00_01 
    183                zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e1t(ji,jj) 
     186               zabe1 = ( ahmt(ji,jj,jk)+rn_ahm_b )   & 
     187                  &     * e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e1t(ji,jj) 
    184188 
    185189               zmskt = 1._wp / MAX(   umask(ji-1,jj,jk  ) + umask(ji,jj,jk+1)     & 
     
    196200         ! j-flux at f-point 
    197201         DO_2D_10_10 
    198             zabe2 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e1f(ji,jj) * e3f(ji,jj,jk) * r1_e2f(ji,jj) 
     202            zabe2 = ( ahmf(ji,jj,jk) + rn_ahm_b )   & 
     203               &     * e1f(ji,jj) * e3f(ji,jj,jk) * r1_e2f(ji,jj) 
    199204 
    200205            zmskf = 1._wp / MAX(   umask(ji,jj+1,jk  )+umask(ji,jj,jk+1)     & 
     
    215220 
    216221         DO_2D_00_10 
    217             zabe1 = ( ahmf(ji,jj,jk) + rn_ahm_b ) * e2f(ji,jj) * e3f(ji,jj,jk) * r1_e1f(ji,jj) 
     222            zabe1 = ( ahmf(ji,jj,jk) + rn_ahm_b )   & 
     223               &     * e2f(ji,jj) * e3f(ji,jj,jk) * r1_e1f(ji,jj) 
    218224 
    219225            zmskf = 1._wp / MAX(  vmask(ji+1,jj,jk  )+vmask(ji,jj,jk+1)     & 
     
    230236         IF( ln_zps ) THEN      ! z-coordinate - partial steps : min(e3u) 
    231237            DO_2D_01_10 
    232                zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * MIN( e3v(ji,jj,jk,Kmm), e3v(ji,jj-1,jk,Kmm) ) * r1_e2t(ji,jj) 
     238               zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj)   & 
     239                  &     * MIN( e3v(ji,jj  ,jk,Kmm),                 & 
     240                  &            e3v(ji,jj-1,jk,Kmm) ) * r1_e2t(ji,jj) 
    233241 
    234242               zmskt = 1._wp / MAX(  vmask(ji,jj-1,jk  )+vmask(ji,jj,jk+1)     & 
     
    243251         ELSE                   ! other coordinate system (zco or sco) : e3t 
    244252            DO_2D_01_10 
    245                zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b ) * e1t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e2t(ji,jj) 
     253               zabe2 = ( ahmt(ji,jj,jk)+rn_ahm_b )   & 
     254                  &     * e1t(ji,jj) * e3t(ji,jj,jk,Kmm) * r1_e2t(ji,jj) 
    246255 
    247256               zmskt = 1./MAX(  vmask(ji,jj-1,jk  )+vmask(ji,jj,jk+1)   & 
     
    261270         DO_2D_00_00 
    262271            puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (  ziut(ji+1,jj) - ziut(ji,jj  )    & 
    263                &                           + zjuf(ji  ,jj) - zjuf(ji,jj-1)  ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
     272               &                           + zjuf(ji  ,jj) - zjuf(ji,jj-1)  ) * r1_e1e2u(ji,jj)   & 
     273               &                           / e3u(ji,jj,jk,Kmm) 
    264274            pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (  zivf(ji,jj  ) - zivf(ji-1,jj)    & 
    265                &                           + zjvt(ji,jj+1) - zjvt(ji,jj  )  ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     275               &                           + zjvt(ji,jj+1) - zjvt(ji,jj  )  ) * r1_e1e2v(ji,jj)   & 
     276               &                           / e3v(ji,jj,jk,Kmm) 
    266277         END_2D 
    267278         !                                             ! =============== 
     
    375386         DO jk = 1, jpkm1 
    376387            DO ji = 2, jpim1 
    377                puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + ( zfuw(ji,jk) - zfuw(ji,jk+1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
    378                pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + ( zfvw(ji,jk) - zfvw(ji,jk+1) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     388               puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + ( zfuw(ji,jk) - zfuw(ji,jk+1) ) * r1_e1e2u(ji,jj)   & 
     389                  &               / e3u(ji,jj,jk,Kmm) 
     390               pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + ( zfvw(ji,jk) - zfvw(ji,jk+1) ) * r1_e1e2v(ji,jj)   & 
     391                  &               / e3v(ji,jj,jk,Kmm) 
    379392            END DO 
    380393         END DO 
  • NEMO/trunk/src/OCE/DYN/dynldf_lap_blp.F90

    r13226 r13237  
    2828   !! * Substitutions 
    2929#  include "do_loop_substitute.h90" 
     30#  include "domzgr_substitute.h90" 
    3031   !!---------------------------------------------------------------------- 
    3132   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/DYN/dynspg_ts.F90

    r13216 r13237  
    8787   !! * Substitutions 
    8888#  include "do_loop_substitute.h90" 
     89#  include "domzgr_substitute.h90" 
    8990   !!---------------------------------------------------------------------- 
    9091   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    161162      REAL(wp), DIMENSION(jpi,jpj) :: zCdU_u, zCdU_v   ! top/bottom stress at u- & v-points 
    162163      REAL(wp), DIMENSION(jpi,jpj) :: zhU, zhV         ! fluxes 
     164      REAL(wp), DIMENSION(jpi, jpj, jpk) :: ze3u, ze3v 
    163165      ! 
    164166      REAL(wp) ::   zwdramp                     ! local scalar - only used if ln_wd_dl = .True.  
     
    227229      !                                   !=  zu_frc =  1/H e3*d/dt(Ua)  =!  (Vertical mean of Ua, the 3D trends) 
    228230      !                                   !  ---------------------------  ! 
    229       zu_frc(:,:) = SUM( e3u(:,:,:,Kmm) * uu(:,:,:,Krhs) * umask(:,:,:) , DIM=3 ) * r1_hu(:,:,Kmm) 
    230       zv_frc(:,:) = SUM( e3v(:,:,:,Kmm) * vv(:,:,:,Krhs) * vmask(:,:,:) , DIM=3 ) * r1_hv(:,:,Kmm) 
     231      DO jk = 1 , jpk 
     232         ze3u(:,:,jk) = e3u(:,:,jk,Kmm) 
     233         ze3v(:,:,jk) = e3v(:,:,jk,Kmm) 
     234      END DO 
     235      ! 
     236      zu_frc(:,:) = SUM( ze3u(:,:,:) * uu(:,:,:,Krhs) * umask(:,:,:) , DIM=3 ) * r1_hu(:,:,Kmm) 
     237      zv_frc(:,:) = SUM( ze3v(:,:,:) * vv(:,:,:,Krhs) * vmask(:,:,:) , DIM=3 ) * r1_hv(:,:,Kmm) 
    231238      ! 
    232239      ! 
     
    250257      zhV(:,:) = pvv_b(:,:,Kmm) * hv(:,:,Kmm) * e1v(:,:)        ! NB: FULL domain : put a value in last row and column 
    251258      ! 
    252       CALL dyn_cor_2d( hu(:,:,Kmm), hv(:,:,Kmm), puu_b(:,:,Kmm), pvv_b(:,:,Kmm), zhU, zhV,  &   ! <<== in 
     259      CALL dyn_cor_2d( ht(:,:), hu(:,:,Kmm), hv(:,:,Kmm), puu_b(:,:,Kmm), pvv_b(:,:,Kmm), zhU, zhV,  &   ! <<== in 
    253260         &                                                                     zu_trd, zv_trd   )   ! ==>> out 
    254261      ! 
     
    568575         ! at each time step. We however keep them constant here for optimization. 
    569576         ! Recall that zhU and zhV hold fluxes at jn+0.5 (extrapolated not backward interpolated) 
    570          CALL dyn_cor_2d( zhup2_e, zhvp2_e, ua_e, va_e, zhU, zhV,    zu_trd, zv_trd   ) 
     577         CALL dyn_cor_2d( zhtp2_e, zhup2_e, zhvp2_e, ua_e, va_e, zhU, zhV,    zu_trd, zv_trd   ) 
    571578         ! 
    572579         ! Add tidal astronomical forcing if defined 
     
    10911098      ! 
    10921099      SELECT CASE( nvor_scheme ) 
    1093       CASE( np_EEN )                != EEN scheme using e3f (energy & enstrophy scheme) 
     1100      CASE( np_EEN )                != EEN scheme using e3f energy & enstrophy scheme 
    10941101         SELECT CASE( nn_een_e3f )              !* ff_f/e3 at F-point 
    10951102         CASE ( 0 )                                   ! original formulation  (masked averaging of e3t divided by 4) 
     
    11181125         END_2D 
    11191126         ! 
    1120       CASE( np_EET )                  != EEN scheme using e3t (energy conserving scheme) 
     1127      CASE( np_EET )                  != EEN scheme using e3t energy conserving scheme 
    11211128         ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 
    11221129         DO_2D_01_01 
     
    11821189 
    11831190 
    1184    SUBROUTINE dyn_cor_2d( phu, phv, punb, pvnb, zhU, zhV,    zu_trd, zv_trd   ) 
     1191   SUBROUTINE dyn_cor_2d( pht, phu, phv, punb, pvnb, zhU, zhV,    zu_trd, zv_trd   ) 
    11851192      !!--------------------------------------------------------------------- 
    11861193      !!                   ***  ROUTINE dyn_cor_2d  *** 
     
    11901197      INTEGER  ::   ji ,jj                             ! dummy loop indices 
    11911198      REAL(wp) ::   zx1, zx2, zy1, zy2, z1_hu, z1_hv   !   -      - 
    1192       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: phu, phv, punb, pvnb, zhU, zhV 
     1199      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pht, phu, phv, punb, pvnb, zhU, zhV 
    11931200      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: zu_trd, zv_trd 
    11941201      !!---------------------------------------------------------------------- 
     
    11991206            z1_hv = ssvmask(ji,jj) / ( phv(ji,jj) + 1._wp - ssvmask(ji,jj) ) 
    12001207            zu_trd(ji,jj) = + r1_4 * r1_e1e2u(ji,jj) * z1_hu                    & 
    1201                &               * (  e1e2t(ji+1,jj)*ht(ji+1,jj)*ff_t(ji+1,jj) * ( pvnb(ji+1,jj) + pvnb(ji+1,jj-1) )   & 
    1202                &                  + e1e2t(ji  ,jj)*ht(ji  ,jj)*ff_t(ji  ,jj) * ( pvnb(ji  ,jj) + pvnb(ji  ,jj-1) )   ) 
     1208               &               * (  e1e2t(ji+1,jj)*pht(ji+1,jj)*ff_t(ji+1,jj) * ( pvnb(ji+1,jj) + pvnb(ji+1,jj-1) )   & 
     1209               &                  + e1e2t(ji  ,jj)*pht(ji  ,jj)*ff_t(ji  ,jj) * ( pvnb(ji  ,jj) + pvnb(ji  ,jj-1) )   ) 
    12031210               ! 
    12041211            zv_trd(ji,jj) = - r1_4 * r1_e1e2v(ji,jj) * z1_hv                    & 
    1205                &               * (  e1e2t(ji,jj+1)*ht(ji,jj+1)*ff_t(ji,jj+1) * ( punb(ji,jj+1) + punb(ji-1,jj+1) )   &  
    1206                &                  + e1e2t(ji,jj  )*ht(ji,jj  )*ff_t(ji,jj  ) * ( punb(ji,jj  ) + punb(ji-1,jj  ) )   )  
     1212               &               * (  e1e2t(ji,jj+1)*pht(ji,jj+1)*ff_t(ji,jj+1) * ( punb(ji,jj+1) + punb(ji-1,jj+1) )   &  
     1213               &                  + e1e2t(ji,jj  )*pht(ji,jj  )*ff_t(ji,jj  ) * ( punb(ji,jj  ) + punb(ji-1,jj  ) )   )  
    12071214         END_2D 
    12081215         !          
  • NEMO/trunk/src/OCE/DYN/dynvor.F90

    r13226 r13237  
    8989   !! * Substitutions 
    9090#  include "do_loop_substitute.h90" 
     91#  include "domzgr_substitute.h90" 
     92 
    9193   !!---------------------------------------------------------------------- 
    9294   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    269271            DO_2D_01_01 
    270272               zwt(ji,jj) = r1_4 * (   zwz(ji-1,jj  ,jk) + zwz(ji,jj  ,jk)   & 
    271                   &                  + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 
     273                  &                  + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) ) & 
     274                  &                  * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 
    272275            END_2D 
    273276         CASE ( np_MET )                           !* metric term 
    274277            DO_2D_01_01 
    275278               zwt(ji,jj) = (   ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj)   & 
    276                   &           - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj)   ) * e3t(ji,jj,jk,Kmm) 
     279                  &           - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj)   ) & 
     280                  &             * e3t(ji,jj,jk,Kmm) 
    277281            END_2D 
    278282         CASE ( np_CRV )                           !* Coriolis + relative vorticity 
    279283            DO_2D_01_01 
    280284               zwt(ji,jj) = (  ff_t(ji,jj) + r1_4 * ( zwz(ji-1,jj  ,jk) + zwz(ji,jj  ,jk)    & 
    281                   &                                 + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) )  ) * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 
     285                  &                                 + zwz(ji-1,jj-1,jk) + zwz(ji,jj-1,jk) )  ) & 
     286                  &                                 * e1e2t(ji,jj)*e3t(ji,jj,jk,Kmm) 
    282287            END_2D 
    283288         CASE ( np_CME )                           !* Coriolis + metric 
     
    285290               zwt(ji,jj) = (  ff_t(ji,jj) * e1e2t(ji,jj)                           & 
    286291                    &        + ( pv(ji,jj,jk) + pv(ji,jj-1,jk) ) * di_e2u_2(ji,jj)  & 
    287                     &        - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj)  ) * e3t(ji,jj,jk,Kmm) 
     292                    &        - ( pu(ji,jj,jk) + pu(ji-1,jj,jk) ) * dj_e1v_2(ji,jj)  ) & 
     293                    &          * e3t(ji,jj,jk,Kmm) 
    288294            END_2D 
    289295         CASE DEFAULT                                             ! error 
     
    545551         CASE ( 0 )                                   ! original formulation  (masked averaging of e3t divided by 4) 
    546552            DO_2D_10_10 
    547                ze3f = (  e3t(ji,jj+1,jk,Kmm)*tmask(ji,jj+1,jk) + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
    548                   &    + e3t(ji,jj  ,jk,Kmm)*tmask(ji,jj  ,jk) + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
     553               ze3f = (  e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)   & 
     554                  &    + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
     555                  &    + e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)   & 
     556                  &    + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
    549557               IF( ze3f /= 0._wp ) THEN   ;   z1_e3f(ji,jj) = 4._wp / ze3f 
    550558               ELSE                       ;   z1_e3f(ji,jj) = 0._wp 
     
    553561         CASE ( 1 )                                   ! new formulation  (masked averaging of e3t divided by the sum of mask) 
    554562            DO_2D_10_10 
    555                ze3f = (  e3t(ji,jj+1,jk,Kmm)*tmask(ji,jj+1,jk) + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
    556                   &    + e3t(ji,jj  ,jk,Kmm)*tmask(ji,jj  ,jk) + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
     563               ze3f = (  e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)   & 
     564                  &    + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
     565                  &    + e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)   & 
     566                  &    + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
    557567               zmsk = (                    tmask(ji,jj+1,jk) +                     tmask(ji+1,jj+1,jk)   & 
    558568                  &                      + tmask(ji,jj  ,jk) +                     tmask(ji+1,jj  ,jk)  ) 
  • NEMO/trunk/src/OCE/DYN/dynzad.F90

    r12377 r13237  
    2929   !! * Substitutions 
    3030#  include "do_loop_substitute.h90" 
     31#  include "domzgr_substitute.h90" 
    3132   !!---------------------------------------------------------------------- 
    3233   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9596      ! 
    9697      DO_3D_00_00( 1, jpkm1 ) 
    97          puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zwuw(ji,jj,jk) + zwuw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
    98          pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zwvw(ji,jj,jk) + zwvw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     98         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zwuw(ji,jj,jk) + zwuw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
     99            &                                      / e3u(ji,jj,jk,Kmm) 
     100         pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zwvw(ji,jj,jk) + zwvw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj)   & 
     101            &                                      / e3v(ji,jj,jk,Kmm) 
    99102      END_3D 
    100103 
  • NEMO/trunk/src/OCE/DYN/dynzdf.F90

    r12489 r13237  
    3838   !! * Substitutions 
    3939#  include "do_loop_substitute.h90" 
     40#  include "domzgr_substitute.h90" 
    4041   !!---------------------------------------------------------------------- 
    4142   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    5556      !! ** Method  :  - Leap-Frog time stepping on all trends but the vertical mixing 
    5657      !!         u(after) =         u(before) + 2*dt *       u(rhs)                vector form or linear free surf. 
    57       !!         u(after) = ( e3u_b*u(before) + 2*dt * e3u_n*u(rhs) ) / e3u(after)   otherwise 
     58      !!         u(after) = ( e3u_b*u(before) + 2*dt * e3u_n*u(rhs) ) / e3u_after   otherwise 
    5859      !!               - update the after velocity with the implicit vertical mixing. 
    5960      !!      This requires to solver the following system:  
    60       !!         u(after) = u(after) + 1/e3u(after) dk+1[ mi(avm) / e3uw(after) dk[ua] ] 
     61      !!         u(after) = u(after) + 1/e3u_after  dk+1[ mi(avm) / e3uw_after dk[ua] ] 
    6162      !!      with the following surface/top/bottom boundary condition: 
    6263      !!      surface: wind stress input (averaged over kt-1/2 & kt+1/2) 
     
    113114         DO jk = 1, jpkm1 
    114115            puu(:,:,jk,Kaa) = (         e3u(:,:,jk,Kbb) * puu(:,:,jk,Kbb)  & 
    115                &          + rDt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs)  ) / e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     116               &            + rDt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs)  )   & 
     117               &                  / e3u(:,:,jk,Kaa) * umask(:,:,jk) 
    116118            pvv(:,:,jk,Kaa) = (         e3v(:,:,jk,Kbb) * pvv(:,:,jk,Kbb)  & 
    117                &          + rDt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs)  ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
     119               &            + rDt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs)  )   & 
     120               &                  / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    118121         END DO 
    119122      ENDIF 
     
    131134            iku = mbku(ji,jj)         ! ocean bottom level at u- and v-points  
    132135            ikv = mbkv(ji,jj)         ! (deepest ocean u- and v-points) 
    133             ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 
    134             ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 
     136            ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm)    & 
     137               &             + r_vvl   * e3u(ji,jj,iku,Kaa) 
     138            ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm)    & 
     139               &             + r_vvl   * e3v(ji,jj,ikv,Kaa) 
    135140            puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
    136141            pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
     
    140145               iku = miku(ji,jj)         ! top ocean level at u- and v-points  
    141146               ikv = mikv(ji,jj)         ! (first wet ocean u- and v-points) 
    142                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 
    143                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 
     147               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm)    & 
     148                  &             + r_vvl   * e3u(ji,jj,iku,Kaa) 
     149               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm)    & 
     150                  &             + r_vvl   * e3v(ji,jj,ikv,Kaa) 
    144151               puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 
    145152               pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 
     
    156163         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
    157164            DO_3D_00_00( 1, jpkm1 ) 
    158                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
     165               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm)    & 
     166                  &             + r_vvl   * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
    159167               zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) + akzu(ji,jj,jk  ) )   & 
    160168                  &         / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
     
    169177         CASE DEFAULT               ! iso-level lateral mixing 
    170178            DO_3D_00_00( 1, jpkm1 ) 
    171                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
    172                zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) ) / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
    173                zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 
     179               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm)    &    ! after scale factor at U-point 
     180                  &             + r_vvl   * e3u(ji,jj,jk,Kaa) 
     181               zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) )   & 
     182                  &         / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
     183               zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) )   & 
     184                  &         / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 
    174185               zWui = ( wi(ji,jj,jk  ) + wi(ji+1,jj,jk  ) ) / ze3ua 
    175186               zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua 
     
    181192         DO_2D_00_00 
    182193            zwi(ji,jj,1) = 0._wp 
    183             ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa) 
    184             zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji  ,jj,2) ) / ( ze3ua * e3uw(ji,jj,2,Kmm) ) * wumask(ji,jj,2) 
     194            ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm)    & 
     195               &             + r_vvl   * e3u(ji,jj,1,Kaa) 
     196            zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji  ,jj,2) )   & 
     197               &         / ( ze3ua * e3uw(ji,jj,2,Kmm) ) * wumask(ji,jj,2) 
    185198            zWus = ( wi(ji  ,jj,2) +  wi(ji+1,jj,2) ) / ze3ua 
    186199            zws(ji,jj,1 ) = zzws - zdt * MAX( zWus, 0._wp ) 
     
    191204         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
    192205            DO_3D_00_00( 1, jpkm1 ) 
    193                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
     206               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm)    & 
     207                  &             + r_vvl   * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
    194208               zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) + akzu(ji,jj,jk  ) )   & 
    195209                  &         / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
     
    202216         CASE DEFAULT               ! iso-level lateral mixing 
    203217            DO_3D_00_00( 1, jpkm1 ) 
    204                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
    205                zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) ) / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
    206                zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 
     218               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm)    & 
     219                  &             + r_vvl   * e3u(ji,jj,jk,Kaa)   ! after scale factor at U-point 
     220               zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) )    & 
     221                  &         / ( ze3ua * e3uw(ji,jj,jk  ,Kmm) ) * wumask(ji,jj,jk  ) 
     222               zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) )    & 
     223                  &         / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 
    207224               zwi(ji,jj,jk) = zzwi 
    208225               zws(ji,jj,jk) = zzws 
     
    226243         DO_2D_00_00 
    227244            iku = mbku(ji,jj)       ! ocean bottom level at u- and v-points 
    228             ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
     245            ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm)    & 
     246               &             + r_vvl   * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
    229247            zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 
    230248         END_2D 
     
    233251               !!gm   top Cd is masked (=0 outside cavities) no need of test on mik>=2  ==>> it has been suppressed 
    234252               iku = miku(ji,jj)       ! ocean top level at u- and v-points  
    235                ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
     253               ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm)    & 
     254                  &             + r_vvl   * e3u(ji,jj,iku,Kaa)   ! after scale factor at T-point 
    236255               zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 
    237256            END_2D 
     
    259278      ! 
    260279      DO_2D_00_00 
    261          ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa)  
     280         ze3ua =  ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm)    & 
     281            &             + r_vvl   * e3u(ji,jj,1,Kaa)  
    262282         puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) )   & 
    263283            &                                      / ( ze3ua * rho0 ) * umask(ji,jj,1)  
     
    282302         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzv) 
    283303            DO_3D_00_00( 1, jpkm1 ) 
    284                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
     304               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm)    & 
     305                  &             + r_vvl   * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
    285306               zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) + akzv(ji,jj,jk  ) )   & 
    286307                  &         / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
     
    295316         CASE DEFAULT               ! iso-level lateral mixing 
    296317            DO_3D_00_00( 1, jpkm1 ) 
    297                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
    298                zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) ) / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
    299                zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 
     318               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm)    & 
     319                  &             + r_vvl   * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
     320               zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) )    & 
     321                  &         / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
     322               zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) )    & 
     323                  &         / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 
    300324               zWvi = ( wi(ji,jj,jk  ) + wi(ji,jj+1,jk  ) ) / ze3va 
    301325               zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va 
     
    307331         DO_2D_00_00 
    308332            zwi(ji,jj,1) = 0._wp 
    309             ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa) 
    310             zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) ) / ( ze3va * e3vw(ji,jj,2,Kmm) ) * wvmask(ji,jj,2) 
     333            ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm)    & 
     334               &             + r_vvl   * e3v(ji,jj,1,Kaa) 
     335            zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) )    & 
     336               &         / ( ze3va * e3vw(ji,jj,2,Kmm) ) * wvmask(ji,jj,2) 
    311337            zWvs = ( wi(ji,jj  ,2) +  wi(ji,jj+1,2) ) / ze3va 
    312338            zws(ji,jj,1 ) = zzws - zdt * MAX( zWvs, 0._wp ) 
     
    317343         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
    318344            DO_3D_00_00( 1, jpkm1 ) 
    319                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
     345               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm)    & 
     346                  &             + r_vvl   * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
    320347               zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) + akzv(ji,jj,jk  ) )   & 
    321348                  &         / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
     
    328355         CASE DEFAULT               ! iso-level lateral mixing 
    329356            DO_3D_00_00( 1, jpkm1 ) 
    330                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
    331                zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) ) / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
    332                zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 
     357               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm)    & 
     358                  &             + r_vvl   * e3v(ji,jj,jk,Kaa)   ! after scale factor at V-point 
     359               zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) )    & 
     360                  &         / ( ze3va * e3vw(ji,jj,jk  ,Kmm) ) * wvmask(ji,jj,jk  ) 
     361               zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) )    & 
     362                  &         / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 
    333363               zwi(ji,jj,jk) = zzwi 
    334364               zws(ji,jj,jk) = zzws 
     
    351381         DO_2D_00_00 
    352382            ikv = mbkv(ji,jj)       ! (deepest ocean u- and v-points) 
    353             ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
     383            ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm)    & 
     384               &             + r_vvl   * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
    354385            zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va            
    355386         END_2D 
     
    357388            DO_2D_00_00 
    358389               ikv = mikv(ji,jj)       ! (first wet ocean u- and v-points) 
    359                ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
     390               ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm)    & 
     391                  &             + r_vvl   * e3v(ji,jj,ikv,Kaa)   ! after scale factor at T-point 
    360392               zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 
    361393            END_2D 
     
    383415      ! 
    384416      DO_2D_00_00 
    385          ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa)  
     417         ze3va =  ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm)    & 
     418            &             + r_vvl   * e3v(ji,jj,1,Kaa)  
    386419         pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) )   & 
    387420            &                                      / ( ze3va * rho0 ) * vmask(ji,jj,1)  
  • NEMO/trunk/src/OCE/DYN/sshwzv.F90

    r13226 r13237  
    5050   !! * Substitutions 
    5151#  include "do_loop_substitute.h90" 
     52#  include "domzgr_substitute.h90" 
     53 
    5254   !!---------------------------------------------------------------------- 
    5355   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    110112      ! 
    111113#if defined key_agrif 
    112       Kbb_a = Kbb; Kmm_a = Kmm; Krhs_a = Kaa; CALL agrif_ssh( kt ) 
     114      Kbb_a = Kbb   ;   Kmm_a = Kmm   ;   Krhs_a = Kaa 
     115      CALL agrif_ssh( kt ) 
    113116#endif 
    114117      ! 
     
    130133 
    131134    
    132    SUBROUTINE wzv( kt, Kbb, Kmm, pww, Kaa ) 
     135   SUBROUTINE wzv( kt, Kbb, Kmm, Kaa, pww ) 
    133136      !!---------------------------------------------------------------------- 
    134137      !!                ***  ROUTINE wzv  *** 
     
    147150      INTEGER                         , INTENT(in)    ::   kt             ! time step 
    148151      INTEGER                         , INTENT(in)    ::   Kbb, Kmm, Kaa  ! time level indices 
    149       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pww            ! now vertical velocity 
     152      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pww            ! vertical velocity at Kmm 
    150153      ! 
    151154      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     
    166169      !                                           !------------------------------! 
    167170      ! 
    168       IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN      ! z_tilde and layer cases 
     171      !                                               !===============================! 
     172      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN      !==  z_tilde and layer cases  ==! 
     173         !                                            !===============================! 
    169174         ALLOCATE( zhdiv(jpi,jpj,jpk) )  
    170175         ! 
     
    181186         DO jk = jpkm1, 1, -1                       ! integrate from the bottom the hor. divergence 
    182187            ! computation of w 
    183             pww(:,:,jk) = pww(:,:,jk+1) - (  e3t(:,:,jk,Kmm) * hdiv(:,:,jk) + zhdiv(:,:,jk)    & 
    184                &                         + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )     ) * tmask(:,:,jk) 
     188            pww(:,:,jk) = pww(:,:,jk+1) - (   e3t(:,:,jk,Kmm) * hdiv(:,:,jk)   & 
     189               &                            +                  zhdiv(:,:,jk)   & 
     190               &                            + r1_Dt * (  e3t(:,:,jk,Kaa)       & 
     191               &                                       - e3t(:,:,jk,Kbb) )   ) * tmask(:,:,jk) 
    185192         END DO 
    186193         !          IF( ln_vvl_layer ) pww(:,:,:) = 0.e0 
    187194         DEALLOCATE( zhdiv )  
    188       ELSE   ! z_star and linear free surface cases 
     195         !                                            !=================================! 
     196      ELSEIF( ln_linssh )   THEN                      !==  linear free surface cases  ==! 
     197         !                                            !=================================! 
     198         DO jk = jpkm1, 1, -1                               ! integrate from the bottom the hor. divergence 
     199            pww(:,:,jk) = pww(:,:,jk+1) - (  e3t(:,:,jk,Kmm) * hdiv(:,:,jk)  ) * tmask(:,:,jk) 
     200         END DO 
     201         !                                            !==========================================! 
     202      ELSE                                            !==  Quasi-Eulerian vertical coordinate  ==!   ('key_qco') 
     203         !                                            !==========================================! 
    189204         DO jk = jpkm1, 1, -1                       ! integrate from the bottom the hor. divergence 
    190             ! computation of w 
    191205            pww(:,:,jk) = pww(:,:,jk+1) - (  e3t(:,:,jk,Kmm) * hdiv(:,:,jk)                 & 
    192                &                         + r1_Dt * ( e3t(:,:,jk,Kaa) - e3t(:,:,jk,Kbb) )  ) * tmask(:,:,jk) 
     206               &                            + r1_Dt * (  e3t(:,:,jk,Kaa)        & 
     207               &                                       - e3t(:,:,jk,Kbb)  )   ) * tmask(:,:,jk) 
    193208         END DO 
    194209      ENDIF 
     
    248263 
    249264 
    250    SUBROUTINE ssh_atf( kt, Kbb, Kmm, Kaa, pssh ) 
     265   SUBROUTINE ssh_atf( kt, Kbb, Kmm, Kaa, pssh, pssh_f ) 
    251266      !!---------------------------------------------------------------------- 
    252267      !!                    ***  ROUTINE ssh_atf  *** 
     
    265280      INTEGER                         , INTENT(in   ) ::   kt             ! ocean time-step index 
    266281      INTEGER                         , INTENT(in   ) ::   Kbb, Kmm, Kaa  ! ocean time level indices 
    267       REAL(wp), DIMENSION(jpi,jpj,jpt), INTENT(inout) ::   pssh           ! SSH field 
     282      REAL(wp), DIMENSION(jpi,jpj,jpt)          , TARGET, INTENT(inout) ::   pssh           ! SSH field 
     283      REAL(wp), DIMENSION(jpi,jpj    ), OPTIONAL, TARGET, INTENT(  out) ::   pssh_f         ! filtered SSH field 
    268284      ! 
    269285      REAL(wp) ::   zcoef   ! local scalar 
     286      REAL(wp), POINTER, DIMENSION(:,:) ::   zssh   ! pointer for filtered SSH  
    270287      !!---------------------------------------------------------------------- 
    271288      ! 
     
    279296      !              !==  Euler time-stepping: no filter, just swap  ==! 
    280297      IF ( .NOT.( l_1st_euler ) ) THEN   ! Only do time filtering for leapfrog timesteps 
     298         IF( PRESENT( pssh_f ) ) THEN   ;   zssh => pssh_f 
     299         ELSE                           ;   zssh => pssh(:,:,Kmm) 
     300         ENDIF 
    281301         !                                                  ! filtered "now" field 
    282302         pssh(:,:,Kmm) = pssh(:,:,Kmm) + rn_atfp * ( pssh(:,:,Kbb) - 2 * pssh(:,:,Kmm) + pssh(:,:,Kaa) ) 
     
    300320   END SUBROUTINE ssh_atf 
    301321 
     322    
    302323   SUBROUTINE wAimp( kt, Kmm ) 
    303324      !!---------------------------------------------------------------------- 
     
    320341      ! 
    321342      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    322       REAL(wp)             ::   zCu, zcff, z1_e3t                     ! local scalars 
     343      REAL(wp)             ::   zCu, zcff, z1_e3t, zdt                ! local scalars 
    323344      REAL(wp) , PARAMETER ::   Cu_min = 0.15_wp                      ! local parameters 
    324345      REAL(wp) , PARAMETER ::   Cu_max = 0.30_wp                      ! local parameters 
     
    337358      ! 
    338359      ! Calculate Courant numbers 
     360      zdt = 2._wp * rn_Dt                            ! 2*rn_Dt and not rDt (for restartability) 
    339361      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 
    340362         DO_3D_00_00( 1, jpkm1 ) 
    341363            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    342             ! 2*rn_Dt and not rDt (for restartability) 
    343             Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )                       &   
    344                &                             + ( MAX( e2u(ji  ,jj)*e3u(ji  ,jj,jk,Kmm)*uu(ji  ,jj,jk,Kmm) + un_td(ji  ,jj,jk), 0._wp ) -   & 
    345                &                                 MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm) + un_td(ji-1,jj,jk), 0._wp ) )   & 
     364            Cu_adv(ji,jj,jk) =   zdt *                                                         & 
     365               &  ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )            & 
     366               &  + ( MAX( e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm)                                  & 
     367               &                        * uu (ji  ,jj,jk,Kmm) + un_td(ji  ,jj,jk), 0._wp ) -   & 
     368               &      MIN( e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm)                                  & 
     369               &                        * uu (ji-1,jj,jk,Kmm) + un_td(ji-1,jj,jk), 0._wp ) )   & 
    346370               &                               * r1_e1e2t(ji,jj)                                                                     & 
    347                &                             + ( MAX( e1v(ji,jj  )*e3v(ji,jj  ,jk,Kmm)*vv(ji,jj  ,jk,Kmm) + vn_td(ji,jj  ,jk), 0._wp ) -   & 
    348                &                                 MIN( e1v(ji,jj-1)*e3v(ji,jj-1,jk,Kmm)*vv(ji,jj-1,jk,Kmm) + vn_td(ji,jj-1,jk), 0._wp ) )   & 
     371               &  + ( MAX( e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm)                                  & 
     372               &                        * vv (ji,jj  ,jk,Kmm) + vn_td(ji,jj  ,jk), 0._wp ) -   & 
     373               &      MIN( e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm)                                  & 
     374               &                        * vv (ji,jj-1,jk,Kmm) + vn_td(ji,jj-1,jk), 0._wp ) )   & 
    349375               &                               * r1_e1e2t(ji,jj)                                                                     & 
    350376               &                             ) * z1_e3t 
     
    353379         DO_3D_00_00( 1, jpkm1 ) 
    354380            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    355             ! 2*rn_Dt and not rDt (for restartability) 
    356             Cu_adv(ji,jj,jk) = 2._wp * rn_Dt * ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )   &  
     381            Cu_adv(ji,jj,jk) =   zdt *                                                      & 
     382               &  ( ( MAX( ww(ji,jj,jk) , 0._wp ) - MIN( ww(ji,jj,jk+1) , 0._wp ) )         & 
    357383               &                             + ( MAX( e2u(ji  ,jj)*e3u(ji  ,jj,jk,Kmm)*uu(ji  ,jj,jk,Kmm), 0._wp ) -   & 
    358384               &                                 MIN( e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kmm)*uu(ji-1,jj,jk,Kmm), 0._wp ) )   & 
  • NEMO/trunk/src/OCE/DYN/wet_dry.F90

    r13226 r13237  
    3333   !! * Substitutions 
    3434#  include "do_loop_substitute.h90" 
     35#  include "domzgr_substitute.h90" 
    3536   !!---------------------------------------------------------------------- 
    3637   !! critical depths,filters, limiters,and masks for  Wetting and Drying 
  • NEMO/trunk/src/OCE/FLO/flo4rk.F90

    r12489 r13237  
    2626   REAL(wp), DIMENSION (3) ::   scoef1 = (/  0.5  ,  0.5  ,  1.0  /)           ! 
    2727 
     28#  include "domzgr_substitute.h90" 
    2829   !!---------------------------------------------------------------------- 
    2930   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/FLO/floblk.F90

    r13216 r13237  
    1919 
    2020   PUBLIC   flo_blk    ! routine called by floats.F90 
     21 
     22#  include "domzgr_substitute.h90" 
    2123 
    2224   !!---------------------------------------------------------------------- 
     
    115117            ! compute the transport across the mesh where the float is.             
    116118!!bug (gm) change e3t into e3. but never checked  
    117             zsurfx(1) = e2u(iiloc(jfl)-1,ijloc(jfl)  ) * e3u(iiloc(jfl)-1,ijloc(jfl)  ,-ikl(jfl),Kmm) 
    118             zsurfx(2) = e2u(iiloc(jfl)  ,ijloc(jfl)  ) * e3u(iiloc(jfl)  ,ijloc(jfl)  ,-ikl(jfl),Kmm) 
    119             zsurfy(1) = e1v(iiloc(jfl)  ,ijloc(jfl)-1) * e3v(iiloc(jfl)  ,ijloc(jfl)-1,-ikl(jfl),Kmm) 
    120             zsurfy(2) = e1v(iiloc(jfl)  ,ijloc(jfl)  ) * e3v(iiloc(jfl)  ,ijloc(jfl)  ,-ikl(jfl),Kmm) 
     119            zsurfx(1) =   & 
     120            &   e2u(iiloc(jfl)-1,ijloc(jfl)  )    & 
     121            & * e3u(iiloc(jfl)-1,ijloc(jfl)  ,-ikl(jfl),Kmm) 
     122            zsurfx(2) =   & 
     123            &   e2u(iiloc(jfl)  ,ijloc(jfl)  )    & 
     124            & * e3u(iiloc(jfl)  ,ijloc(jfl)  ,-ikl(jfl),Kmm) 
     125            zsurfy(1) =   & 
     126            &   e1v(iiloc(jfl)  ,ijloc(jfl)-1)    & 
     127            & * e3v(iiloc(jfl)  ,ijloc(jfl)-1,-ikl(jfl),Kmm) 
     128            zsurfy(2) =   & 
     129            &   e1v(iiloc(jfl)  ,ijloc(jfl)  )    & 
     130            & * e3v(iiloc(jfl)  ,ijloc(jfl)  ,-ikl(jfl),Kmm) 
    121131 
    122132            ! for a isobar float zsurfz is put to zero. The vertical velocity will be zero too. 
  • NEMO/trunk/src/OCE/IOM/restart.F90

    r12489 r13237  
    295295         vv   (:,:,:  ,Kbb) = vv   (:,:,:  ,Kmm) 
    296296         ssh  (:,:    ,Kbb) = ssh  (:,:    ,Kmm) 
    297          ! 
    298          IF( .NOT.ln_linssh ) THEN 
    299             DO jk = 1, jpk 
    300                e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 
    301             END DO 
    302          ENDIF 
    303          ! 
    304297      ENDIF 
    305298      ! 
  • NEMO/trunk/src/OCE/ISF/isfcavgam.F90

    r12077 r13237  
    3030   PUBLIC   isfcav_gammats 
    3131 
     32#  include "domzgr_substitute.h90" 
    3233   !!---------------------------------------------------------------------- 
    3334   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/ISF/isfcpl.F90

    r13226 r13237  
    1515   USE isfutils, ONLY : debug 
    1616   USE lib_mpp , ONLY: mpp_sum, mpp_max ! mpp routine 
     17#if ! defined key_qco 
    1718   USE domvvl  , ONLY: dom_vvl_zgr      ! vertical scale factor interpolation 
     19#else 
     20   USE domqco   , ONLY: dom_qco_zgr      ! vertical scale factor interpolation 
     21#endif 
    1822   USE domngb  , ONLY: dom_ngb          ! find the closest grid point from a given lon/lat position 
    1923   ! 
     
    4347   !! * Substitutions 
    4448#  include "do_loop_substitute.h90" 
     49#  include "domzgr_substitute.h90" 
    4550   !!---------------------------------------------------------------------- 
    4651   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    112117      vv   (:,:,:,Kbb)   = vv   (:,:,:,Kmm) 
    113118      ssh (:,:,Kbb)     = ssh (:,:,Kmm) 
     119#if ! defined key_qco 
    114120      e3t(:,:,:,Kbb)   = e3t(:,:,:,Kmm) 
    115   
     121#endif  
    116122      ! prepare writing restart 
    117123      IF( lwxios ) THEN 
     
    135141      INTEGER, INTENT(in) :: Kmm    ! ocean time level index 
    136142      !!---------------------------------------------------------------------- 
     143      INTEGER :: jk                               ! loop index 
     144      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t, ze3u, ze3v, zgdepw  ! e3t , e3u, e3v !!st patch to use substitution 
     145      !!---------------------------------------------------------------------- 
     146      ! 
     147      DO jk = 1, jpk 
     148         ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     149         ze3u(:,:,jk) = e3u(:,:,jk,Kmm) 
     150         ze3v(:,:,jk) = e3v(:,:,jk,Kmm) 
     151         ! 
     152         zgdepw(:,:,jk) = gdepw(:,:,jk,Kmm) 
     153      END DO  
    137154      ! 
    138155      IF( lwxios ) CALL iom_swap( cwxios_context ) 
    139156      CALL iom_rstput( kt, nitrst, numrow, 'tmask'  , tmask , ldxios = lwxios ) 
    140157      CALL iom_rstput( kt, nitrst, numrow, 'ssmask' , ssmask, ldxios = lwxios ) 
    141       CALL iom_rstput( kt, nitrst, numrow, 'e3t_n'  , e3t(:,:,:,Kmm) , ldxios = lwxios ) 
    142       CALL iom_rstput( kt, nitrst, numrow, 'e3u_n'  , e3u(:,:,:,Kmm) , ldxios = lwxios ) 
    143       CALL iom_rstput( kt, nitrst, numrow, 'e3v_n'  , e3v(:,:,:,Kmm) , ldxios = lwxios ) 
    144       CALL iom_rstput( kt, nitrst, numrow, 'gdepw_n', gdepw(:,:,:,Kmm) , ldxios = lwxios ) 
     158      CALL iom_rstput( kt, nitrst, numrow, 'e3t_n'  , ze3t , ldxios = lwxios ) 
     159      CALL iom_rstput( kt, nitrst, numrow, 'e3u_n'  , ze3u , ldxios = lwxios ) 
     160      CALL iom_rstput( kt, nitrst, numrow, 'e3v_n'  , ze3v , ldxios = lwxios ) 
     161      CALL iom_rstput( kt, nitrst, numrow, 'gdepw_n', zgdepw , ldxios = lwxios ) 
    145162      IF( lwxios ) CALL iom_swap( cxios_context ) 
    146163      ! 
     
    209226      IF(lwp) write(numout,*) 'isfcpl_ssh : recompute scale factor from ssh (new wet cell,Kmm)' 
    210227      IF(lwp) write(numout,*) '~~~~~~~~~~~' 
     228#if ! defined key_qco 
    211229      DO jk = 1, jpk 
    212          e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & 
    213              &                          / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    214              &          + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk)) 
     230         e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + (ht_0(:,:) + ssh(:,:,Kmm)) * r1_ht_0(:,:) ) 
    215231      END DO 
    216232      e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    217233      CALL dom_vvl_zgr(Kbb, Kmm, Kaa) 
     234#else 
     235      CALL dom_qco_zgr(Kbb, Kmm, Kaa) 
     236#endif 
    218237      ! 
    219238   END SUBROUTINE isfcpl_ssh 
     
    400419         ! 1.1: get volume flux before coupling (>0 out) 
    401420         DO_2D_00_00 
    402             zqvolb(ji,jj,jk) =  (   e2u(ji,jj) * ze3u_b(ji,jj,jk) * uu(ji,jj,jk,Kmm) - e2u(ji-1,jj  ) * ze3u_b(ji-1,jj  ,jk) * uu(ji-1,jj  ,jk,Kmm)    & 
    403                &                  + e1v(ji,jj) * ze3v_b(ji,jj,jk) * vv(ji,jj,jk,Kmm) - e1v(ji  ,jj-1) * ze3v_b(ji  ,jj-1,jk) * vv(ji  ,jj-1,jk,Kmm)  ) & 
     421            zqvolb(ji,jj,jk) =    & 
     422               &  (   e2u(ji  ,jj  ) * ze3u_b(ji  ,jj  ,jk) * uu(ji  ,jj  ,jk,Kmm)      & 
     423               &    - e2u(ji-1,jj  ) * ze3u_b(ji-1,jj  ,jk) * uu(ji-1,jj  ,jk,Kmm)      & 
     424               &    + e1v(ji  ,jj  ) * ze3v_b(ji  ,jj  ,jk) * vv(ji  ,jj  ,jk,Kmm)      & 
     425               &    - e1v(ji  ,jj-1) * ze3v_b(ji  ,jj-1,jk) * vv(ji  ,jj-1,jk,Kmm)  )   & 
    404426               &                * ztmask_b(ji,jj,jk) 
    405427         END_2D 
     
    412434         ! compute volume flux divergence after coupling 
    413435         DO_2D_00_00 
    414             zqvoln(ji,jj,jk) = (   e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) - e2u(ji-1,jj  ) * e3u(ji-1,jj  ,jk,Kmm) * uu(ji-1,jj  ,jk,Kmm)    & 
    415                &                 + e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) - e1v(ji  ,jj-1) * e3v(ji  ,jj-1,jk,Kmm) * vv(ji  ,jj-1,jk,Kmm)  ) & 
     436            zqvoln(ji,jj,jk) =   & 
     437               &  (   e2u(ji  ,jj  ) * e3u(ji  ,jj  ,jk,Kmm) * uu(ji  ,jj  ,jk,Kmm)    & 
     438               &    - e2u(ji-1,jj  ) * e3u(ji-1,jj  ,jk,Kmm) * uu(ji-1,jj  ,jk,Kmm)    & 
     439               &    + e1v(ji  ,jj  ) * e3v(ji  ,jj  ,jk,Kmm) * vv(ji  ,jj  ,jk,Kmm)    & 
     440               &    - e1v(ji  ,jj-1) * e3v(ji  ,jj-1,jk,Kmm) * vv(ji  ,jj-1,jk,Kmm)  ) & 
    416441               &               * tmask(ji,jj,jk) 
    417442         END_2D 
     
    523548 
    524549               ! volume diff 
    525                zdvol = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) - ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
     550               zdvol =   e3t  (ji,jj,jk,Kmm) *  tmask  (ji,jj,jk)   & 
     551                  &   - ze3t_b(ji,jj,jk    ) * ztmask_b(ji,jj,jk) 
    526552 
    527553               ! heat diff 
     
    555581            DO ji = nldi,nlei 
    556582               jip1=MIN(ji+1,jpi) ; jim1=MAX(ji-1,1) ; jjp1=MIN(jj+1,jpj) ; jjm1=MAX(jj-1,1) ; 
    557                IF ( tmask(ji,jj,jk) == 0._wp .AND. ztmask_b(ji,jj,jk) == 1._wp ) nisfl(narea) = nisfl(narea) + MAX(SUM(tmask(jim1:jip1,jjm1:jjp1,jk)),1._wp) 
     583               IF ( tmask(ji,jj,jk) == 0._wp .AND. ztmask_b(ji,jj,jk) == 1._wp ) THEN  
     584                  nisfl(narea) = nisfl(narea) + MAX(SUM(tmask(jim1:jip1,jjm1:jjp1,jk)),1._wp) 
     585               ENDIF 
    558586            ENDDO 
    559587         ENDDO 
  • NEMO/trunk/src/OCE/ISF/isfdiags.F90

    r12905 r13237  
    2626   !! * Substitutions 
    2727#  include "do_loop_substitute.h90" 
     28#  include "domzgr_substitute.h90" 
    2829   !!---------------------------------------------------------------------- 
    2930   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/ISF/isfdynatf.F90

    r12489 r13237  
    1414 
    1515   USE phycst , ONLY: r1_rho0         ! physical constant 
    16    USE dom_oce, ONLY: tmask, ssmask, ht, e3t, r1_e1e2t   ! time and space domain 
     16   USE dom_oce                        ! time and space domain 
     17   USE oce, ONLY : ssh                ! sea-surface height !!st needed for substitution 
    1718 
    1819   USE in_out_manager 
     
    2526   !! * Substitutions 
    2627#  include "do_loop_substitute.h90" 
     28#  include "domzgr_substitute.h90" 
    2729 
    2830CONTAINS 
     
    8183      ! add the increment 
    8284      DO jk = 1, jpkm1 
    83          pe3t_f(:,:,jk) = pe3t_f(:,:,jk) - tmask(:,:,jk) * zfwfinc(:,:) * e3t(:,:,jk,Kmm) 
     85         pe3t_f(:,:,jk) = pe3t_f(:,:,jk) - tmask(:,:,jk) * zfwfinc(:,:)   & 
     86            &                              * e3t(:,:,jk,Kmm) 
    8487      END DO 
    8588      ! 
  • NEMO/trunk/src/OCE/ISF/isfhdiv.F90

    r12489 r13237  
    2626   !! * Substitutions 
    2727#  include "do_loop_substitute.h90" 
     28#  include "domzgr_substitute.h90" 
    2829 
    2930CONTAINS 
     
    134135      ! 
    135136      DO jk=1,jpk  
    136          phdiv(:,:,jk) =  phdiv(:,:,jk) + pqvol(:,:,jk) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) 
     137         phdiv(:,:,jk) =  phdiv(:,:,jk) + pqvol(:,:,jk) * r1_e1e2t(:,:)   & 
     138            &                             / e3t(:,:,jk,Kmm) 
    137139      END DO 
    138140      ! 
  • NEMO/trunk/src/OCE/ISF/isfload.F90

    r12340 r13237  
    1313   USE isf_oce, ONLY: cn_isfload, rn_isfload_T, rn_isfload_S ! ice shelf variables 
    1414 
    15    USE dom_oce, ONLY: e3w, gdept, risfdep, mikt     ! vertical scale factor 
     15   USE dom_oce                                      ! vertical scale factor 
    1616   USE eosbn2 , ONLY: eos                           ! eos routine 
    1717 
     
    2626   !! * Substitutions 
    2727#  include "do_loop_substitute.h90" 
     28#  include "domzgr_substitute.h90" 
    2829 
    2930CONTAINS 
     
    99100            ! 
    100101            ! top layer of the ice shelf 
    101             pisfload(ji,jj) = pisfload(ji,jj) + (znad + zrhd(ji,jj,1) ) * e3w(ji,jj,1,Kmm) 
     102            pisfload(ji,jj) = pisfload(ji,jj) + (znad + zrhd(ji,jj,1) )   & 
     103               &                                * e3w(ji,jj,1,Kmm) 
    102104            ! 
    103105            ! core layers of the ice shelf 
    104106            DO jk = 2, ikt-1 
    105                pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) * e3w(ji,jj,jk,Kmm) 
     107               pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk))   & 
     108                  &                                * e3w(ji,jj,jk,Kmm) 
    106109            END DO 
    107110            ! 
  • NEMO/trunk/src/OCE/ISF/isfstp.F90

    r12242 r13237  
    1313   !!   isfstp       : compute iceshelf melt and heat flux 
    1414   !!---------------------------------------------------------------------- 
    15    ! 
    1615   USE isf_oce                                      ! isf variables 
    1716   USE isfload, ONLY: isf_load                      ! ice shelf load 
     
    2120   USE isfcpl , ONLY: isfcpl_rst_write, isfcpl_init ! isf variables 
    2221 
    23    USE dom_oce, ONLY: ht, e3t, ln_isfcav, ln_linssh     ! ocean space and time domain 
     22   USE dom_oce        ! ocean space and time domain 
     23   USE oce      , ONLY: ssh                           ! sea surface height 
    2424   USE domvvl,  ONLY: ln_vvl_zstar                      ! zstar logical 
    2525   USE zdfdrg,  ONLY: r_Cdmin_top, r_ke0_top            ! vertical physics: top/bottom drag coef. 
     
    3131 
    3232   IMPLICIT NONE 
    33  
    3433   PRIVATE 
    3534 
    3635   PUBLIC   isf_stp, isf_init, isf_nam  ! routine called in sbcmod and divhor 
    3736 
     37   !! * Substitutions 
     38#  include "domzgr_substitute.h90" 
    3839   !!---------------------------------------------------------------------- 
    3940   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    4142   !! Software governed by the CeCILL license (see ./LICENSE) 
    4243   !!---------------------------------------------------------------------- 
     44 
    4345CONTAINS 
    4446  
     
    6062      INTEGER, INTENT(in) ::   kt   ! ocean time step 
    6163      INTEGER, INTENT(in) ::   Kmm  ! ocean time level index 
     64      !!---------------------------------------------------------------------- 
     65      INTEGER :: jk                               ! loop index 
     66      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t    ! e3t  
    6267      !!--------------------------------------------------------------------- 
    6368      ! 
     
    7883         ! 1.2: compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 
    7984         rhisf_tbl_cav(:,:) = rn_htbl * mskisf_cav(:,:) 
    80          CALL isf_tbl_lvl(ht, e3t(:,:,:,Kmm), misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav) 
     85         DO jk = 1, jpk 
     86            ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     87         END DO  
     88         CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav) 
    8189         ! 
    8290         ! 1.3: compute ice shelf melt 
     
    100108         ! by simplicity, we assume the top level where param applied do not change with time (done in init part) 
    101109         rhisf_tbl_par(:,:) = rhisf0_tbl_par(:,:) 
    102          CALL isf_tbl_lvl(ht, e3t(:,:,:,Kmm), misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par) 
     110         DO jk = 1, jpk 
     111            ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     112         END DO 
     113         CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par) 
    103114         ! 
    104115         ! 2.3: compute ice shelf melt 
  • NEMO/trunk/src/OCE/ISF/isftbl.F90

    r12340 r13237  
    2525   !! * Substitutions 
    2626#  include "do_loop_substitute.h90" 
     27#  include "domzgr_substitute.h90" 
    2728 
    2829CONTAINS 
     
    5657      REAL(wp), DIMENSION(jpi,jpj) :: zhtbl   ! thickness of the tbl 
    5758      REAL(wp), DIMENSION(jpi,jpj) :: zfrac   ! thickness of the tbl 
     59      INTEGER :: jk                            ! loop index 
     60      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t,ze3u,ze3v ! e3  
    5861      !!-------------------------------------------------------------------- 
    5962      !  
     
    6467         zhtbl = phtbl 
    6568         ! 
     69         DO jk = 1, jpk 
     70            ze3u(:,:,jk) = e3u(:,:,jk,Kmm) 
     71         END DO  
    6672         ! compute tbl lvl and thickness 
    67          CALL isf_tbl_lvl( hu(:,:,Kmm), e3u(:,:,:,Kmm), ktop, ikbot, zhtbl, zfrac ) 
     73         CALL isf_tbl_lvl( hu(:,:,Kmm), ze3u, ktop, ikbot, zhtbl, zfrac ) 
    6874         ! 
    6975         ! compute tbl property at U point 
    70          CALL isf_tbl_avg( miku, ikbot, zhtbl, zfrac, e3u(:,:,:,Kmm), pvarin, zvarout ) 
     76         CALL isf_tbl_avg( miku, ikbot, zhtbl, zfrac, ze3u, pvarin, zvarout ) 
    7177         ! 
    7278         ! compute tbl property at T point 
     
    8288         zhtbl = phtbl 
    8389         ! 
     90         DO jk = 1, jpk 
     91            ze3v(:,:,jk) = e3v(:,:,jk,Kmm) 
     92         END DO  
    8493         ! compute tbl lvl and thickness 
    85          CALL isf_tbl_lvl( hv(:,:,Kmm), e3v(:,:,:,Kmm), ktop, ikbot, zhtbl, zfrac ) 
     94         CALL isf_tbl_lvl( hv(:,:,Kmm), ze3v, ktop, ikbot, zhtbl, zfrac ) 
    8695         ! 
    8796         ! compute tbl property at V point 
    88          CALL isf_tbl_avg( mikv, ikbot, zhtbl, zfrac, e3v(:,:,:,Kmm), pvarin, zvarout ) 
     97         CALL isf_tbl_avg( mikv, ikbot, zhtbl, zfrac, ze3v, pvarin, zvarout ) 
    8998         ! 
    9099         ! pvarout is an averaging of wet point 
     
    98107         ! 
    99108         ! compute tbl property at T point 
    100          CALL isf_tbl_avg( ktop, kbot, phtbl, pfrac, e3t(:,:,:,Kmm), pvarin, pvarout ) 
     109         DO jk = 1, jpk 
     110            ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 
     111         END DO  
     112         CALL isf_tbl_avg( ktop, kbot, phtbl, pfrac, ze3t, pvarin, pvarout ) 
    101113         ! 
    102114      END SELECT 
     
    212224      ! phtbl need to be bounded by water column thickness before 
    213225      ! test: if htbl = water column thickness, should return mbathy 
    214       ! test: if htbl = 0 should return ktop (phtbl cap to e3t(ji,jj,1)) 
     226      ! test: if htbl = 0 should return ktop (phtbl cap to pe3t(ji,jj,1)) 
    215227      ! 
    216228      ! get ktbl 
  • NEMO/trunk/src/OCE/LDF/ldfslp.F90

    r13226 r13237  
    7575   !! * Substitutions 
    7676#  include "do_loop_substitute.h90" 
     77#  include "domzgr_substitute.h90" 
    7778   !!---------------------------------------------------------------------- 
    7879   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    198199         !                                      !              max slope = 1/2 * e3 / e1 
    199200         IF (ln_zps .AND. jk==mbku(ji,jj)) & 
    200             zbu = MIN(  zbu, - z1_slpmax * ABS( zau ) , - 2._wp * e1u(ji,jj) / e3u(ji,jj,jk,Kmm)* ABS( zau )  ) 
     201            zbu = MIN(  zbu, - z1_slpmax * ABS( zau ) ,   & 
     202               &                - 2._wp * e1u(ji,jj) / e3u(ji,jj,jk,Kmm)* ABS( zau )  ) 
    201203         IF (ln_zps .AND. jk==mbkv(ji,jj)) & 
    202             zbv = MIN(  zbv, - z1_slpmax * ABS( zav ) , - 2._wp * e2v(ji,jj) / e3v(ji,jj,jk,Kmm)* ABS( zav )  ) 
     204            zbv = MIN(  zbv, - z1_slpmax * ABS( zav ) ,   & 
     205               &                - 2._wp * e2v(ji,jj) / e3v(ji,jj,jk,Kmm)* ABS( zav )  ) 
    203206         !                                      ! uslp and vslp output in zwz and zww, resp. 
    204207         zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 
     
    206209         ! thickness of water column between surface and level k at u/v point 
    207210         zdepu = 0.5_wp * ( ( gdept (ji,jj,jk,Kmm) + gdept (ji+1,jj,jk,Kmm) )                            & 
    208                           - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj) ) - e3u(ji,jj,miku(ji,jj),Kmm)   ) 
     211            &              - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj) )        & 
     212            &              - e3u(ji,jj,miku(ji,jj),Kmm)   ) 
    209213         zdepv = 0.5_wp * ( ( gdept (ji,jj,jk,Kmm) + gdept (ji,jj+1,jk,Kmm) )                            & 
    210                           - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) ) - e3v(ji,jj,mikv(ji,jj),Kmm)   ) 
     214            &              - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) )        & 
     215            &              - e3v(ji,jj,mikv(ji,jj),Kmm)   ) 
    211216         ! 
    212217         zwz(ji,jj,jk) = ( ( 1._wp - zfi) * zau / ( zbu - zeps )                                     & 
     
    293298!               !                                         ! jk must be >= ML level for zfk=1. otherwise  zfk=0. 
    294299!               zfk = REAL( 1 - 1/(1 + jk / nmln(ji+1,jj)), wp ) 
    295 !               zck = gdepw(ji,jj,jk)    / MAX( hmlp(ji,jj), 10. ) 
     300!               zck = gdepw(ji,jj,jk,Kmm)    / MAX( hmlp(ji,jj), 10. ) 
    296301!               zwz(ji,jj,jk) = ( zfk * zai / ( zbi - zeps ) + ( 1._wp - zfk ) * wslpiml(ji,jj) * zck ) * tmask(ji,jj,jk) 
    297302!               zww(ji,jj,jk) = ( zfk * zaj / ( zbj - zeps ) + ( 1._wp - zfk ) * wslpjml(ji,jj) * zck ) * tmask(ji,jj,jk) 
  • NEMO/trunk/src/OCE/LDF/ldftra.F90

    r13226 r13237  
    9595   !! * Substitutions 
    9696#  include "do_loop_substitute.h90" 
     97#  include "domzgr_substitute.h90" 
    9798   !!---------------------------------------------------------------------- 
    9899   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/SBC/fldread.F90

    r13226 r13237  
    127127   !! * Substitutions 
    128128#  include "do_loop_substitute.h90" 
     129#  include "domzgr_substitute.h90" 
    129130   !!---------------------------------------------------------------------- 
    130131   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    632633               zdhalf(jk) = zdhalf(jk-1) + e3v(ji,jj,jk-1,Kmm) 
    633634               zdepth(jk) =          zcoef  * ( zdhalf(jk  ) + 0.5_wp * e3vw(ji,jj,jk,Kmm))  & 
    634                   &         + (1._wp-zcoef) * ( zdepth(jk-1) +          e3vw(ji,jj,jk,Kmm)) 
     635                     + (1._wp-zcoef) * ( zdepth(jk-1) +          e3vw(ji,jj,jk,Kmm)) 
    635636            END DO 
    636637         END SELECT 
  • NEMO/trunk/src/OCE/SBC/sbccpl.F90

    r13226 r13237  
    199199   !! Substitution 
    200200#  include "do_loop_substitute.h90" 
     201#  include "domzgr_substitute.h90" 
    201202   !!---------------------------------------------------------------------- 
    202203   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/SBC/sbcice_cice.F90

    r13226 r13237  
    1212   USE oce             ! ocean dynamics and tracers 
    1313   USE dom_oce         ! ocean space and time domain 
     14# if ! defined key_qco 
    1415   USE domvvl 
     16# else 
     17   USE domqco 
     18# endif 
    1519   USE phycst, only : rcp, rho0, r1_rho0, rhos, rhoi 
    1620   USE in_out_manager  ! I/O manager 
     
    233237!!gm This should be put elsewhere....   (same remark for limsbc) 
    234238!!gm especially here it is assumed zstar coordinate, but it can be ztilde.... 
     239#if defined key_qco 
     240            IF( .NOT.ln_linssh )   CALL dom_qco_zgr( Kbb, Kmm, Kaa )   ! interpolation scale factor, depth and water column 
     241#else 
    235242            IF( .NOT.ln_linssh ) THEN 
    236243               ! 
    237244               DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
    238                   e3t(:,:,jk,Kmm) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kmm)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
    239                   e3t(:,:,jk,Kbb) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kbb)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
     245                  e3t(:,:,jk,Kmm) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kmm)*r1_ht_0(:,:)*tmask(:,:,jk) ) 
     246                  e3t(:,:,jk,Kbb) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kbb)*r1_ht_0(:,:)*tmask(:,:,jk) ) 
    240247               ENDDO 
    241248               e3t(:,:,:,Krhs) = e3t(:,:,:,Kbb) 
     
    267274               END DO 
    268275            ENDIF 
     276#endif 
    269277         ENDIF 
    270278      ENDIF 
  • NEMO/trunk/src/OCE/SBC/sbcrnf.F90

    r12489 r13237  
    7272   !! * Substitutions 
    7373#  include "do_loop_substitute.h90" 
     74#  include "domzgr_substitute.h90" 
    7475   !!---------------------------------------------------------------------- 
    7576   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/SBC/sbcssm.F90

    r12377 r13237  
    3232   LOGICAL, SAVE ::   l_ssm_mean = .FALSE.   ! keep track of whether means have been read from restart file 
    3333    
     34#  include "domzgr_substitute.h90" 
    3435   !!---------------------------------------------------------------------- 
    3536   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/SBC/sbcwave.F90

    r13226 r13237  
    7373   !! * Substitutions 
    7474#  include "do_loop_substitute.h90" 
     75#  include "domzgr_substitute.h90" 
    7576   !!---------------------------------------------------------------------- 
    7677   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    207208            &                 - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * usd(ji-1,jj,jk)    & 
    208209            &                 + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vsd(ji,jj  ,jk)    & 
    209             &                 - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj) 
     210            &                 - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vsd(ji,jj-1,jk)  ) & 
     211            &                * r1_e1e2t(ji,jj) 
    210212      END_3D 
    211213      ! 
  • NEMO/trunk/src/OCE/TRA/eosbn2.F90

    r12489 r13237  
    180180   !! * Substitutions 
    181181#  include "do_loop_substitute.h90" 
     182#  include "domzgr_substitute.h90" 
    182183   !!---------------------------------------------------------------------- 
    183184   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/traadv.F90

    r12489 r13237  
    6666   INTEGER, PARAMETER ::   np_QCK     = 5   ! QUICK scheme 
    6767    
     68#  include "domzgr_substitute.h90" 
    6869   !!---------------------------------------------------------------------- 
    6970   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9899      IF( ln_wave .AND. ln_sdw )  THEN 
    99100         DO jk = 1, jpkm1                                                       ! eulerian transport + Stokes Drift 
    100             zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) ) 
    101             zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) ) 
    102             zww(:,:,jk) = e1e2t(:,:)                 * ( ww(:,:,jk) + wsd(:,:,jk) ) 
     101            zuu(:,:,jk) =   & 
     102               &  e2u  (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) ) 
     103            zvv(:,:,jk) =   &  
     104               &  e1v  (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) ) 
     105            zww(:,:,jk) =   &  
     106               &  e1e2t(:,:)                 * ( ww(:,:,jk) + wsd(:,:,jk) ) 
    103107         END DO 
    104108      ELSE 
  • NEMO/trunk/src/OCE/TRA/traadv_cen.F90

    r13226 r13237  
    3737   !! * Substitutions 
    3838#  include "do_loop_substitute.h90" 
     39#  include "domzgr_substitute.h90" 
    3940   !!---------------------------------------------------------------------- 
    4041   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    161162               &             - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )    & 
    162163               &                + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  )    & 
    163                &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     164               &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) & 
     165               &                * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    164166         END_3D 
    165167         !                             ! trend diagnostics 
  • NEMO/trunk/src/OCE/TRA/traadv_fct.F90

    r13226 r13237  
    4646   !! * Substitutions 
    4747#  include "do_loop_substitute.h90" 
     48#  include "domzgr_substitute.h90" 
    4849   !!---------------------------------------------------------------------- 
    4950   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    139140         ALLOCATE(zwdia(jpi,jpj,jpk), zwinf(jpi,jpj,jpk),zwsup(jpi,jpj,jpk)) 
    140141         DO_3D_00_00( 1, jpkm1 ) 
    141             zwdia(ji,jj,jk) =  1._wp + p2dt * ( MAX( wi(ji,jj,jk  ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) / e3t(ji,jj,jk,Krhs) 
     142            zwdia(ji,jj,jk) =  1._wp + p2dt * ( MAX( wi(ji,jj,jk) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) )   & 
     143            &                               / e3t(ji,jj,jk,Krhs) 
    142144            zwinf(ji,jj,jk) =  p2dt * MIN( wi(ji,jj,jk  ) , 0._wp ) / e3t(ji,jj,jk,Krhs) 
    143145            zwsup(ji,jj,jk) = -p2dt * MAX( wi(ji,jj,jk+1) , 0._wp ) / e3t(ji,jj,jk,Krhs) 
     
    180182               &      + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1) ) * r1_e1e2t(ji,jj) 
    181183            !                             ! update and guess with monotonic sheme 
    182             pt(ji,jj,jk,jn,Krhs) =                     pt(ji,jj,jk,jn,Krhs) +        ztra   / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 
    183             zwi(ji,jj,jk)    = ( e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * ztra ) / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk) 
     184            pt(ji,jj,jk,jn,Krhs) =                   pt(ji,jj,jk,jn,Krhs) +       ztra   & 
     185               &                                  / e3t(ji,jj,jk,Kmm ) * tmask(ji,jj,jk) 
     186            zwi(ji,jj,jk)    = ( e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * ztra ) & 
     187               &                                  / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk) 
    184188         END_3D 
    185189          
  • NEMO/trunk/src/OCE/TRA/traadv_mus.F90

    r13226 r13237  
    4747   !! * Substitutions 
    4848#  include "do_loop_substitute.h90" 
     49#  include "domzgr_substitute.h90" 
    4950   !!---------------------------------------------------------------------- 
    5051   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    227228         ! 
    228229         DO_3D_00_00( 1, jpkm1 ) 
    229             pt(ji,jj,jk,jn,Krhs) =  pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     230            pt(ji,jj,jk,jn,Krhs) =  pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) )   & 
     231               &                                      * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    230232         END_3D 
    231233         !                                ! send trends for diagnostic 
  • NEMO/trunk/src/OCE/TRA/traadv_qck.F90

    r13226 r13237  
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/traadv_ubs.F90

    r13226 r13237  
    3939   !! * Substitutions 
    4040#  include "do_loop_substitute.h90" 
     41#  include "domzgr_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    4243   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    8283      !! 
    8384      !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404.  
    84       !!             Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731Ð1741.  
     85      !!             Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731�1741. 
    8586      !!---------------------------------------------------------------------- 
    8687      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     
    158159               pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs)                        & 
    159160                  &             - (  ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk)    & 
    160                   &                + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     161                  &                + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) & 
     162                  &                * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    161163            END_2D 
    162164            !                                              
     
    202204            ! 
    203205            DO_3D_00_00( 1, jpkm1 ) 
    204                ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     206               ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) )    & 
     207                  &     * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    205208               pt(ji,jj,jk,jn,Krhs) =   pt(ji,jj,jk,jn,Krhs) +  ztak  
    206209               zti(ji,jj,jk)    = ( pt(ji,jj,jk,jn,Kbb) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk) 
     
    228231         ! 
    229232         DO_3D_00_00( 1, jpkm1 ) 
    230             pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     233            pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) )    & 
     234               &                                        * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    231235         END_3D 
    232236         ! 
  • NEMO/trunk/src/OCE/TRA/traatf.F90

    r13226 r13237  
    5858   !! * Substitutions 
    5959#  include "do_loop_substitute.h90" 
     60#  include "domzgr_substitute.h90" 
    6061   !!---------------------------------------------------------------------- 
    6162   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    162163      ! 
    163164      IF( l_trdtra .AND. ln_linssh ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt      
    164          zfact = 1._wp / rDt              
    165165         DO jk = 1, jpkm1 
    166             ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kmm) - ztrdt(:,:,jk) ) * zfact 
    167             ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kmm) - ztrds(:,:,jk) ) * zfact 
     166            ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kmm) - ztrdt(:,:,jk) ) * r1_Dt 
     167            ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kmm) - ztrds(:,:,jk) ) * r1_Dt 
    168168         END DO 
    169169         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_atf, ztrdt ) 
     
    229229      !!  
    230230      !! ** Method  : - Apply a thickness weighted Asselin time filter on now fields. 
    231       !!             pt(Kmm)  = ( e3t(Kmm)*pt(Kmm) + rn_atfp*[ e3t(Kbb)*pt(Kbb) - 2 e3t(Kmm)*pt(Kmm) + e3t_a*pt(Kaa) ] ) 
    232       !!                       /( e3t(Kmm)         + rn_atfp*[ e3t(Kbb)         - 2 e3t(Kmm)         + e3t(Kaa)    ] ) 
     231      !!             pt(Kmm)  = ( e3t_Kmm*pt(Kmm) + rn_atfp*[ e3t_Kbb*pt(Kbb) - 2 e3t_Kmm*pt(Kmm) + e3t_Kaa*pt(Kaa) ] ) 
     232      !!                       /( e3t_Kmm         + rn_atfp*[ e3t_Kbb         - 2 e3t_Kmm         + e3t_Kaa    ] ) 
    233233      !! 
    234234      !! ** Action  : - pt(Kmm) ready for the next time step 
  • NEMO/trunk/src/OCE/TRA/trabbc.F90

    r13226 r13237  
    4646   !! * Substitutions 
    4747#  include "do_loop_substitute.h90" 
     48#  include "domzgr_substitute.h90" 
    4849   !!---------------------------------------------------------------------- 
    4950   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9192      !                             !  Add the geothermal trend on temperature 
    9293      DO_2D_00_00 
    93          pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs) = pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs) + qgh_trd0(ji,jj) / e3t(ji,jj,mbkt(ji,jj),Kmm) 
     94         pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs) = pts(ji,jj,mbkt(ji,jj),jp_tem,Krhs)   & 
     95            &             + qgh_trd0(ji,jj) / e3t(ji,jj,mbkt(ji,jj),Kmm) 
    9496      END_2D 
    9597      ! 
  • NEMO/trunk/src/OCE/TRA/trabbl.F90

    r13226 r13237  
    6868   !! * Substitutions 
    6969#  include "do_loop_substitute.h90" 
     70#  include "domzgr_substitute.h90" 
    7071   !!---------------------------------------------------------------------- 
    7172   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/traisf.F90

    r12377 r13237  
    1111   !!---------------------------------------------------------------------- 
    1212   USE isf_oce                                     ! Ice shelf variables 
    13    USE dom_oce , ONLY : e3t, r1_e1e2t            ! ocean space domain variables 
     13   USE dom_oce                                     ! ocean space domain variables 
    1414   USE isfutils, ONLY : debug                      ! debug option 
    1515   USE timing  , ONLY : timing_start, timing_stop  ! Timing 
     
    2323   !! * Substitutions 
    2424#  include "do_loop_substitute.h90" 
     25#  include "domzgr_substitute.h90" 
    2526   !!---------------------------------------------------------------------- 
    2627   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    140141      ! 
    141142      DO jk = 1,jpk 
    142          ptsa(:,:,jk,jp_tem) = ptsa(:,:,jk,jp_tem) + ptsc(:,:,jk,jp_tem) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) 
    143          ptsa(:,:,jk,jp_sal) = ptsa(:,:,jk,jp_sal) + ptsc(:,:,jk,jp_sal) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) 
     143         ptsa(:,:,jk,jp_tem) =   & 
     144            &  ptsa(:,:,jk,jp_tem) + ptsc(:,:,jk,jp_tem) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) 
     145         ptsa(:,:,jk,jp_sal) =   & 
     146            &  ptsa(:,:,jk,jp_sal) + ptsc(:,:,jk,jp_sal) * r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) 
    144147      END DO 
    145148      ! 
  • NEMO/trunk/src/OCE/TRA/traldf_iso.F90

    r12489 r13237  
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    167168            IF( ln_traldf_blp ) THEN                ! bilaplacian operator 
    168169               DO_3D_10_10( 2, jpkm1 ) 
    169                   akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   & 
    170                      &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
     170                  akz(ji,jj,jk) = 16._wp   & 
     171                     &   * ah_wslp2   (ji,jj,jk)   & 
     172                     &   * (  akz     (ji,jj,jk)   & 
     173                     &      + ah_wslp2(ji,jj,jk)   & 
     174                     &        / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
    171175               END_3D 
    172176            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator 
     
    247251            ! 
    248252            DO_2D_00_00 
    249                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      & 
    250                   &                                                 + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
     253               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
     254                  &       + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk) + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
    251255                  &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    252256            END_2D 
     
    294298            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2 
    295299               DO_3D_00_00( 2, jpkm1 ) 
    296                   ztfw(ji,jj,jk) = ztfw(ji,jj,jk)                       & 
    297                      &           + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj)   & 
     300                  ztfw(ji,jj,jk) =   & 
     301                     &  ztfw(ji,jj,jk) + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj)   & 
    298302                     &           * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk) 
    299303               END_3D 
     
    308312         !          
    309313         DO_3D_00_00( 1, jpkm1 ) 
    310             pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   & 
    311                &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     314            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  ) * r1_e1e2t(ji,jj)   & 
     315               &                                             / e3t(ji,jj,jk,Kmm) 
    312316         END_3D 
    313317         ! 
  • NEMO/trunk/src/OCE/TRA/traldf_lap_blp.F90

    r13226 r13237  
    3838   !! * Substitutions 
    3939#  include "do_loop_substitute.h90" 
     40#  include "domzgr_substitute.h90" 
    4041   !!---------------------------------------------------------------------- 
    4142   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/traldf_triad.F90

    r12489 r13237  
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    179180            IF( ln_traldf_blp ) THEN                ! bilaplacian operator 
    180181               DO_3D_10_10( 2, jpkm1 ) 
    181                   akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   & 
    182                      &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
     182                  akz(ji,jj,jk) = 16._wp           & 
     183                     &   * ah_wslp2   (ji,jj,jk)   & 
     184                     &   * (  akz     (ji,jj,jk)   & 
     185                     &      + ah_wslp2(ji,jj,jk)   & 
     186                     &        / ( e3w (ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
    183187               END_3D 
    184188            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator 
     
    326330            !                             !==  horizontal divergence and add to the general trend  ==! 
    327331            DO_2D_00_00 
    328                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji-1,jj,jk) - zftu(ji,jj,jk)       & 
     332               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
     333                  &                       + zsign * (  zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)       & 
    329334                  &                                           + zftv(ji,jj-1,jk) - zftv(ji,jj,jk)   )   & 
    330335                  &                                        / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
     
    357362         ! 
    358363         DO_3D_00_00( 1, jpkm1 ) 
    359             pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk)  )   & 
     364            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
     365            &                                  + zsign * (  ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk)  )   & 
    360366               &                                              / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) ) 
    361367         END_3D 
  • NEMO/trunk/src/OCE/TRA/tramle.F90

    r13226 r13237  
    4949   !! * Substitutions 
    5050#  include "do_loop_substitute.h90" 
     51#  include "domzgr_substitute.h90" 
    5152   !!---------------------------------------------------------------------- 
    5253   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/tranpc.F90

    r13226 r13237  
    3535   !! * Substitutions 
    3636#  include "do_loop_substitute.h90" 
     37#  include "domzgr_substitute.h90" 
    3738   !!---------------------------------------------------------------------- 
    3839   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRA/traqsr.F90

    r13205 r13237  
    6868   !! * Substitutions 
    6969#  include "do_loop_substitute.h90" 
     70#  include "domzgr_substitute.h90" 
    7071   !!---------------------------------------------------------------------- 
    7172   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    266267      DO_3D_00_00( 1, nksr ) 
    267268         pts(ji,jj,jk,jp_tem,Krhs) = pts(ji,jj,jk,jp_tem,Krhs)   & 
    268             &                      + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) / e3t(ji,jj,jk,Kmm) 
     269            &                      + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) )   & 
     270            &                             / e3t(ji,jj,jk,Kmm) 
    269271      END_3D 
    270272      ! 
  • NEMO/trunk/src/OCE/TRA/trasbc.F90

    r12489 r13237  
    4343   !! * Substitutions 
    4444#  include "do_loop_substitute.h90" 
     45#  include "domzgr_substitute.h90" 
    4546   !!---------------------------------------------------------------------- 
    4647   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    138139      DO jn = 1, jpts               !==  update tracer trend  ==! 
    139140         DO_2D_01_00 
    140             pts(ji,jj,1,jn,Krhs) = pts(ji,jj,1,jn,Krhs) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) / e3t(ji,jj,1,Kmm) 
     141            pts(ji,jj,1,jn,Krhs) = pts(ji,jj,1,jn,Krhs) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) )    & 
     142               &                                                / e3t(ji,jj,1,Kmm) 
    141143         END_2D 
    142144      END DO 
  • NEMO/trunk/src/OCE/TRA/trazdf.F90

    r13226 r13237  
    3737   !! * Substitutions 
    3838#  include "do_loop_substitute.h90" 
     39#  include "domzgr_substitute.h90" 
    3940   !!---------------------------------------------------------------------- 
    4041   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    8485      IF( l_trdtra )   THEN                      ! save the vertical diffusive trends for further diagnostics 
    8586         DO jk = 1, jpkm1 
    86             ztrdt(:,:,jk) = ( ( pts(:,:,jk,jp_tem,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_tem,Kbb)*e3t(:,:,jk,Kbb) ) & 
    87                &          / (e3t(:,:,jk,Kmm)*rDt) ) - ztrdt(:,:,jk) 
    88             ztrds(:,:,jk) = ( ( pts(:,:,jk,jp_sal,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_sal,Kbb)*e3t(:,:,jk,Kbb) ) & 
    89               &           / (e3t(:,:,jk,Kmm)*rDt) ) - ztrds(:,:,jk) 
     87            ztrdt(:,:,jk) = (   (  pts(:,:,jk,jp_tem,Kaa)*e3t(:,:,jk,Kaa)     & 
     88               &                 - pts(:,:,jk,jp_tem,Kbb)*e3t(:,:,jk,Kbb)  )  & 
     89               &              / (  e3t(:,:,jk,Kmm)*rDt  )   )                 & 
     90               &          - ztrdt(:,:,jk) 
     91            ztrds(:,:,jk) = (   (  pts(:,:,jk,jp_sal,Kaa)*e3t(:,:,jk,Kaa)     & 
     92               &                 - pts(:,:,jk,jp_sal,Kbb)*e3t(:,:,jk,Kbb)  )  & 
     93               &             / (   e3t(:,:,jk,Kmm)*rDt  )   )                 & 
     94               &          - ztrds(:,:,jk) 
    9095         END DO 
    9196!!gm this should be moved in trdtra.F90 and done on all trends 
     
    213218         !          
    214219         DO_2D_00_00 
    215             pt(ji,jj,1,jn,Kaa) = e3t(ji,jj,1,Kbb) * pt(ji,jj,1,jn,Kbb) + p2dt * e3t(ji,jj,1,Kmm) * pt(ji,jj,1,jn,Krhs) 
     220            pt(ji,jj,1,jn,Kaa) =        e3t(ji,jj,1,Kbb) * pt(ji,jj,1,jn,Kbb)    & 
     221               &               + p2dt * e3t(ji,jj,1,Kmm) * pt(ji,jj,1,jn,Krhs) 
    216222         END_2D 
    217223         DO_3D_00_00( 2, jpkm1 ) 
    218             zrhs = e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * e3t(ji,jj,jk,Kmm) * pt(ji,jj,jk,jn,Krhs)   ! zrhs=right hand side 
     224            zrhs =        e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb)    &  
     225               & + p2dt * e3t(ji,jj,jk,Kmm) * pt(ji,jj,jk,jn,Krhs)   ! zrhs=right hand side 
    219226            pt(ji,jj,jk,jn,Kaa) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pt(ji,jj,jk-1,jn,Kaa) 
    220227         END_3D 
  • NEMO/trunk/src/OCE/TRA/zpshde.F90

    r13226 r13237  
    3232   !! * Substitutions 
    3333#  include "do_loop_substitute.h90" 
     34#  include "domzgr_substitute.h90" 
    3435   !!---------------------------------------------------------------------- 
    3536   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6566      !!              ___ |   |   |           ___  |   |   | 
    6667      !!                   
    67       !!      case 1->   e3w(i+1) >= e3w(i) ( and e3w(j+1) >= e3w(j) ) then 
    68       !!          t~ = t(i+1,j  ,k) + (e3w(i+1) - e3w(i)) * dk(Ti+1)/e3w(i+1) 
    69       !!        ( t~ = t(i  ,j+1,k) + (e3w(j+1) - e3w(j)) * dk(Tj+1)/e3w(j+1)  ) 
     68      !!      case 1->   e3w(i+1,:,:,Kmm) >= e3w(i,:,:,Kmm) ( and e3w(:,j+1,:,Kmm) >= e3w(:,j,:,Kmm) ) then 
     69      !!          t~ = t(i+1,j  ,k) + (e3w(i+1,j,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Ti+1)/e3w(i+1,j,k,Kmm) 
     70      !!        ( t~ = t(i  ,j+1,k) + (e3w(i,j+1,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Tj+1)/e3w(i,j+1,k,Kmm)  ) 
    7071      !!          or 
    71       !!      case 2->   e3w(i+1) <= e3w(i) ( and e3w(j+1) <= e3w(j) ) then 
    72       !!          t~ = t(i,j,k) + (e3w(i) - e3w(i+1)) * dk(Ti)/e3w(i ) 
    73       !!        ( t~ = t(i,j,k) + (e3w(j) - e3w(j+1)) * dk(Tj)/e3w(j ) ) 
     72      !!      case 2->   e3w(i+1,:,:,Kmm) <= e3w(i,:,:,Kmm) ( and e3w(:,j+1,:,Kmm) <= e3w(:,j,:,Kmm) ) then 
     73      !!          t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i+1,j,k,Kmm)) * dk(Ti)/e3w(i,j,k,Kmm) 
     74      !!        ( t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i,j+1,k,Kmm)) * dk(Tj)/e3w(i,j,k,Kmm) ) 
    7475      !!          Idem for di(s) and dj(s)           
    7576      !! 
     
    109110            iku = mbku(ji,jj)   ;   ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    110111            ikv = mbkv(ji,jj)   ;   ikvm1 = MAX( ikv - 1 , 1 )    ! if level first is a p-step, ik.m1=1 
    111 !!gm BUG ? when applied to before fields, e3w(:,:,:,Kbb) should be used.... 
     112!!gm BUG ? when applied to before fields, e3w(:,:,k,Kbb) should be used.... 
    112113            ze3wu = e3w(ji+1,jj  ,iku,Kmm) - e3w(ji,jj,iku,Kmm) 
    113114            ze3wv = e3w(ji  ,jj+1,ikv,Kmm) - e3w(ji,jj,ikv,Kmm) 
     
    214215      !!              ___ |   |   |           ___  |   |   | 
    215216      !!                   
    216       !!      case 1->   e3w(i+1) >= e3w(i) ( and e3w(j+1) >= e3w(j) ) then 
    217       !!          t~ = t(i+1,j  ,k) + (e3w(i+1) - e3w(i)) * dk(Ti+1)/e3w(i+1) 
    218       !!        ( t~ = t(i  ,j+1,k) + (e3w(j+1) - e3w(j)) * dk(Tj+1)/e3w(j+1)  ) 
     217      !!      case 1->   e3w(i+1,j,k,Kmm) >= e3w(i,j,k,Kmm) ( and e3w(i,j+1,k,Kmm) >= e3w(i,j,k,Kmm) ) then 
     218      !!          t~ = t(i+1,j  ,k) + (e3w(i+1,j  ,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Ti+1)/e3w(i+1,j  ,k,Kmm) 
     219      !!        ( t~ = t(i  ,j+1,k) + (e3w(i  ,j+1,k,Kmm) - e3w(i,j,k,Kmm)) * dk(Tj+1)/e3w(i  ,j+1,k,Kmm)  ) 
    219220      !!          or 
    220       !!      case 2->   e3w(i+1) <= e3w(i) ( and e3w(j+1) <= e3w(j) ) then 
    221       !!          t~ = t(i,j,k) + (e3w(i) - e3w(i+1)) * dk(Ti)/e3w(i ) 
    222       !!        ( t~ = t(i,j,k) + (e3w(j) - e3w(j+1)) * dk(Tj)/e3w(j ) ) 
     221      !!      case 2->   e3w(i+1,j,k,Kmm) <= e3w(i,j,k,Kmm) ( and e3w(i,j+1,k,Kmm) <= e3w(i,j,k,Kmm) ) then 
     222      !!          t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i+1,j  ,k,Kmm)) * dk(Ti)/e3w(i,j,k,Kmm) 
     223      !!        ( t~ = t(i,j,k) + (e3w(i,j,k,Kmm) - e3w(i  ,j+1,k,Kmm)) * dk(Tj)/e3w(i,j,k,Kmm) ) 
    223224      !!          Idem for di(s) and dj(s)           
    224225      !! 
     
    356357            ! (ISF) case partial step top and bottom in adjacent cell in vertical 
    357358            ! cannot used e3w because if 2 cell water column, we have ps at top and bottom 
    358             ! in this case e3w(i,j) - e3w(i,j+1) is not the distance between Tj~ and Tj 
     359            ! in this case e3w(i,j,k,Kmm) - e3w(i,j+1,k,Kmm) is not the distance between Tj~ and Tj 
    359360            ! the only common depth between cells (i,j) and (i,j+1) is gdepw_0 
    360361            ze3wu  =  gdept(ji,jj,iku,Kmm) - gdept(ji+1,jj,iku,Kmm) 
  • NEMO/trunk/src/OCE/TRD/trddyn.F90

    r13226 r13237  
    3737   !! * Substitutions 
    3838#  include "do_loop_substitute.h90" 
     39#  include "domzgr_substitute.h90" 
    3940   !!---------------------------------------------------------------------- 
    4041   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRD/trdglo.F90

    r12489 r13237  
    5252   !! * Substitutions 
    5353#  include "do_loop_substitute.h90" 
     54#  include "domzgr_substitute.h90" 
    5455   !!---------------------------------------------------------------------- 
    5556   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    211212         zcof   = 0.5_wp / rho0           ! Density flux at u and v-points 
    212213         DO_3D_10_10( 1, jpkm1 ) 
    213             zkx(ji,jj,jk) = zcof * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * uu(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) ) 
    214             zky(ji,jj,jk) = zcof * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) * vv(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) ) 
     214            zkx(ji,jj,jk) = zcof * e2u(ji,jj) * e3u(ji,jj,jk,Kmm)   & 
     215               &                              *  uu(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) ) 
     216            zky(ji,jj,jk) = zcof * e1v(ji,jj) * e3v(ji,jj,jk,Kmm)   & 
     217               &                              *  vv(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) ) 
    215218         END_3D 
    216219          
     
    226229         peke = 0._wp 
    227230         DO jk = 1, jpkm1 
    228             peke = peke + SUM( zkepe(:,:,jk) * gdept(:,:,jk,Kmm) * e1e2t(:,:) * e3t(:,:,jk,Kmm) ) 
     231            peke = peke + SUM( zkepe(:,:,jk) * gdept(:,:,jk,Kmm) * e1e2t(:,:)   & 
     232               &                               * e3t(:,:,jk,Kmm) ) 
    229233         END DO 
    230234         peke = grav * peke 
     
    524528 
    525529      DO_3D_00_00( 1, jpk ) 
    526          tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk) 
    527          tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) 
     530         tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm)   & 
     531            &                                       * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk) 
     532         tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * e3v(ji,jj,jk,Kmm)   & 
     533            &                                       * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) 
    528534      END_3D 
    529535      CALL mpp_sum( 'trdglo', tvolu )   ! sums over the global domain 
  • NEMO/trunk/src/OCE/TRD/trdken.F90

    r13226 r13237  
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/TRD/trdmxl.F90

    r13226 r13237  
    7070   !! * Substitutions 
    7171#  include "do_loop_substitute.h90" 
     72#  include "domzgr_substitute.h90" 
    7273   !!---------------------------------------------------------------------- 
    7374   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    120121         wkx(:,:,:) = 0._wp         !==  now ML weights for vertical averaging  ==! 
    121122         DO_3D_11_11( 1, jpktrd ) 
    122             IF( jk - kmxln(ji,jj) < 0 )   wkx(ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 
     123            IF( jk - kmxln(ji,jj) < 0 )   THEN 
     124               wkx(ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 
     125            ENDIF 
    123126         END_3D 
    124127         hmxl(:,:) = 0._wp               ! NOW mixed-layer depth 
  • NEMO/trunk/src/OCE/TRD/trdpen.F90

    r12377 r13237  
    3535   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   rab_pe   ! partial derivatives of PE anomaly with respect to T and S 
    3636 
     37   !! * Substitutions 
     38#  include "domzgr_substitute.h90" 
    3739   !!---------------------------------------------------------------------- 
    3840   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    4042   !! Software governed by the CeCILL license (see ./LICENSE) 
    4143   !!---------------------------------------------------------------------- 
     44 
    4245CONTAINS 
    4346 
  • NEMO/trunk/src/OCE/TRD/trdtra.F90

    r13084 r13237  
    4242   !! * Substitutions 
    4343#  include "do_loop_substitute.h90" 
     44#  include "domzgr_substitute.h90" 
    4445   !!---------------------------------------------------------------------- 
    4546   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    131132            zwt(:,:,jpk) = 0._wp   ;   zws(:,:,jpk) = 0._wp 
    132133            DO jk = 2, jpk 
    133                zwt(:,:,jk) = avt(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
    134                zws(:,:,jk) = avs(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
     134               zwt(:,:,jk) = avt(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) )   & 
     135                  &        / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
     136               zws(:,:,jk) = avs(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) )   & 
     137                  &        / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
    135138            END DO 
    136139            ! 
     
    145148            zwt(:,:,:) = 0._wp   ;   zws(:,:,:) = 0._wp            ! vertical diffusive fluxes 
    146149            DO jk = 2, jpk 
    147                zwt(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
    148                zws(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
     150               zwt(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) )   & 
     151                  &            / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
     152               zws(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) )   & 
     153                  &            / e3w(:,:,jk,Kmm) * tmask(:,:,jk) 
    149154            END DO 
    150155            ! 
  • NEMO/trunk/src/OCE/TRD/trdvor.F90

    r13226 r13237  
    5757   !! * Substitutions 
    5858#  include "do_loop_substitute.h90" 
     59#  include "domzgr_substitute.h90" 
    5960   !!---------------------------------------------------------------------- 
    6061   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    192193         DO jj = 1, jpjm1 
    193194            vortrd(ji,jj,ktrd) = (    zvdpvor(ji+1,jj) - zvdpvor(ji,jj)       & 
    194                  &                - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) )   ) / ( e1f(ji,jj) * e2f(ji,jj) ) 
     195                 &                - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) )   ) & 
     196                 &                  / ( e1f(ji,jj) * e2f(ji,jj) ) 
    195197         END DO 
    196198      END DO 
     
    268270            DO jj = 1, jpjm1 
    269271               vortrd(ji,jj,jpvor_bev) = (    zvbet(ji+1,jj) - zvbet(ji,jj)     & 
    270                   &                       - ( zubet(ji,jj+1) - zubet(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) 
     272                  &                       - ( zubet(ji,jj+1) - zubet(ji,jj) ) ) & 
     273                  &                           / ( e1f(ji,jj) * e2f(ji,jj) ) 
    271274            END DO 
    272275         END DO 
     
    283286         DO jj=1,jpjm1 
    284287            vortrd(ji,jj,ktrd) = (    zvdpvor(ji+1,jj) - zvdpvor(ji,jj)     & 
    285                &                  - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) 
     288               &                  - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) & 
     289               &                         / ( e1f(ji,jj) * e2f(ji,jj) ) 
    286290         END DO 
    287291      END DO 
     
    345349         DO jj = 1, jpjm1 
    346350            vor_avr(ji,jj) = (  ( zvv(ji+1,jj) - zvv(ji,jj) )    & 
    347                &              - ( zuu(ji,jj+1) - zuu(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1) 
     351               &              - ( zuu(ji,jj+1) - zuu(ji,jj) ) )  & 
     352               &             / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1) 
    348353         END DO 
    349354      END DO 
  • NEMO/trunk/src/OCE/ZDF/zdfddm.F90

    r13226 r13237  
    9696         DO_2D_11_11 
    9797            zrw =   ( gdepw(ji,jj,jk  ,Kmm) - gdept(ji,jj,jk,Kmm) )   & 
    98 !!gm please, use e3w(:,:,:,Kmm) below  
     98!!gm please, use e3w at Kmm below  
    9999               &  / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) )  
    100100            ! 
  • NEMO/trunk/src/OCE/ZDF/zdfdrg.F90

    r12489 r13237  
    7474   !! * Substitutions 
    7575#  include "do_loop_substitute.h90" 
     76#  include "domzgr_substitute.h90" 
    7677   !!---------------------------------------------------------------------- 
    7778   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/OCE/ZDF/zdfgls.F90

    r12489 r13237  
    105105   !! * Substitutions 
    106106#  include "do_loop_substitute.h90" 
     107#  include "domzgr_substitute.h90" 
    107108   !!---------------------------------------------------------------------- 
    108109   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    263264         zcof = rfact_tke * tmask(ji,jj,jk) 
    264265         !                                        ! lower diagonal, in fact not used for jk = 2 (see surface conditions) 
    265          zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t(ji,jj,jk-1,Kmm) * e3w(ji,jj,jk,Kmm) ) 
     266         zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) )   & 
     267            &                 / ( e3t(ji,jj,jk-1,Kmm) * e3w(ji,jj,jk,Kmm) ) 
    266268         !                                        ! upper diagonal, in fact not used for jk = ibotm1 (see bottom conditions) 
    267          zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t(ji,jj,jk  ,Kmm) * e3w(ji,jj,jk,Kmm) ) 
     269         zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) )   & 
     270            &                 / ( e3t(ji,jj,jk  ,Kmm) * e3w(ji,jj,jk,Kmm) ) 
    268271         !                                        ! diagonal 
    269272         zdiag(ji,jj,jk) = 1._wp - zd_lw(ji,jj,jk) - zd_up(ji,jj,jk)  + rn_Dt * zdiss * wmask(ji,jj,jk)  
     
    473476         zcof = rfact_psi * zwall_psi(ji,jj,jk) * tmask(ji,jj,jk) 
    474477         !                                               ! lower diagonal 
    475          zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) ) / ( e3t(ji,jj,jk-1,Kmm) * e3w(ji,jj,jk,Kmm) ) 
     478         zd_lw(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) )   & 
     479            &            / ( e3t(ji,jj,jk-1,Kmm) * e3w(ji,jj,jk,Kmm) ) 
    476480         !                                               ! upper diagonal 
    477          zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) ) / ( e3t(ji,jj,jk  ,Kmm) * e3w(ji,jj,jk,Kmm) ) 
     481         zd_up(ji,jj,jk) = zcof * ( p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) )   & 
     482            &            / ( e3t(ji,jj,jk  ,Kmm) * e3w(ji,jj,jk,Kmm) ) 
    478483         !                                               ! diagonal 
    479484         zdiag(ji,jj,jk) = 1._wp - zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) + rn_Dt * zdiss * wmask(ji,jj,jk) 
     
    11001105   !!====================================================================== 
    11011106END MODULE zdfgls 
    1102  
  • NEMO/trunk/src/OCE/ZDF/zdfiwm.F90

    r12510 r13237  
    5151   !! * Substitutions 
    5252#  include "do_loop_substitute.h90" 
     53#  include "domzgr_substitute.h90" 
    5354   !!---------------------------------------------------------------------- 
    5455   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9495      !!                 2. Pycnocline-intensified low-mode dissipation 
    9596      !!                     zemx_iwm(z) = ( epyc_iwm / rho0 ) * ( sqrt(rn2(z))^nn_zpyc ) 
    96       !!                                   / SUM( sqrt(rn2(z))^nn_zpyc * e3w(z) ) 
     97      !!                                   / SUM( sqrt(rn2(z))^nn_zpyc * e3w[z) ) 
    9798      !!              where epyc_iwm is a map of available power, and nn_zpyc 
    9899      !!              is the chosen stratification-dependence of the internal wave 
     
    100101      !!                 3. WKB-height dependent high mode dissipation 
    101102      !!                     zemx_iwm(z) = ( ebot_iwm / rho0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 
    102       !!                                   / SUM( rn2(z) * EXP(-z_wkb(z)/hbot_iwm) * e3w(z) ) 
     103      !!                                   / SUM( rn2(z) * EXP(-z_wkb(z)/hbot_iwm) * e3w[z) ) 
    103104      !!              where hbot_iwm is the characteristic length scale of the WKB bottom  
    104105      !!              intensification, ebot_iwm is a map of available power, and z_wkb is the 
    105106      !!              WKB-stretched height above bottom defined as 
    106       !!                    z_wkb(z) = H * SUM( sqrt(rn2(z'>=z)) * e3w(z'>=z) ) 
    107       !!                                 / SUM( sqrt(rn2(z'))    * e3w(z')    ) 
     107      !!                    z_wkb(z) = H * SUM( sqrt(rn2(z'>=z)) * e3w[z'>=z) ) 
     108      !!                                 / SUM( sqrt(rn2(z'))    * e3w[z')    ) 
    108109      !! 
    109110      !!              - update the model vertical eddy viscosity and diffusivity:  
     
    178179         zfact(:,:) = 0._wp 
    179180         DO jk = 2, jpkm1              ! part independent of the level 
    180             zfact(:,:) = zfact(:,:) + e3w(:,:,jk,Kmm) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
     181            zfact(:,:) =   & 
     182               &  zfact(:,:) +   & 
     183               &  e3w(:,:,jk,Kmm) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
    181184         END DO 
    182185         ! 
  • NEMO/trunk/src/OCE/ZDF/zdfmxl.F90

    r12489 r13237  
    3838   !! * Substitutions 
    3939#  include "do_loop_substitute.h90" 
     40#  include "domzgr_substitute.h90" 
    4041   !!---------------------------------------------------------------------- 
    4142   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    100101      DO_3D_11_11( nlb10, jpkm1 ) 
    101102         ikt = mbkt(ji,jj) 
    102          hmlp(ji,jj) = hmlp(ji,jj) + MAX( rn2b(ji,jj,jk) , 0._wp ) * e3w(ji,jj,jk,Kmm) 
     103         hmlp(ji,jj) =   & 
     104            & hmlp(ji,jj) + MAX( rn2b(ji,jj,jk) , 0._wp ) * e3w(ji,jj,jk,Kmm) 
    103105         IF( hmlp(ji,jj) < zN2_c )   nmln(ji,jj) = MIN( jk , ikt ) + 1   ! Mixed layer level 
    104106      END_3D 
  • NEMO/trunk/src/OCE/ZDF/zdfosm.F90

    r13226 r13237  
    103103   INTEGER :: idebug = 236 
    104104   INTEGER :: jdebug = 228 
     105    
    105106   !! * Substitutions 
    106107#  include "do_loop_substitute.h90" 
     108#  include "domzgr_substitute.h90" 
    107109   !!---------------------------------------------------------------------- 
    108110   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    503505                       & - zbeta * ( zs_bl(ji,jj) - ts(ji,jj,jm,jp_sal,Kmm) ) ), 0.0 ) + zvel_max 
    504506 
    505                   zhbl_s = zhbl_s + MIN( - zwb_ent(ji,jj) / zdb * rn_Dt / FLOAT(ibld(ji,jj)-imld(ji,jj) ), e3w(ji,jj,jk,Kmm) ) 
     507                  zhbl_s = zhbl_s + MIN( - zwb_ent(ji,jj) / zdb * rn_Dt / FLOAT(ibld(ji,jj)-imld(ji,jj) ),   & 
     508                     &                     e3w(ji,jj,jk,Kmm) ) 
    506509                  zhbl_s = MIN(zhbl_s, ht(ji,jj)) 
    507510 
     
    594597                     zwb_ent(ji,jj) = 0._wp 
    595598                  ENDIF 
    596                   inhml = MAX( INT( zari * zhbl(ji,jj) / e3t(ji,jj,ibld(ji,jj),Kmm) ) , 1 ) 
     599                  inhml = MAX( INT( zari * zhbl(ji,jj)   & 
     600                     &              / e3t(ji,jj,ibld(ji,jj),Kmm) ), 1 ) 
    597601                  imld(ji,jj) = MAX( ibld(ji,jj) - inhml, 1) 
    598602                  zhml(ji,jj) = gdepw(ji,jj,imld(ji,jj),Kmm) 
     
    610614                     zari = MIN( 4.5 * ( zvstr(ji,jj)**2 ) & 
    611615                       & / ( zdb_bl(ji,jj) * zhbl(ji,jj) ) + 0.01  , 0.2 ) 
    612                      inhml = MAX( INT( zari * zhbl(ji,jj) / e3t(ji,jj,ibld(ji,jj),Kmm) ) , 1 ) 
     616                     inhml = MAX( INT( zari * zhbl(ji,jj)   & 
     617                        &             / e3t(ji,jj,ibld(ji,jj),Kmm) ), 1 ) 
    613618                     imld(ji,jj) = MAX( ibld(ji,jj) - inhml, 1) 
    614619                     zhml(ji,jj) = gdepw(ji,jj,imld(ji,jj),Kmm) 
  • NEMO/trunk/src/OCE/ZDF/zdfsh2.F90

    r12377 r13237  
    2424   !! * Substitutions 
    2525#  include "do_loop_substitute.h90" 
     26#  include "domzgr_substitute.h90" 
    2627   !!---------------------------------------------------------------------- 
    2728   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6263            zsh2u(ji,jj) = ( p_avm(ji+1,jj,jk) + p_avm(ji,jj,jk) ) & 
    6364               &         * (   uu(ji,jj,jk-1,Kmm) -   uu(ji,jj,jk,Kmm) ) & 
    64                &         * (   uu(ji,jj,jk-1,Kbb) -   uu(ji,jj,jk,Kbb) ) / ( e3uw(ji,jj,jk,Kmm) * e3uw(ji,jj,jk,Kbb) ) * wumask(ji,jj,jk) 
     65               &         * (   uu(ji,jj,jk-1,Kbb) -   uu(ji,jj,jk,Kbb) ) &  
     66               &         / ( e3uw(ji,jj,jk  ,Kmm) * e3uw(ji,jj,jk,Kbb) ) & 
     67               &         * wumask(ji,jj,jk) 
    6568            zsh2v(ji,jj) = ( p_avm(ji,jj+1,jk) + p_avm(ji,jj,jk) ) & 
    6669               &         * (   vv(ji,jj,jk-1,Kmm) -   vv(ji,jj,jk,Kmm) ) & 
    67                &         * (   vv(ji,jj,jk-1,Kbb) -   vv(ji,jj,jk,Kbb) ) / ( e3vw(ji,jj,jk,Kmm) * e3vw(ji,jj,jk,Kbb) ) * wvmask(ji,jj,jk) 
     70               &         * (   vv(ji,jj,jk-1,Kbb) -   vv(ji,jj,jk,Kbb) ) & 
     71               &         / ( e3vw(ji,jj,jk  ,Kmm) * e3vw(ji,jj,jk,Kbb) ) & 
     72               &         * wvmask(ji,jj,jk) 
    6873         END_2D 
    6974         DO_2D_00_00 
  • NEMO/trunk/src/OCE/ZDF/zdftke.F90

    r13226 r13237  
    9898   !! * Substitutions 
    9999#  include "do_loop_substitute.h90" 
     100#  include "domzgr_substitute.h90" 
    100101   !!---------------------------------------------------------------------- 
    101102   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    267268         zpelc(:,:,1) =  MAX( rn2b(:,:,1), 0._wp ) * gdepw(:,:,1,Kmm) * e3w(:,:,1,Kmm) 
    268269         DO jk = 2, jpk 
    269             zpelc(:,:,jk)  = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0._wp ) * gdepw(:,:,jk,Kmm) * e3w(:,:,jk,Kmm) 
     270            zpelc(:,:,jk)  = zpelc(:,:,jk-1) +   & 
     271               &        MAX( rn2b(:,:,jk), 0._wp ) * gdepw(:,:,jk,Kmm) * e3w(:,:,jk,Kmm) 
    270272         END DO 
    271273         !                        !* finite Langmuir Circulation depth 
     
    325327         !                                   ! eddy coefficient (ensure numerical stability) 
    326328         zzd_up = zcof * MAX(  p_avm(ji,jj,jk+1) + p_avm(ji,jj,jk  ) , 2.e-5_wp  )   &  ! upper diagonal 
    327             &          /    (  e3t(ji,jj,jk  ,Kmm) * e3w(ji,jj,jk  ,Kmm)  ) 
     329            &          /    (  e3t(ji,jj,jk  ,Kmm)   & 
     330            &                * e3w(ji,jj,jk  ,Kmm)  ) 
    328331         zzd_lw = zcof * MAX(  p_avm(ji,jj,jk  ) + p_avm(ji,jj,jk-1) , 2.e-5_wp  )   &  ! lower diagonal 
    329             &          /    (  e3t(ji,jj,jk-1,Kmm) * e3w(ji,jj,jk  ,Kmm)  ) 
     332            &          /    (  e3t(ji,jj,jk-1,Kmm)   & 
     333            &                * e3w(ji,jj,jk  ,Kmm)  ) 
    330334         ! 
    331335         zd_up(ji,jj,jk) = zzd_up            ! Matrix (zdiag, zd_up, zd_lw) 
     
    515519            &            gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) - gdepw(ji,jj,jk,Kmm) ) 
    516520            ! wmask prevent zmxlm = 0 if jk = mikt(ji,jj) 
    517             zmxlm(ji,jj,jk) = zemxl * wmask(ji,jj,jk) + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 
    518             zmxld(ji,jj,jk) = zemxl * wmask(ji,jj,jk) + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 
     521            zmxlm(ji,jj,jk) = zemxl * wmask(ji,jj,jk)   & 
     522               &            + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 
     523            zmxld(ji,jj,jk) = zemxl * wmask(ji,jj,jk)   & 
     524               &            + MIN( zmxlm(ji,jj,jk) , e3w(ji,jj,jk,Kmm) ) * (1 - wmask(ji,jj,jk)) 
    519525         END_3D 
    520526         ! 
     
    528534      CASE ( 2 )           ! |dk[xml]| bounded by e3t : 
    529535         DO_3D_00_00( 2, jpkm1 ) 
    530             zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 
     536            zmxlm(ji,jj,jk) =   & 
     537               &    MIN( zmxlm(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 
    531538         END_3D 
    532539         DO_3DS_00_00( jpkm1, 2, -1 ) 
     
    538545      CASE ( 3 )           ! lup and ldown, |dk[xml]| bounded by e3t : 
    539546         DO_3D_00_00( 2, jpkm1 ) 
    540             zmxld(ji,jj,jk) = MIN( zmxld(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 
     547            zmxld(ji,jj,jk) =    & 
     548               &    MIN( zmxld(ji,jj,jk-1) + e3t(ji,jj,jk-1,Kmm), zmxlm(ji,jj,jk) ) 
    541549         END_3D 
    542550         DO_3DS_00_00( jpkm1, 2, -1 ) 
    543             zmxlm(ji,jj,jk) = MIN( zmxlm(ji,jj,jk+1) + e3t(ji,jj,jk+1,Kmm), zmxlm(ji,jj,jk) ) 
     551            zmxlm(ji,jj,jk) =   & 
     552               &    MIN( zmxlm(ji,jj,jk+1) + e3t(ji,jj,jk+1,Kmm), zmxlm(ji,jj,jk) ) 
    544553         END_3D 
    545554         DO_3D_00_00( 2, jpkm1 ) 
  • NEMO/trunk/src/OCE/nemogcm.F90

    r13226 r13237  
    6060   USE diacfl         ! CFL diagnostics               (dia_cfl_init routine) 
    6161   USE diamlr         ! IOM context management for multiple-linear-regression analysis 
     62#if defined key_qco 
     63   USE stepMLF        ! NEMO time-stepping               (stp_MLF   routine) 
     64#else 
    6265   USE step           ! NEMO time-stepping                 (stp     routine) 
     66#endif 
    6367   USE isfstp         ! ice shelf                     (isf_stp_init routine) 
    6468   USE icbini         ! handle bergs, initialisation 
     
    178182      ! 
    179183      DO WHILE( istp <= nitend .AND. nstop == 0 ) 
     184#if defined key_qco 
     185         CALL stp_MLF 
     186#else 
    180187         CALL stp 
     188#endif 
    181189         istp = istp + 1 
    182190      END DO 
     
    195203            ENDIF 
    196204             
     205#if defined key_qco 
     206            CALL stp_MLF      ( istp ) 
     207#else 
    197208            CALL stp        ( istp )  
     209#endif 
    198210            istp = istp + 1 
    199211 
     
    426438#endif 
    427439                           CALL     dom_init( Nbb, Nnn, Naa, "OPA") ! Domain 
    428  
    429  
    430  
    431440      IF( ln_crs       )   CALL     crs_init(      Nnn )       ! coarsened grid: domain initialization  
    432441      IF( sn_cfctl%l_prtctl )   & 
     
    704713   !!====================================================================== 
    705714END MODULE nemogcm 
    706  
  • NEMO/trunk/src/OCE/oce.F90

    r12489 r13237  
    3232   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   Cu_adv                   !: vertical Courant number (adaptive-implicit) 
    3333 
    34    !! free surface                                      !  before  ! now    ! after  ! 
    35    !! ------------                                      !  fields  ! fields ! fields ! 
     34   !! free surface 
     35   !! ------------ 
    3636   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ssh, uu_b,  vv_b   !: SSH [m] and barotropic velocities [m/s] 
    3737 
  • NEMO/trunk/src/OCE/step.F90

    r12933 r13237  
    3333   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rewrite in preparation for new timestepping scheme 
    3434   !!---------------------------------------------------------------------- 
    35  
     35#if defined key_qco 
     36   !!---------------------------------------------------------------------- 
     37   !!   'key_qco'      EMPTY MODULE      Quasi-Eulerian vertical coordonate 
     38   !!---------------------------------------------------------------------- 
     39#else 
    3640   !!---------------------------------------------------------------------- 
    3741   !!   stp             : OPA system time-stepping 
     
    179183                            CALL ssh_nxt       ( kstp, Nbb, Nnn, ssh, Naa )    ! after ssh (includes call to div_hor) 
    180184      IF( .NOT.ln_linssh )  CALL dom_vvl_sf_nxt( kstp, Nbb, Nnn,      Naa )    ! after vertical scale factors  
    181                             CALL wzv           ( kstp, Nbb, Nnn, ww,  Naa )    ! now cross-level velocity  
     185                            CALL wzv           ( kstp, Nbb, Nnn, Naa, ww )    ! now cross-level velocity  
    182186      IF( ln_zad_Aimp )     CALL wAimp         ( kstp,      Nnn           )  ! Adaptive-implicit vertical advection partitioning 
    183187                            CALL eos    ( ts(:,:,:,:,Nnn), rhd, rhop, gdept(:,:,:,Nnn) )  ! now in situ density for hpg computation 
     
    208212                            CALL dyn_zdf    ( kstp, Nbb, Nnn, Nrhs, uu, vv, Naa  )  ! vertical diffusion 
    209213      IF( ln_dynspg_ts ) THEN                                                       ! vertical scale factors and vertical velocity need to be updated 
    210                             CALL wzv        ( kstp, Nbb, Nnn, ww, Naa )             ! now cross-level velocity  
     214                            CALL wzv        ( kstp, Nbb, Nnn, Naa, ww )             ! now cross-level velocity  
    211215         IF( ln_zad_Aimp )  CALL wAimp      ( kstp,      Nnn )                      ! Adaptive-implicit vertical advection partitioning 
    212216      ENDIF 
     
    364368   END SUBROUTINE stp 
    365369   ! 
     370#endif 
    366371   !!====================================================================== 
    367372END MODULE step 
  • NEMO/trunk/src/OFF/dtadyn.F90

    r12489 r13237  
    2323   USE c1d             ! 1D configuration: lk_c1d 
    2424   USE dom_oce         ! ocean domain: variables 
     25#if ! defined key_qco  
    2526   USE domvvl          ! variable volume 
     27#else 
     28   USE domqco 
     29#endif 
    2630   USE zdf_oce         ! ocean vertical physics: variables 
    2731   USE sbc_oce         ! surface module: variables 
     
    5256   PUBLIC   dta_dyn_sed        ! called by nemo_gcm 
    5357   PUBLIC   dta_dyn_atf        ! called by nemo_gcm 
     58#if ! defined key_qco 
    5459   PUBLIC   dta_dyn_sf_interp  ! called by nemo_gcm 
     60#endif 
    5561 
    5662   CHARACTER(len=100) ::   cn_dir          !: Root directory for location of ssr files 
     
    149155         emp_b  (:,:)   = sf_dyn(jf_empb)%fnow(:,:,1) * tmask(:,:,1)    ! E-P 
    150156         zemp   (:,:)   = ( 0.5_wp * ( emp(:,:) + emp_b(:,:) ) + rnf(:,:) + fwbcorr ) * tmask(:,:,1) 
    151          CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa), e3t(:,:,:,Kaa) )  !=  ssh, vertical scale factor & vertical transport 
     157#if defined key_qco 
     158         CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa) ) 
     159         CALL dom_qco_r3c( ssh(:,:,Kaa), r3t(:,:,Kaa), r3u(:,:,Kaa), r3v(:,:,Kaa) ) 
     160#else 
     161         CALL dta_dyn_ssh( kt, zhdivtr, ssh(:,:,Kbb), zemp, ssh(:,:,Kaa), e3t(:,:,:,Kaa) )  !=  ssh, vertical scale factor 
     162#endif 
    152163         DEALLOCATE( zemp , zhdivtr ) 
    153164         !                                           Write in the tracer restart file 
     
    329340        ENDIF 
    330341        ! 
     342#if defined key_qco 
     343        CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb) ) 
     344        CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm) ) 
     345#else 
    331346        DO jk = 1, jpkm1 
    332            e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) ) 
     347           e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * r1_ht_0(:,:) * tmask(:,:,jk) ) 
    333348        ENDDO 
    334349        e3t(:,:,jpk,Kaa) = e3t_0(:,:,jpk) 
     
    342357        ! ------------------------------------ 
    343358        CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w(:,:,:,Kmm), 'W' ) 
    344    
     359!!gm this should be computed from ssh(Kbb)   
    345360        e3t(:,:,:,Kbb)  = e3t(:,:,:,Kmm) 
    346361        e3u(:,:,:,Kbb)  = e3u(:,:,:,Kmm) 
     
    367382        ! 
    368383      ENDIF 
     384#endif 
    369385      ! 
    370386      IF( ln_dynrnf .AND. ln_dynrnf_depth ) THEN       ! read depht over which runoffs are distributed 
     
    389405            ENDIF 
    390406         END_2D 
     407!!st pourquoi on n'utilise pas le gde3w ici plutôt que de faire une boucle ?  
    391408         DO_2D_11_11 
    392409            h_rnf(ji,jj) = 0._wp 
     
    413430   END SUBROUTINE dta_dyn_init 
    414431 
     432    
    415433   SUBROUTINE dta_dyn_sed( kt, Kmm ) 
    416434      !!---------------------------------------------------------------------- 
     
    529547   END SUBROUTINE dta_dyn_sed_init 
    530548 
     549    
    531550   SUBROUTINE dta_dyn_atf( kt, Kbb, Kmm, Kaa ) 
    532551     !!--------------------------------------------------------------------- 
     
    552571   END SUBROUTINE dta_dyn_atf 
    553572    
     573    
     574#if ! defined key_qco     
    554575   SUBROUTINE dta_dyn_sf_interp( kt, Kmm ) 
    555576      !!--------------------------------------------------------------------- 
     
    588609      ! 
    589610   END SUBROUTINE dta_dyn_sf_interp 
     611#endif 
     612 
    590613 
    591614   SUBROUTINE dta_dyn_ssh( kt, phdivtr, psshb,  pemp, pssha, pe3ta ) 
     
    606629      !!          The boundary conditions are w=0 at the bottom (no flux) 
    607630      !! 
    608       !! ** action  :   ssh(:,:,Kaa) / e3t(:,:,:,Kaa) / ww 
     631      !! ** action  :   ssh(:,:,Kaa) / e3t(:,:,k,Kaa) / ww 
    609632      !! 
    610633      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling. 
     
    630653      !                                                ! Sea surface  elevation time-stepping 
    631654      pssha(:,:) = ( psshb(:,:) - z2dt * ( r1_rho0 * pemp(:,:)  + zhdiv(:,:) ) ) * ssmask(:,:) 
    632       !                                                 !  
    633       !                                                 ! After acale factors at t-points ( z_star coordinate ) 
     655      ! 
     656      IF( PRESENT( pe3ta ) ) THEN                      ! After acale factors at t-points ( z_star coordinate ) 
    634657      DO jk = 1, jpkm1 
    635         pe3ta(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + pssha(:,:) * tmask(:,:,1) / ( ht_0(:,:) + 1.0 - tmask(:,:,1) ) ) 
     658            pe3ta(:,:,jk) = e3t_0(:,:,jk) * ( 1._wp + pssha(:,:) * r1_ht_0(:,:) * tmask(:,:,jk) ) 
    636659      END DO 
     660      ENDIF 
    637661      ! 
    638662   END SUBROUTINE dta_dyn_ssh 
     
    657681      !!---------------------------------------------------------------------- 
    658682      ! 
     683!!st code dupliqué même remarque que plus haut pourquoi ne pas utiliser gdepw ? 
    659684      DO_2D_11_11 
    660685         h_rnf(ji,jj) = 0._wp 
  • NEMO/trunk/src/OFF/nemogcm.F90

    r13011 r13237  
    2828   USE usrdef_nam     ! user defined configuration 
    2929   USE eosbn2         ! equation of state            (eos bn2 routine) 
     30#if defined key_qco 
     31   USE domqco         ! tools for scale factor         (dom_qco_r3c  routine) 
     32#endif 
    3033   USE bdy_oce,  ONLY : ln_bdy 
    3134   USE bdyini         ! open boundary cond. setting       (bdy_init routine) 
     
    119122                                CALL dta_dyn    ( istp, Nbb, Nnn, Naa )       ! Interpolation of the dynamical fields 
    120123#endif 
     124#if ! defined key_sed_off 
     125         IF( .NOT.ln_linssh ) THEN 
     126                                CALL dta_dyn_atf( istp, Nbb, Nnn, Naa )       ! time filter of sea  surface height and vertical scale factors 
     127# if defined key_qco 
     128                                CALL dom_qco_r3c( ssh(:,:,Kmm), r3t_f, r3u_f, r3v_f ) 
     129# endif 
     130         ENDIF 
    121131                                CALL trc_stp    ( istp, Nbb, Nnn, Nrhs, Naa ) ! time-stepping 
    122 #if ! defined key_sed_off 
    123          IF( .NOT.ln_linssh )   CALL dta_dyn_atf( istp, Nbb, Nnn, Naa )       ! time filter of sea  surface height and vertical scale factors 
     132# if defined key_qco 
     133                                !r3t(:,:,Kmm) = r3t_f(:,:)                     ! update ssh to h0 ratio 
     134                                !r3u(:,:,Kmm) = r3u_f(:,:) 
     135                                !r3v(:,:,Kmm) = r3v_f(:,:) 
     136# endif 
    124137#endif 
    125138         ! Swap time levels 
     
    129142         Naa = Nrhs 
    130143         ! 
     144#if ! defined key_qco 
    131145#if ! defined key_sed_off 
    132146         IF( .NOT.ln_linssh )   CALL dta_dyn_sf_interp( istp, Nnn )  ! calculate now grid parameters 
    133147#endif 
     148#endif          
    134149                                CALL stp_ctl    ( istp )             ! Time loop: control and print 
    135150         istp = istp + 1 
  • NEMO/trunk/src/TOP/C14/trcsms_c14.F90

    r12489 r13237  
    2828   !! * Substitutions 
    2929#  include "do_loop_substitute.h90" 
     30#  include "domzgr_substitute.h90" 
    3031   !!---------------------------------------------------------------------- 
    3132   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/CFC/trcsms_cfc.F90

    r12489 r13237  
    4949   !! * Substitutions 
    5050#  include "do_loop_substitute.h90" 
     51#  include "domzgr_substitute.h90" 
    5152   !!---------------------------------------------------------------------- 
    5253   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P2Z/p2zbio.F90

    r13226 r13237  
    5858   !! * Substitutions 
    5959#  include "do_loop_substitute.h90" 
     60#  include "domzgr_substitute.h90" 
    6061   !!---------------------------------------------------------------------- 
    6162   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P2Z/p2zexp.F90

    r13226 r13237  
    3939   !! * Substitutions 
    4040#  include "do_loop_substitute.h90" 
     41#  include "domzgr_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    4243   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P2Z/p2zopt.F90

    r12377 r13237  
    4040   !! * Substitutions 
    4141#  include "do_loop_substitute.h90" 
     42#  include "domzgr_substitute.h90" 
    4243   !!---------------------------------------------------------------------- 
    4344   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P2Z/p2zsed.F90

    r12377 r13237  
    3333   !! * Substitutions 
    3434#  include "do_loop_substitute.h90" 
     35#  include "domzgr_substitute.h90" 
    3536   !!---------------------------------------------------------------------- 
    3637   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zbc.F90

    r13226 r13237  
    4848   !! * Substitutions 
    4949#  include "do_loop_substitute.h90" 
     50#  include "domzgr_substitute.h90" 
    5051   !!---------------------------------------------------------------------- 
    5152   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zbio.F90

    r12377 r13237  
    4040   !! * Substitutions 
    4141#  include "do_loop_substitute.h90" 
     42#  include "domzgr_substitute.h90" 
    4243   !!---------------------------------------------------------------------- 
    4344   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zche.F90

    r12377 r13237  
    132132   !! * Substitutions 
    133133#  include "do_loop_substitute.h90" 
     134#  include "domzgr_substitute.h90" 
    134135   !!---------------------------------------------------------------------- 
    135136   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zfechem.F90

    r12377 r13237  
    3333   !! * Substitutions 
    3434#  include "do_loop_substitute.h90" 
     35#  include "domzgr_substitute.h90" 
    3536   !!---------------------------------------------------------------------- 
    3637   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zflx.F90

    r12377 r13237  
    5454   !! * Substitutions 
    5555#  include "do_loop_substitute.h90" 
     56#  include "domzgr_substitute.h90" 
    5657   !!---------------------------------------------------------------------- 
    5758   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zopt.F90

    r13226 r13237  
    4444   !! * Substitutions 
    4545#  include "do_loop_substitute.h90" 
     46#  include "domzgr_substitute.h90" 
    4647   !!---------------------------------------------------------------------- 
    4748   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zpoc.F90

    r12377 r13237  
    3939   !! * Substitutions 
    4040#  include "do_loop_substitute.h90" 
     41#  include "domzgr_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    4243   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zprod.F90

    r12377 r13237  
    4848   !! * Substitutions 
    4949#  include "do_loop_substitute.h90" 
     50#  include "domzgr_substitute.h90" 
    5051   !!---------------------------------------------------------------------- 
    5152   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zrem.F90

    r12377 r13237  
    4444   !! * Substitutions 
    4545#  include "do_loop_substitute.h90" 
     46#  include "domzgr_substitute.h90" 
    4647   !!---------------------------------------------------------------------- 
    4748   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsed.F90

    r12377 r13237  
    3939   !! * Substitutions 
    4040#  include "do_loop_substitute.h90" 
     41#  include "domzgr_substitute.h90" 
    4142   !!---------------------------------------------------------------------- 
    4243   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsink.F90

    r12377 r13237  
    4040   !! * Substitutions 
    4141#  include "do_loop_substitute.h90" 
     42#  include "domzgr_substitute.h90" 
    4243   !!---------------------------------------------------------------------- 
    4344   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p4zsms.F90

    r13030 r13237  
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/P4Z/p5zprod.F90

    r12377 r13237  
    5252   !! * Substitutions 
    5353#  include "do_loop_substitute.h90" 
     54#  include "domzgr_substitute.h90" 
    5455   !!---------------------------------------------------------------------- 
    5556   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/PISCES/SED/oce_sed.F90

    r12489 r13237  
    1313   USE dom_oce , ONLY :   glamt     =>   glamt          !: longitude of t-point (degre) 
    1414   USE dom_oce , ONLY :   gphit     =>   gphit          !: latitude  of t-point (degre) 
     15!!st  
     16#if ! defined key_qco 
    1517   USE dom_oce , ONLY :   e3t       =>   e3t            !: latitude  of t-point (degre) 
     18#endif 
    1619   USE dom_oce , ONLY :   e3t_1d    =>   e3t_1d         !: reference depth of t-points (m) 
    1720   USE dom_oce , ONLY :   gdepw_0   =>   gdepw_0        !: reference depth of t-points (m) 
     
    5356 
    5457END MODULE oce_sed 
    55  
    56  
  • NEMO/trunk/src/TOP/PISCES/SED/seddta.F90

    r12489 r13237  
    2424   !! * Substitutions 
    2525#  include "do_loop_substitute.h90" 
     26#  include "domzgr_substitute.h90" 
    2627   !! $Id$ 
    2728CONTAINS 
     
    164165      CALL pack_arr ( jpoce,  rainrm_dta(1:jpoce,jscal), trc_data(1:jpi,1:jpj,14), iarroce(1:jpoce) ) 
    165166      rainrm_dta(1:jpoce,jscal) = rainrm_dta(1:jpoce,jscal) * 1e-4 
    166       ! vector temperature [°C] and salinity  
     167      ! vector temperature [C] and salinity  
    167168      CALL pack_arr ( jpoce,  temp(1:jpoce), trc_data(1:jpi,1:jpj,15), iarroce(1:jpoce) ) 
    168169      CALL pack_arr ( jpoce,  salt(1:jpoce), trc_data(1:jpi,1:jpj,16), iarroce(1:jpoce) ) 
  • NEMO/trunk/src/TOP/PISCES/trcwri_pisces.F90

    r12377 r13237  
    2121   !! * Substitutions 
    2222#  include "do_loop_substitute.h90" 
     23#  include "domzgr_substitute.h90" 
    2324   !!---------------------------------------------------------------------- 
    2425   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/TRP/trcadv.F90

    r12489 r13237  
    5959   INTEGER, PARAMETER ::   np_QCK     = 5   ! QUICK scheme 
    6060    
     61#  include "domzgr_substitute.h90" 
    6162   !!---------------------------------------------------------------------- 
    6263   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/TRP/trcatf.F90

    r12489 r13237  
    3131   USE trd_oce 
    3232   USE trdtra 
     33# if defined key_qco 
     34   USE traatfqco 
     35# else 
    3336   USE traatf 
     37# endif 
    3438   USE bdy_oce   , ONLY: ln_bdy 
    3539   USE trcbdy          ! BDY open boundaries 
     
    5054   !! * Substitutions 
    5155#  include "do_loop_substitute.h90" 
     56#  include "domzgr_substitute.h90" 
    5257   !!---------------------------------------------------------------------- 
    5358   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    112117         ! total trend for the non-time-filtered variables.  
    113118         zfact = 1.0 / rn_Dt 
    114          ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from ts(Kmm) terms 
     119         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3ta*Ta)/e3tn; e3tn cancel from ts(Kmm) terms 
    115120         IF( ln_linssh ) THEN       ! linear sea surface height only 
    116121            DO jn = 1, jptra 
     
    151156      ELSE      
    152157         IF( .NOT. l_offline ) THEN ! Leap-Frog + Asselin filter time stepping 
     158# if defined key_qco 
     159            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix_lf( kt, Kbb, Kmm, Kaa, nittrc000,        'TRC', ptr, jptra )                     !     linear ssh 
     160            ELSE                   ;   CALL tra_atf_qco_lf( kt, Kbb, Kmm, Kaa, nittrc000, rn_Dt, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh 
     161# else 
    153162            IF( ln_linssh ) THEN   ;   CALL tra_atf_fix( kt, Kbb, Kmm, Kaa, nittrc000,         'TRC', ptr, jptra )                     !     linear ssh 
    154163            ELSE                   ;   CALL tra_atf_vvl( kt, Kbb, Kmm, Kaa, nittrc000, rn_Dt, 'TRC', ptr, sbc_trc, sbc_trc_b, jptra ) ! non-linear ssh 
     164# endif 
    155165            ENDIF 
    156166         ELSE 
     
    182192   END SUBROUTINE trc_atf 
    183193 
    184  
     194# if ! defined key_qco 
    185195   SUBROUTINE trc_atf_off( kt, Kbb, Kmm, Kaa, ptr ) 
    186196      !!---------------------------------------------------------------------- 
     
    198208      !!                This can be summurized for tempearture as: 
    199209      !!             ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )   ln_dynhpg_imp = T 
    200       !!                  /( e3t(:,:,:,Kmm)    + rbcp*[ e3t(:,:,:,Kbb)    - 2 e3t(:,:,:,Kmm)    + e3t(:,:,:,Kaa)    ] )    
     210      !!                  /( e3t(:,:,jk,Kmm)    + rbcp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )    
    201211      !!             ztm = 0                                                       otherwise 
    202212      !!             tb  = ( e3t_n*tn + rn_atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) 
    203       !!                  /( e3t(:,:,:,Kmm)    + rn_atfp*[ e3t(:,:,:,Kbb)    - 2 e3t(:,:,:,Kmm)    + e3t(:,:,:,Kaa)    ] ) 
     213      !!                  /( e3t(:,:,jk,Kmm)    + rn_atfp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] ) 
    204214      !!             tn  = ta  
    205215      !!             ta  = zt        (NB: reset to 0 after eos_bn2 call) 
     
    257267      ! 
    258268   END SUBROUTINE trc_atf_off 
     269# else 
     270   SUBROUTINE trc_atf_off( kt, Kbb, Kmm, Kaa, ptr ) 
     271      !!---------------------------------------------------------------------- 
     272      !!                   ***  ROUTINE tra_atf_off  *** 
     273      !! 
     274      !!          !!!!!!!!!!!!!!!!! REWRITE HEADER COMMENTS !!!!!!!!!!!!!! 
     275      !! 
     276      !! ** Purpose :   Time varying volume: apply the Asselin time filter   
     277      !!  
     278      !! ** Method  : - Apply a thickness weighted Asselin time filter on now fields. 
     279      !!              - save in (ta,sa) a thickness weighted average over the three  
     280      !!             time levels which will be used to compute rdn and thus the semi- 
     281      !!             implicit hydrostatic pressure gradient (ln_dynhpg_imp = T) 
     282      !!              - swap tracer fields to prepare the next time_step. 
     283      !!                This can be summurized for tempearture as: 
     284      !!             ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )   ln_dynhpg_imp = T 
     285      !!                  /( e3t(:,:,jk,Kmm)    + rbcp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] )    
     286      !!             ztm = 0                                                       otherwise 
     287      !!             tb  = ( e3t_n*tn + rn_atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] ) 
     288      !!                  /( e3t(:,:,jk,Kmm)    + rn_atfp*[ e3t(:,:,jk,Kbb)    - 2 e3t(:,:,jk,Kmm)    + e3t(:,:,jk,Kaa)    ] ) 
     289      !!             tn  = ta  
     290      !!             ta  = zt        (NB: reset to 0 after eos_bn2 call) 
     291      !! 
     292      !! ** Action  : - (tb,sb) and (tn,sn) ready for the next time step 
     293      !!              - (ta,sa) time averaged (t,s)   (ln_dynhpg_imp = T) 
     294      !!---------------------------------------------------------------------- 
     295      INTEGER                                   , INTENT(in   ) ::  kt            ! ocean time-step index 
     296      INTEGER                                   , INTENT(in   ) ::  Kbb, Kmm, Kaa ! time level indices 
     297      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) ::  ptr           ! passive tracers 
     298      !!      
     299      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices 
     300      REAL(wp) ::   ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar 
     301      REAL(wp) ::   ze3t_b, ze3t_n, ze3t_a, ze3t_f           !   -      - 
     302      !!---------------------------------------------------------------------- 
     303      ! 
     304      IF( kt == nittrc000 )  THEN 
     305         IF(lwp) WRITE(numout,*) 
     306         IF(lwp) WRITE(numout,*) 'trc_atf_off : Asselin time filtering' 
     307         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 
     308         IF( .NOT. ln_linssh ) THEN 
     309            rfact1 = rn_atfp * rn_Dt 
     310            rfact2 = rfact1 / rho0 
     311         ENDIF 
     312        !   
     313      ENDIF 
     314      ! 
     315      DO jn = 1, jptra       
     316         DO_3D_11_11( 1, jpkm1 ) 
     317            ze3t_b = 1._wp + r3t(ji,jj,Kbb) * tmask(ji,jj,jk) 
     318            ze3t_n = 1._wp + r3t(ji,jj,Kmm) * tmask(ji,jj,jk) 
     319            ze3t_a = 1._wp + r3t(ji,jj,Kaa) * tmask(ji,jj,jk) 
     320            !                                         ! tracer content at Before, now and after 
     321            ztc_b  = ptr(ji,jj,jk,jn,Kbb) * ze3t_b 
     322            ztc_n  = ptr(ji,jj,jk,jn,Kmm) * ze3t_n 
     323            ztc_a  = ptr(ji,jj,jk,jn,Kaa) * ze3t_a 
     324            ! 
     325            ztc_d  = ztc_a  - 2. * ztc_n  + ztc_b 
     326            ! 
     327            ze3t_f = 1._wp + r3t_f(ji,jj)*tmask(ji,jj,jk) 
     328            ztc_f  = ztc_n  + rn_atfp * ztc_d 
     329            ! 
     330            IF( .NOT. ln_linssh .AND. jk == mikt(ji,jj) ) THEN           ! first level  
     331               ztc_f  = ztc_f  - rfact1 * ( sbc_trc(ji,jj,jn) - sbc_trc_b(ji,jj,jn) ) 
     332            ENDIF 
     333 
     334            ze3t_f = 1.e0 / ze3t_f 
     335            ptr(ji,jj,jk,jn,Kmm) = ztc_f * ze3t_f     ! time filtered "now" field 
     336            ! 
     337         END_3D 
     338         !  
     339      END DO 
     340      ! 
     341   END SUBROUTINE trc_atf_off 
     342# endif 
    259343 
    260344#else 
  • NEMO/trunk/src/TOP/TRP/trcdmp.F90

    r12377 r13237  
    4545   !! * Substitutions 
    4646#  include "do_loop_substitute.h90" 
     47#  include "domzgr_substitute.h90" 
    4748   !!---------------------------------------------------------------------- 
    4849   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/TRP/trcldf.F90

    r12377 r13237  
    4444   !! * Substitutions 
    4545#  include "do_loop_substitute.h90" 
     46#  include "domzgr_substitute.h90" 
    4647   !!---------------------------------------------------------------------- 
    4748   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/TRP/trcsbc.F90

    r13226 r13237  
    3030   !! * Substitutions 
    3131#  include "do_loop_substitute.h90" 
     32#  include "domzgr_substitute.h90" 
    3233   !!---------------------------------------------------------------------- 
    3334   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    4950      !!            The surface freshwater flux modify the ocean volume 
    5051      !!         and thus the concentration of a tracer as : 
    51       !!            tr(Krhs) = tr(Krhs) + emp * tr(Kmm) / e3t   for k=1 
     52      !!            tr(Krhs) = tr(Krhs) + emp * tr(Kmm) / e3t_   for k=1 
    5253      !!         where emp, the surface freshwater budget (evaporation minus 
    5354      !!         precipitation ) given in kg/m2/s is divided 
  • NEMO/trunk/src/TOP/TRP/trcsink.F90

    r13226 r13237  
    2626   !! * Substitutions 
    2727#  include "do_loop_substitute.h90" 
     28#  include "domzgr_substitute.h90" 
    2829   !!---------------------------------------------------------------------- 
    2930   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/TRP/trdmxl_trc.F90

    r12489 r13237  
    5151   !! * Substitutions 
    5252#  include "do_loop_substitute.h90" 
     53#  include "domzgr_substitute.h90" 
    5354   !!---------------------------------------------------------------------- 
    5455   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/trcbc.F90

    r12852 r13237  
    4848   !! * Substitutions 
    4949#  include "do_loop_substitute.h90" 
     50#  include "domzgr_substitute.h90" 
    5051   !!---------------------------------------------------------------------- 
    5152   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/trcdta.F90

    r12377 r13237  
    4141   !! Substitutions 
    4242#include "do_loop_substitute.h90" 
     43#include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
     
    206207                     ztp(jk) = ptrcdta(ji,jj,jpkm1) 
    207208                  ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1 
    208                      DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1) 
     209                     DO jkk = 1, jpkm1                                  ! when  gdept_1d(jkk) < zl < gdept_1d(jkk+1) 
    209210                        IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 
    210211                           zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 
  • NEMO/trunk/src/TOP/trcini.F90

    r12377 r13237  
    3131   PUBLIC   trc_init   ! called by opa 
    3232 
     33#  include "domzgr_substitute.h90" 
    3334   !!---------------------------------------------------------------------- 
    3435   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/trcrst.F90

    r12489 r13237  
    3333   PUBLIC   trc_rst_cal 
    3434 
     35#  include "domzgr_substitute.h90" 
    3536   !!---------------------------------------------------------------------- 
    3637   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/trcstp.F90

    r12620 r13237  
    3737   REAL(wp), DIMENSION(:,:,:), SAVE, ALLOCATABLE ::   qsr_arr   ! save qsr during TOP time-step 
    3838 
     39#  include "domzgr_substitute.h90" 
    3940   !!---------------------------------------------------------------------- 
    4041   !! NEMO/TOP 4.0 , NEMO Consortium (2018) 
  • NEMO/trunk/src/TOP/trcwri.F90

    r12377 r13237  
    6060       CALL iom_put( "e3v_0", e3v_0(:,:,:) ) 
    6161       ! 
     62#if ! defined key_qco 
    6263       CALL iom_put( "e3t" , e3t(:,:,:,Kmm) ) 
    6364       CALL iom_put( "e3u" , e3u(:,:,:,Kmm) ) 
    6465       CALL iom_put( "e3v" , e3v(:,:,:,Kmm) ) 
     66#endif  
    6567       ! 
    6668      ENDIF 
Note: See TracChangeset for help on using the changeset viewer.