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 13540 for NEMO/branches/2020/r12377_ticket2386/src/OCE/DOM/domvvl.F90 – NEMO

Ignore:
Timestamp:
2020-09-29T12:41:06+02:00 (4 years ago)
Author:
andmirek
Message:

Ticket #2386: update to latest trunk

Location:
NEMO/branches/2020/r12377_ticket2386
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/r12377_ticket2386

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
        88 
        99# SETTE 
        10 ^/utils/CI/sette@HEAD         sette 
         10^/utils/CI/sette@13507        sette 
  • NEMO/branches/2020/r12377_ticket2386/src/OCE/DOM/domvvl.F90

    r12511 r13540  
    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      !!---------------------------------------------------------------------- 
     
    190202      gdept(:,:,1,Kbb) = 0.5_wp * e3w(:,:,1,Kbb) 
    191203      gdepw(:,:,1,Kbb) = 0.0_wp 
    192       DO_3D_11_11( 2, jpk ) 
     204      DO_3D( 1, 1, 1, 1, 2, jpk )                     ! vertical sum 
    193205         !    zcoef = tmask - wmask    ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
    194206         !                             ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 
     
    238250         ENDIF 
    239251         IF ( ln_vvl_zstar_at_eqtor ) THEN   ! use z-star in vicinity of the Equator 
    240             DO_2D_11_11 
     252            DO_2D( 1, 1, 1, 1 ) 
    241253!!gm  case |gphi| >= 6 degrees is useless   initialized just above by default 
    242254               IF( ABS(gphit(ji,jj)) >= 6.) THEN 
     
    261273            IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN 
    262274               IF( nn_cfg == 3 ) THEN   ! ORCA2: Suppress ztilde in the Foxe Basin for ORCA2 
    263                   ii0 = 103   ;   ii1 = 111        
    264                   ij0 = 128   ;   ij1 = 135   ;    
     275                  ii0 = 103 + nn_hls - 1   ;   ii1 = 111 + nn_hls - 1       
     276                  ij0 = 128 + nn_hls       ;   ij1 = 135 + nn_hls 
    265277                  frq_rst_e3t( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  0.0_wp 
    266278                  frq_rst_hdv( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  1.e0_wp / rn_Dt 
     
    322334      LOGICAL                ::   ll_do_bclinic         ! local logical 
    323335      REAL(wp), DIMENSION(jpi,jpj)     ::   zht, z_scale, zwu, zwv, zhdiv 
    324       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ze3t 
     336      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ze3t 
     337      LOGICAL , DIMENSION(:,:,:), ALLOCATABLE ::   llmsk 
    325338      !!---------------------------------------------------------------------- 
    326339      ! 
     
    407420         zwu(:,:) = 0._wp 
    408421         zwv(:,:) = 0._wp 
    409          DO_3D_10_10( 1, jpkm1 ) 
     422         DO_3D( 1, 0, 1, 0, 1, jpkm1 )   ! a - first derivative: diffusive fluxes 
    410423            un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj)           & 
    411424               &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
     
    415428            zwv(ji,jj) = zwv(ji,jj) + vn_td(ji,jj,jk) 
    416429         END_3D 
    417          DO_2D_11_11 
     430         DO_2D( 1, 1, 1, 1 )             ! b - correction for last oceanic u-v points 
    418431            un_td(ji,jj,mbku(ji,jj)) = un_td(ji,jj,mbku(ji,jj)) - zwu(ji,jj) 
    419432            vn_td(ji,jj,mbkv(ji,jj)) = vn_td(ji,jj,mbkv(ji,jj)) - zwv(ji,jj) 
    420433         END_2D 
    421          DO_3D_00_00( 1, jpkm1 ) 
     434         DO_3D( 0, 0, 0, 0, 1, jpkm1 )   ! c - second derivative: divergence of diffusive fluxes 
    422435            tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   un_td(ji-1,jj  ,jk) - un_td(ji,jj,jk)    & 
    423436               &                                          +     vn_td(ji  ,jj-1,jk) - vn_td(ji,jj,jk)    & 
    424437               &                                            ) * r1_e1e2t(ji,jj) 
    425438         END_3D 
    426          !                       ! d - thickness diffusion transport: boundary conditions 
     439         !                               ! d - thickness diffusion transport: boundary conditions 
    427440         !                             (stored for tracer advction and continuity equation) 
    428441         CALL lbc_lnk_multi( 'domvvl', un_td , 'U' , -1._wp, vn_td , 'V' , -1._wp) 
     
    435448         ! Maximum deformation control 
    436449         ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~ 
    437          ze3t(:,:,jpk) = 0._wp 
    438          DO jk = 1, jpkm1 
    439             ze3t(:,:,jk) = tilde_e3t_a(:,:,jk) / e3t_0(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) 
    440          END DO 
    441          z_tmax = MAXVAL( ze3t(:,:,:) ) 
    442          CALL mpp_max( 'domvvl', z_tmax )                 ! max over the global domain 
    443          z_tmin = MINVAL( ze3t(:,:,:) ) 
    444          CALL mpp_min( 'domvvl', z_tmin )                 ! min over the global domain 
     450         ALLOCATE( ze3t(jpi,jpj,jpk), llmsk(jpi,jpj,jpk) ) 
     451         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     452            ze3t(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) / e3t_0(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) 
     453         END_3D 
     454         ! 
     455         llmsk(   1:Nis1,:,:) = .FALSE.   ! exclude halos from the checked region 
     456         llmsk(Nie1: jpi,:,:) = .FALSE. 
     457         llmsk(:,   1:Njs1,:) = .FALSE. 
     458         llmsk(:,Nje1: jpj,:) = .FALSE. 
     459         ! 
     460         llmsk(Nis0:Nie0,Njs0:Nje0,:) = tmask(Nis0:Nie0,Njs0:Nje0,:) == 1._wp                  ! define only the inner domain 
     461         z_tmax = MAXVAL( ze3t(:,:,:), mask = llmsk )   ;   CALL mpp_max( 'domvvl', z_tmax )   ! max over the global domain 
     462         z_tmin = MINVAL( ze3t(:,:,:), mask = llmsk )   ;   CALL mpp_min( 'domvvl', z_tmin )   ! min over the global domain 
    445463         ! - ML - test: for the moment, stop simulation for too large e3_t variations 
    446464         IF( ( z_tmax >  rn_zdef_max ) .OR. ( z_tmin < - rn_zdef_max ) ) THEN 
    447             IF( lk_mpp ) THEN 
    448                CALL mpp_maxloc( 'domvvl', ze3t, tmask, z_tmax, ijk_max ) 
    449                CALL mpp_minloc( 'domvvl', ze3t, tmask, z_tmin, ijk_min ) 
    450             ELSE 
    451                ijk_max = MAXLOC( ze3t(:,:,:) ) 
    452                ijk_max(1) = ijk_max(1) + nimpp - 1 
    453                ijk_max(2) = ijk_max(2) + njmpp - 1 
    454                ijk_min = MINLOC( ze3t(:,:,:) ) 
    455                ijk_min(1) = ijk_min(1) + nimpp - 1 
    456                ijk_min(2) = ijk_min(2) + njmpp - 1 
    457             ENDIF 
     465            CALL mpp_maxloc( 'domvvl', ze3t, llmsk, z_tmax, ijk_max ) 
     466            CALL mpp_minloc( 'domvvl', ze3t, llmsk, z_tmin, ijk_min ) 
    458467            IF (lwp) THEN 
    459468               WRITE(numout, *) 'MAX( tilde_e3t_a(:,:,:) / e3t_0(:,:,:) ) =', z_tmax 
     
    464473            ENDIF 
    465474         ENDIF 
     475         DEALLOCATE( ze3t, llmsk ) 
    466476         ! - ML - end test 
    467477         ! - ML - Imposing these limits will cause a baroclinicity error which is corrected for below 
     
    647657      gdepw(:,:,1,Kmm) = 0.0_wp 
    648658      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    649       DO_3D_11_11( 2, jpk ) 
     659      DO_3D( 1, 1, 1, 1, 2, jpk ) 
    650660        !    zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))   ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
    651661                                                           ! 1 for jk = mikt 
     
    702712         ! 
    703713      CASE( 'U' )                   !* from T- to U-point : hor. surface weighted mean 
    704          DO_3D_10_10( 1, jpk ) 
     714         DO_3D( 1, 0, 1, 0, 1, jpk ) 
    705715            pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2u(ji,jj)   & 
    706716               &                       * (   e1e2t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
     
    711721         ! 
    712722      CASE( 'V' )                   !* from T- to V-point : hor. surface weighted mean 
    713          DO_3D_10_10( 1, jpk ) 
     723         DO_3D( 1, 0, 1, 0, 1, jpk ) 
    714724            pe3_out(ji,jj,jk) = 0.5_wp * ( vmask(ji,jj,jk)  * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2v(ji,jj)   & 
    715725               &                       * (   e1e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
     
    720730         ! 
    721731      CASE( 'F' )                   !* from U-point to F-point : hor. surface weighted mean 
    722          DO_3D_10_10( 1, jpk ) 
     732         DO_3D( 1, 0, 1, 0, 1, jpk ) 
    723733            pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * umask(ji,jj+1,jk) * (1.0_wp - zlnwd) + zlnwd ) & 
    724734               &                       *    r1_e1e2f(ji,jj)                                                  & 
     
    793803         IF( ln_rstart ) THEN                   !* Read the restart file 
    794804            CALL rst_read_open                  !  open the restart file if necessary 
    795             CALL iom_get( numror, jpdom_autoglo, 'sshn'   , ssh(:,:,Kmm), ldxios = lrxios    ) 
     805            CALL iom_get( numror, jpdom_auto, 'sshn'   , ssh(:,:,Kmm), ldxios = lrxios    ) 
    796806            ! 
    797807            id1 = iom_varid( numror, 'e3t_b', ldstop = .FALSE. ) 
     
    806816            ! 
    807817            IF( MIN( id1, id2 ) > 0 ) THEN       ! all required arrays exist 
    808                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
    809                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
     818               CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     819               CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
    810820               ! needed to restart if land processor not computed  
    811821               IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files' 
     
    821831               IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.' 
    822832               IF(lwp) write(numout,*) 'l_1st_euler is forced to true' 
    823                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     833               CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
    824834               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
    825835               l_1st_euler = .true. 
     
    828838               IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.' 
    829839               IF(lwp) write(numout,*) 'l_1st_euler is forced to true' 
    830                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
     840               CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
    831841               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    832842               l_1st_euler = .true. 
     
    853863               !                          ! ----------------------- ! 
    854864               IF( MIN( id3, id4 ) > 0 ) THEN  ! all required arrays exist 
    855                   CALL iom_get( numror, jpdom_autoglo, 'tilde_e3t_b', tilde_e3t_b(:,:,:), ldxios = lrxios ) 
    856                   CALL iom_get( numror, jpdom_autoglo, 'tilde_e3t_n', tilde_e3t_n(:,:,:), ldxios = lrxios ) 
     865                  CALL iom_get( numror, jpdom_auto, 'tilde_e3t_b', tilde_e3t_b(:,:,:), ldxios = lrxios ) 
     866                  CALL iom_get( numror, jpdom_auto, 'tilde_e3t_n', tilde_e3t_n(:,:,:), ldxios = lrxios ) 
    857867               ELSE                            ! one at least array is missing 
    858868                  tilde_e3t_b(:,:,:) = 0.0_wp 
     
    863873                  !                       ! ------------ ! 
    864874                  IF( id5 > 0 ) THEN  ! required array exists 
    865                      CALL iom_get( numror, jpdom_autoglo, 'hdiv_lf', hdiv_lf(:,:,:), ldxios = lrxios ) 
     875                     CALL iom_get( numror, jpdom_auto, 'hdiv_lf', hdiv_lf(:,:,:), ldxios = lrxios ) 
    866876                  ELSE                ! array is missing 
    867877                     hdiv_lf(:,:,:) = 0.0_wp 
     
    887897                  ssh(:,:,Kbb) = -ssh_ref 
    888898 
    889                   DO_2D_11_11 
     899                  DO_2D( 1, 1, 1, 1 ) 
    890900                     IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth 
    891901                        ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 
     
    903913               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    904914 
    905                DO ji = 1, jpi 
    906                   DO jj = 1, jpj 
    907                      IF ( ht_0(ji,jj) .LE. 0.0 .AND. NINT( ssmask(ji,jj) ) .EQ. 1) THEN 
    908                        CALL ctl_stop( 'dom_vvl_rst: ht_0 must be positive at potentially wet points' ) 
    909                      ENDIF 
    910                   END DO  
    911                END DO  
     915               DO_2D( 1, 1, 1, 1 ) 
     916                  IF ( ht_0(ji,jj) .LE. 0.0 .AND. NINT( ssmask(ji,jj) ) .EQ. 1) THEN 
     917                     CALL ctl_stop( 'dom_vvl_rst: ht_0 must be positive at potentially wet points' ) 
     918                  ENDIF 
     919               END_2D 
    912920               ! 
    913921            ELSE 
     
    10311039   END SUBROUTINE dom_vvl_ctl 
    10321040 
     1041#endif 
     1042 
    10331043   !!====================================================================== 
    10341044END MODULE domvvl 
Note: See TracChangeset for help on using the changeset viewer.