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 10164 for NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE – NEMO

Ignore:
Timestamp:
2018-10-02T16:02:24+02:00 (6 years ago)
Author:
jchanut
Message:

ztilde update, #2126: Add higher order filtering + alternative scale factor decomposition

Location:
NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE/cfgs/SHARED/field_def_nemo-oce.xml

    r9990 r10164  
    2121        <field id="e3t"          long_name="T-cell thickness"                    standard_name="cell_thickness"     unit="m"   grid_ref="grid_T_3D" /> 
    2222        <field id="e3t_surf"     long_name="T-cell thickness"   field_ref="e3t"  standard_name="cell_thickness"     unit="m"   grid_ref="grid_T_SFC"/> 
     23        <field id="e3t_0"        long_name="Initial T-cell thickness"            standard_name="ref_cell_thickness" unit="m"   grid_ref="grid_T_3D" /> 
     24        <field id="e3t_star"     long_name="Barotropic T-cell thickness anomaly" standard_name="barotropic_cell_thickness_anomaly" unit="m"   grid_ref="grid_T_3D" /> 
     25        <field id="e3t_tilde"    long_name="Baroclinic T-cell thickness anomaly" standard_name="baroclinic_cell_thickness_anomaly" unit="m"   grid_ref="grid_T_3D" /> 
    2326        <field id="e3t_0"        long_name="Initial T-cell thickness"            standard_name="ref_cell_thickness" unit="m"   grid_ref="grid_T_3D" /> 
    2427 
  • NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE/src/OCE/DIA/diawri.F90

    r10126 r10164  
    5858   USE diurnal_bulk    ! diurnal warm layer 
    5959   USE cool_skin       ! Cool skin 
     60   USE domvvl 
    6061 
    6162   IMPLICIT NONE 
     
    124125 
    125126      ! Output of initial vertical scale factor 
    126       CALL iom_put("e3t_0", e3t_0(:,:,:) ) 
     127      CALL iom_put("e3t_0", e3t_0(:,:,:)) 
    127128      CALL iom_put("e3u_0", e3u_0(:,:,:) ) 
    128129      CALL iom_put("e3v_0", e3v_0(:,:,:) ) 
     
    132133      CALL iom_put( "e3v" , e3v_n(:,:,:) ) 
    133134      CALL iom_put( "e3w" , e3w_n(:,:,:) ) 
     135      ! 
     136      IF ( (ln_vvl_ztilde).OR.(ln_vvl_layer) ) THEN 
     137         IF( iom_use("e3t_star") ) THEN    ! Barotropic cell thickness anomaly 
     138            z3d(:,:,:) = (e3t_n(:,:,:)-tilde_e3t_n(:,:,:)-e3t_0(:,:,:))*tmask(:,:,:)  
     139            CALL iom_put( "e3t_star" , z3d(:,:,:) ) 
     140         ENDIF 
     141         IF( iom_use("e3t_tilde") )  THEN  ! Baroclinic cell thickness anomaly 
     142            CALL iom_put( "e3t_tilde" , tilde_e3t_n(:,:,:) ) 
     143         ENDIF 
     144      ENDIF 
     145 
    134146      IF( iom_use("e3tdef") )   & 
    135147         CALL iom_put( "e3tdef"  , ( ( e3t_n(:,:,:) - e3t_0(:,:,:) ) / e3t_0(:,:,:) * 100 * tmask(:,:,:) ) ** 2 ) 
  • NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE/src/OCE/DOM/domvvl.F90

    r10128 r10164  
    6161   LOGICAL          :: ln_vvl_kepe            = .FALSE.    ! kinetic/potential energy transfer 
    6262   !                                                       ! conservation: not used yet 
     63   INTEGER          :: nn_filt_order=1 
    6364   REAL(wp)         :: rn_ahe3_lap               ! thickness diffusion coefficient (Laplacian) 
    6465   REAL(wp)         :: rn_ahe3_blp               ! thickness diffusion coefficient (Bilaplacian) 
     
    6970 
    7071   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_td, vn_td                ! thickness diffusion transport 
    71    REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: un_lf, vn_lf, hdivn_lf      ! low frequency fluxes and divergence 
     72   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: un_lf, vn_lf, hdivn_lf    ! low frequency fluxes and divergence 
    7273   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_b, tilde_e3t_n    ! baroclinic scale factors 
    73    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a, dtilde_e3t_a   ! baroclinic scale factors 
     74   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tilde_e3t_a                 ! baroclinic scale factors 
    7475   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: tildemask                   ! mask tilde tendency 
    7576   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_e3t                 ! restoring period for scale factors 
     
    9495      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 
    9596         ALLOCATE( tilde_e3t_b(jpi,jpj,jpk)  , tilde_e3t_n(jpi,jpj,jpk) , tilde_e3t_a(jpi,jpj,jpk) ,   & 
    96             &      dtilde_e3t_a(jpi,jpj,jpk) , un_td  (jpi,jpj,jpk)     , vn_td  (jpi,jpj,jpk)     ,   & 
     97            &      un_td  (jpi,jpj,jpk)      , vn_td  (jpi,jpj,jpk)     ,                              & 
    9798            &      tildemask(jpi,jpj) , hsm(jpi,jpj) , dsm(jpi,jpj) , i_int_bot(jpi,jpj), STAT = dom_vvl_alloc ) 
    9899         IF( lk_mpp             )   CALL mpp_sum ( dom_vvl_alloc ) 
     
    102103      ENDIF 
    103104      IF( ln_vvl_ztilde ) THEN 
    104          ALLOCATE( frq_rst_e3t(jpi,jpj) , frq_rst_hdv(jpi,jpj) , hdivn_lf(jpi,jpj,jpk)             ,    & 
    105             &      un_lf(jpi,jpj,jpk), vn_lf(jpi,jpj,jpk)      , STAT= dom_vvl_alloc ) 
     105         ALLOCATE( frq_rst_e3t(jpi,jpj) , frq_rst_hdv(jpi,jpj) , hdivn_lf(jpi,jpj,jpk,nn_filt_order),  & 
     106            &      un_lf(jpi,jpj,jpk,nn_filt_order), vn_lf(jpi,jpj,jpk,nn_filt_order), STAT= dom_vvl_alloc ) 
    106107         IF( lk_mpp             )   CALL mpp_sum ( dom_vvl_alloc ) 
    107108         IF( dom_vvl_alloc /= 0 )   CALL ctl_warn('dom_vvl_alloc: failed to allocate arrays') 
     
    263264         ! 
    264265         IF ( ln_vvl_zstar_on_shelf ) THEN 
    265             zhmin =  50._wp 
     266            zhmin = 50._wp 
    266267            zhmax = 100._wp 
    267268            DO jj = 1, jpj 
     
    284285         IF ( (ztmp*rdt) > 1._wp) CALL ctl_stop( 'dom_vvl_init: rn_lf_cuttoff is too small' ) 
    285286         ! 
     287      ENDIF 
     288 
     289      IF( ln_vvl_layer ) THEN 
     290         IF ( ln_vvl_zstar_on_shelf ) THEN 
     291            zhmin = 100._wp 
     292            zhmax = 150._wp 
     293            DO jj = 1, jpj 
     294               DO ji = 1, jpi 
     295                  zwgt = 1._wp 
     296                  IF(( ht_0(ji,jj)>zhmin).AND.(ht_0(ji,jj) <=zhmax)) THEN 
     297                     zwgt = (ht_0(ji,jj)-zhmin)/(zhmax-zhmin) 
     298                  ELSEIF ( ht_0(ji,jj)<=zhmin) THEN 
     299                     zwgt = 0._wp 
     300                  ENDIF 
     301                  tildemask(ji,jj)   = MIN(tildemask(ji,jj), zwgt) 
     302               END DO 
     303            END DO 
     304         ENDIF 
     305         IF ( ln_vvl_zstar_at_eqtor ) THEN 
     306            DO jj = 1, jpj 
     307               DO ji = 1, jpi 
     308!!gm  case |gphi| >= 6 degrees is useless   initialized just above by default   
     309                  IF( ABS(gphit(ji,jj)) >= 6.) THEN 
     310                     ! values outside the equatorial band and transition zone (ztilde) 
     311               
     312                  ELSEIF( ABS(gphit(ji,jj)) <= 2.5) THEN 
     313                     ! values inside the equatorial band (ztilde as zstar) 
     314                     tildemask(ji,jj) = 0._wp 
     315                  ELSE 
     316                     tildemask(ji,jj) = 0.5_wp * (  1._wp  - COS( rad*(ABS(gphit(ji,jj))-2.5_wp)  & 
     317                        &                                                 * 180._wp / 3.5_wp ) ) 
     318                  ENDIF 
     319               END DO 
     320            END DO 
     321         ENDIF 
    286322      ENDIF 
    287323      ! 
     
    335371      ! 
    336372      LOGICAL                ::   ll_do_bclinic         ! local logical 
    337       INTEGER                ::   ji, jj, jk            ! dummy loop indices 
     373      INTEGER                ::   ji, jj, jk, jo        ! dummy loop indices 
    338374      INTEGER                ::   ib, ib_bdy, ip, jp    !   "     "     " 
    339375      INTEGER , DIMENSION(3) ::   ijk_max, ijk_min      ! temporary integers 
     
    341377      REAL(wp)               ::   z2dt  , z_tmin, z_tmax! local scalars                 
    342378      REAL(wp)               ::   zalpha, zwgt          ! temporary scalars 
    343       REAL(wp)               ::   zdu, zdv, zramp 
     379      REAL(wp)               ::   zdu, zdv, zramp, zmet 
    344380      REAL(wp)               ::   ztra, zbtr, ztout, ztin, zfac, zmsku, zmskv 
    345381      REAL(wp), DIMENSION(jpi,jpj)     ::   zht, z_scale, zwu, zwv, zhdiv 
     
    357393      ENDIF 
    358394 
     395      zmet = 1._wp 
     396 
    359397      ll_do_bclinic = .TRUE. 
    360398      ncall = 1 
     
    362400      IF( PRESENT(kcall) ) THEN 
    363401         ! comment line below if tilda coordinate has to be computed at each call 
    364          ! This is not working yet 
    365 !         IF (kcall == 2 .AND. (ln_vvl_ztilde.OR.ln_vvl_layer) ) ll_do_bclinic = .FALSE. 
     402         IF (kcall == 2 .AND. (ln_vvl_ztilde.OR.ln_vvl_layer) ) ll_do_bclinic = .FALSE. 
    366403         ncall = kcall   
    367404      ENDIF 
     
    386423      END DO 
    387424      ! 
     425      IF ((ln_vvl_ztilde.OR.ln_vvl_layer).AND.(zmet==1._wp)) THEN 
     426         DO jk = 1, jpkm1 
     427            e3t_a(:,:,jk) = e3t_a(:,:,jk) - tilde_e3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
     428         END DO   
     429      ENDIF 
     430 
    388431      IF( (ln_vvl_ztilde.OR.ln_vvl_layer) .AND. ll_do_bclinic ) THEN   ! z_tilde or layer coordinate ! 
    389          !                                                               ! ------baroclinic part------ ! 
     432         !                                                             ! ------baroclinic part------ ! 
    390433         tilde_e3t_a(:,:,:) = 0.0_wp  ! tilde_e3t_a used to store tendency terms 
    391434         un_td(:,:,:) = 0.0_wp        ! Transport corrections 
     
    415458                  ztu(:,:,jk) = un(:,:,jk) * e3u_n(:,:,jk) * e2u(:,:) 
    416459                  ztv(:,:,jk) = vn(:,:,jk) * e3v_n(:,:,jk) * e1v(:,:) 
    417                   ze3t(:,:,jk) = -tilde_e3t_a(:,:,jk) - e3t_n(:,:,jk) * zhdiv(:,:) 
     460                  ze3t(:,:,jk) = -tilde_e3t_a(:,:,jk) - (e3t_n(:,:,jk)-zmet*tilde_e3t_n(:,:,jk)) * zhdiv(:,:) 
    418461               END DO 
    419462               ! 
    420                un_lf(:,:,:)  =    un_lf(:,:,:) * (1._wp - zalpha) + zalpha * ztu(:,:,:) 
    421                vn_lf(:,:,:)  =    vn_lf(:,:,:) * (1._wp - zalpha) + zalpha * ztv(:,:,:) 
    422             hdivn_lf(:,:,:)  = hdivn_lf(:,:,:) * (1._wp - zalpha) + zalpha * ze3t(:,:,:) 
    423 ! 2nd order filtering: 
    424 !               un_lf2(:,:,:)  =    un_lf2(:,:,:) * (1._wp - zalpha) + zalpha * ztu(:,:,:) 
    425 !               vn_lf2(:,:,:)  =    vn_lf2(:,:,:) * (1._wp - zalpha) + zalpha * ztv(:,:,:) 
    426 !            hdivn_lf2(:,:,:)  = hdivn_lf2(:,:,:) * (1._wp - zalpha) + zalpha * ze3t(:,:,:) 
    427 !               un_lf(:,:,:)   =     un_lf(:,:,:) * (1._wp - zalpha) + zalpha * un_lf2(:,:,:) 
    428 !               vn_lf(:,:,:)   =     vn_lf(:,:,:) * (1._wp - zalpha) + zalpha * vn_lf2(:,:,:) 
    429 !            hdivn_lf(:,:,:)   =  hdivn_lf(:,:,:) * (1._wp - zalpha) + zalpha * hdivn_lf2(:,:,:) 
     463                  un_lf(:,:,:,nn_filt_order)  =    un_lf(:,:,:,nn_filt_order) * (1._wp - zalpha) + zalpha * ztu(:,:,:) 
     464                  vn_lf(:,:,:,nn_filt_order)  =    vn_lf(:,:,:,nn_filt_order) * (1._wp - zalpha) + zalpha * ztv(:,:,:) 
     465               hdivn_lf(:,:,:,nn_filt_order)  = hdivn_lf(:,:,:,nn_filt_order) * (1._wp - zalpha) + zalpha * ze3t(:,:,:) 
     466 
     467               DO jo = nn_filt_order-1,1,-1 
     468                     un_lf(:,:,:,jo)  =    un_lf(:,:,:,jo) * (1._wp - zalpha) + zalpha *    un_lf(:,:,:,jo+1) 
     469                     vn_lf(:,:,:,jo)  =    vn_lf(:,:,:,jo) * (1._wp - zalpha) + zalpha *    vn_lf(:,:,:,jo+1) 
     470                  hdivn_lf(:,:,:,jo)  = hdivn_lf(:,:,:,jo) * (1._wp - zalpha) + zalpha * hdivn_lf(:,:,:,jo+1) 
     471               END DO 
    430472            ENDIF 
    431473            ! 
    432474            DO jk = 1, jpkm1 
    433                ztu(:,:,jk) = (un(:,:,jk)-un_lf(:,:,jk)/e3u_n(:,:,jk)*r1_e2u(:,:))*umask(:,:,jk) 
    434                ztv(:,:,jk) = (vn(:,:,jk)-vn_lf(:,:,jk)/e3v_n(:,:,jk)*r1_e1v(:,:))*vmask(:,:,jk) 
    435                tilde_e3t_a(:,:,jk) =  tilde_e3t_a(:,:,jk) + hdivn_lf(:,:,jk) - frq_rst_e3t(:,:) * tilde_e3t_b(:,:,jk) 
    436             END DO 
    437  
     475               ztu(:,:,jk) = (un(:,:,jk)-un_lf(:,:,jk,1)/e3u_n(:,:,jk)*r1_e2u(:,:))*umask(:,:,jk) 
     476               ztv(:,:,jk) = (vn(:,:,jk)-vn_lf(:,:,jk,1)/e3v_n(:,:,jk)*r1_e1v(:,:))*vmask(:,:,jk) 
     477               tilde_e3t_a(:,:,jk) =  tilde_e3t_a(:,:,jk) + hdivn_lf(:,:,jk,1) - frq_rst_e3t(:,:) * tilde_e3t_b(:,:,jk) 
     478            END DO 
    438479            ! 
    439480         ELSEIF ( ln_vvl_layer ) THEN 
     
    447488         ! 
    448489         ! Block fluxes through small layers: 
    449          DO jk=1,jpkm1 
    450             DO ji = 1, jpi 
    451                DO jj= 1, jpj 
    452                   zmsku = 0.5_wp * ( 1._wp + SIGN(1._wp, e3u_n(ji,jj,jk) - hsmall) ) 
    453                   un_td(ji,jj,jk) = un_td(ji,jj,jk) - (1. - zmsku) * un(ji,jj,jk) * e3u_n(ji,jj,jk) * e2u(ji,jj) 
    454                   ztu(ji,jj,jk) = zmsku * ztu(ji,jj,jk) 
    455                   IF ( ln_vvl_ztilde ) un_lf(ji,jj,jk) = zmsku * un_lf(ji,jj,jk) 
    456                   ! 
    457                   zmskv = 0.5_wp * ( 1._wp + SIGN(1._wp, e3v_n(ji,jj,jk) - hsmall) ) 
    458                   vn_td(ji,jj,jk) = vn_td(ji,jj,jk) - (1. - zmskv) * vn(ji,jj,jk) * e3v_n(ji,jj,jk) * e1v(ji,jj) 
    459                   ztv(ji,jj,jk) = zmskv * ztv(ji,jj,jk) 
    460                   IF ( ln_vvl_ztilde ) vn_lf(ji,jj,jk) = zmskv * vn_lf(ji,jj,jk) 
    461                END DO 
    462             END DO 
    463          END DO       
     490!         DO jk=1,jpkm1 
     491!            DO ji = 1, jpi 
     492!               DO jj= 1, jpj 
     493!                  zmsku = 0.5_wp * ( 1._wp + SIGN(1._wp, e3u_n(ji,jj,jk) - hsmall) ) 
     494!                  un_td(ji,jj,jk) = un_td(ji,jj,jk) - (1. - zmsku) * un(ji,jj,jk) * e3u_n(ji,jj,jk) * e2u(ji,jj) 
     495!                  ztu(ji,jj,jk) = zmsku * ztu(ji,jj,jk) 
     496!                  IF ( ln_vvl_ztilde ) un_lf(ji,jj,jk) = zmsku * un_lf(ji,jj,jk) 
     497!                  ! 
     498!                  zmskv = 0.5_wp * ( 1._wp + SIGN(1._wp, e3v_n(ji,jj,jk) - hsmall) ) 
     499!                  vn_td(ji,jj,jk) = vn_td(ji,jj,jk) - (1. - zmskv) * vn(ji,jj,jk) * e3v_n(ji,jj,jk) * e1v(ji,jj) 
     500!                  ztv(ji,jj,jk) = zmskv * ztv(ji,jj,jk) 
     501!                  IF ( ln_vvl_ztilde ) vn_lf(ji,jj,jk) = zmskv * vn_lf(ji,jj,jk) 
     502!               END DO 
     503!            END DO 
     504!         END DO       
    464505         ! 
    465506         ! Do advection 
     
    481522         !  
    482523         ! Thickness anomaly diffusion: 
    483          ! ----------------------------- 
     524         ! ---------------------------- 
    484525         ztu(:,:,:) = 0.0_wp 
    485526         ztv(:,:,:) = 0.0_wp 
     
    498539               DO jj = 1, jpjm1                 ! First derivative (gradient) 
    499540                  DO ji = 1, fs_jpim1   ! vector opt.                   
    500 !                     ztu(ji,jj,jk) =  umask(ji,jj,jk) * e2_e1u(ji,jj) & 
    501 !                                     &  * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
    502 !                     ztv(ji,jj,jk) =  vmask(ji,jj,jk) * e1_e2v(ji,jj) &  
    503 !                                     &  * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
    504541                     ztu(ji,jj,jk) =  umask(ji,jj,jk) * e2_e1u(ji,jj) & 
    505542                                     &  * ( ze3t(ji,jj,jk) - ze3t(ji+1,jj  ,jk) ) 
     
    542579               DO jj = 1, jpjm1 
    543580                  DO ji = 1, fs_jpim1   ! vector opt.                   
    544 !                     zdu =  rn_ahe3_lap * umask(ji,jj,jk) * e2_e1u(ji,jj) & 
    545 !                         &  * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
    546 !                     zdv =  rn_ahe3_lap * vmask(ji,jj,jk) * e1_e2v(ji,jj) &  
    547 !                         &  * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
    548581                     zdu =  rn_ahe3_lap * umask(ji,jj,jk) * e2_e1u(ji,jj) & 
    549582                         &  * ( ze3t(ji,jj,jk) - ze3t(ji+1,jj  ,jk) ) 
     
    561594            DO jj = 2, jpjm1 
    562595               DO ji = fs_2, fs_jpim1   ! vector opt. 
    563 !                  tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   ztu(ji-1,jj  ,jk) - ztu(ji,jj,jk)    & 
    564 !                     &                                          +     ztv(ji  ,jj-1,jk) - ztv(ji,jj,jk)    & 
    565 !                     &                                            ) * r1_e1e2t(ji,jj) 
    566596                  un_td(ji,jj,jk) = un_td(ji,jj,jk) + ztu(ji,jj,jk+1) - ztu(ji,jj,jk  )  
    567597                  vn_td(ji,jj,jk) = vn_td(ji,jj,jk) + ztv(ji,jj,jk+1) - ztv(ji,jj,jk  )  
     
    575605         END DO 
    576606 
    577   
    578 !         un_td(:,:,:) = un_td(:,:,:) + ztu(:,:,:) 
    579 !         vn_td(:,:,:) = vn_td(:,:,:) + ztv(:,:,:) 
    580  
    581607         CALL lbc_lnk_multi( un_td, 'U', -1., vn_td, 'V', -1. )     !* local domain boundaries 
    582608         ! 
     
    591617         ! Remove "external thickness" tendency: 
    592618         DO jk = 1, jpkm1 
    593             tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) +   e3t_n(:,:,jk) * zhdiv(:,:) 
     619            tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) +   (e3t_n(:,:,jk)-zmet*tilde_e3t_n(:,:,jk)) * zhdiv(:,:) * tmask(:,:,jk) 
    594620         END DO  
    595621                    
     
    646672         END DO 
    647673         CALL lbc_lnk( ze3t(:,:,:), 'T', 1. ) 
    648          hdivn_lf(:,:,:)  =  hdivn_lf(:,:,:)  + zalpha * ze3t(:,:,:)  
    649 ! 2nd order filtering: 
    650 !         hdivn_lf2(:,:,:)  = hdivn_lf2(:,:,:) + zalpha * ze3t(:,:,:) 
    651 !         hdivn_lf(:,:,:)  =  hdivn_lf(:,:,:)  + zalpha * zalpha * ze3t(:,:,:) 
     674 
     675         DO jo = nn_filt_order,1,-1 
     676            hdivn_lf(:,:,:,jo) =  hdivn_lf(:,:,:,jo)  + zalpha**(nn_filt_order-jo+1) * ze3t(:,:,:) 
     677         END DO 
    652678      ENDIF 
    653679 
    654680      IF( ln_vvl_ztilde .OR. ln_vvl_layer )  THEN   ! z_tilde or layer coordinate ! 
    655681      !                                             ! ---baroclinic part--------- ! 
     682 
     683         IF ( (ncall==2).AND.(.NOT.ll_do_bclinic) ) THEN 
     684 
     685            DO jk = 1, jpkm1 
     686               ztu(:,:,jk) = (un_adv(:,:)*r1_hu_n(:,:) - un_b(:,:) ) * e3u_n(:,:,jk) * e2u(:,:) * umask(:,:,jk) 
     687               ztv(:,:,jk) = (vn_adv(:,:)*r1_hv_n(:,:) - vn_b(:,:) ) * e3v_n(:,:,jk) * e1v(:,:) * vmask(:,:,jk) 
     688               DO jj = 2, jpjm1 
     689                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     690                     ze3t(ji,jj,jk) = (  ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk) & 
     691                                    &  + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk) & 
     692                                    & ) / ( e1t(ji,jj) * e2t(ji,jj) ) 
     693                  END DO 
     694               END DO 
     695            END DO 
     696            ! 
     697            zhdiv(:,:) = 0. 
     698            DO jk = 1, jpkm1 
     699               zhdiv(:,:) = zhdiv(:,:) + ze3t(:,:,jk) * tmask(:,:,jk) 
     700            END DO 
     701            zhdiv(:,:) = zhdiv(:,:) / ( ht_0(:,:) + sshn(:,:) + 1._wp - tmask_i(:,:) )  
     702            ! 
     703            DO jk = 1, jpkm1 
     704               tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - z2dt * (ze3t(:,:,jk) &  
     705                                   & - zhdiv(:,:)*(e3t_n(:,:,jk)-zmet*tilde_e3t_n(:,:,jk))*tmask(:,:,jk)) 
     706            END DO 
     707            CALL lbc_lnk( tilde_e3t_a(:,:,:), 'T', 1. ) 
     708         ENDIF 
     709         ! 
    656710         DO jk = 1, jpkm1 
    657             e3t_a(:,:,jk) = e3t_a(:,:,jk) + (tilde_e3t_a(:,:,jk) - tilde_e3t_b(:,:,jk))                      
    658          END DO 
     711            e3t_a(:,:,jk) = e3t_a(:,:,jk) + tilde_e3t_a(:,:,jk) - tilde_e3t_b(:,:,jk) 
     712         END DO 
     713!         IF( ln_vvl_ztilde ) THEN !  Relax barotropic component: 
     714!            DO jk = 1, jpkm1 
     715!               e3t_a(:,:,jk) = e3t_a(:,:,jk)  & 
     716!                             & - z2dt * frq_rst_e3t(:,:) * (e3t_b(:,:,jk) - tilde_e3t_b(:,:,jk)   &  
     717!                             & - e3t_0(:,:,jk) * (ht_0(:,:) + sshb(:,:))/ (ht_0(:,:)*tmask(:,:,1) + 1._wp - tmask(:,:,1)))                  
     718!            END DO 
     719!         ENDIF 
    659720      ENDIF 
    660721 
     
    686747            ENDIF 
    687748            IF (lwp) THEN 
     749               ji = ijk_min(1) ; jj = ijk_min(2) ; jk = ijk_min(3) 
    688750               WRITE(numout, *) 'Negative scale factor, e3t_n =', z_tmin 
    689751               WRITE(numout, *) 'at i, j, k=', ijk_min             
     
    12701332      ! 
    12711333      INTEGER ::   ji, jj, jk 
    1272       INTEGER ::   id1, id2, id3, id4, id5     ! local integers 
     1334      INTEGER ::   id1, id2, id3, id4, id5, id6, id7     ! local integers 
    12731335      !!---------------------------------------------------------------------- 
    12741336      ! 
     
    12841346            id4 = iom_varid( numror, 'tilde_e3t_n', ldstop = .FALSE. ) 
    12851347            id5 = iom_varid( numror, 'hdivn_lf', ldstop = .FALSE. ) 
     1348            id6 = iom_varid( numror, 'un_lf', ldstop = .FALSE. ) 
     1349            id7 = iom_varid( numror, 'vn_lf', ldstop = .FALSE. ) 
    12861350            !                             ! --------- ! 
    12871351            !                             ! all cases ! 
     
    13441408               IF( ln_vvl_ztilde ) THEN   ! z_tilde case ! 
    13451409                  !                       ! ------------ ! 
    1346                   IF( id5 > 0 ) THEN  ! required array exists 
    1347                      CALL iom_get( numror, jpdom_autoglo, 'hdivn_lf', hdivn_lf(:,:,:), ldxios = lrxios ) 
     1410                  IF( MIN(id5, id6, id7) > 0 ) THEN  ! required arrays exist 
     1411                     CALL iom_get( numror, jpdom_autoglo, 'hdivn_lf', hdivn_lf(:,:,:,1), ldxios = lrxios ) 
     1412                     CALL iom_get( numror, jpdom_autoglo, 'un_lf', un_lf(:,:,:,1), ldxios = lrxios ) 
     1413                     CALL iom_get( numror, jpdom_autoglo, 'vn_lf', vn_lf(:,:,:,1), ldxios = lrxios ) 
    13481414                  ELSE                ! array is missing 
    1349                      hdivn_lf(:,:,:) = 0.0_wp 
     1415                     hdivn_lf(:,:,:,:) = 0.0_wp 
     1416                     un_lf(:,:,:,:) = 0.0_wp 
     1417                     vn_lf(:,:,:,:) = 0.0_wp 
    13501418                  ENDIF 
    13511419               ENDIF 
     
    14181486               tilde_e3t_n(:,:,:) = 0._wp 
    14191487               IF( ln_vvl_ztilde ) THEN  
    1420                   hdivn_lf(:,:,:) = 0._wp  
    1421                   un_lf(:,:,:) = 0._wp  
    1422                   vn_lf(:,:,:) = 0._wp  
     1488                  hdivn_lf(:,:,:,:) = 0._wp  
     1489                  un_lf(:,:,:,:) = 0._wp  
     1490                  vn_lf(:,:,:,:) = 0._wp  
    14231491               ENDIF 
    14241492            END IF 
     
    14431511         IF( ln_vvl_ztilde ) THEN                    ! z_tilde case ! 
    14441512            !                                        ! ------------ ! 
    1445             CALL iom_rstput( kt, nitrst, numrow, 'hdivn_lf', hdivn_lf(:,:,:), ldxios = lwxios) 
     1513            CALL iom_rstput( kt, nitrst, numrow, 'hdivn_lf', hdivn_lf(:,:,:,1), ldxios = lwxios) 
     1514            CALL iom_rstput( kt, nitrst, numrow, 'un_lf', un_lf(:,:,:,1), ldxios = lwxios) 
     1515            CALL iom_rstput( kt, nitrst, numrow, 'vn_lf', vn_lf(:,:,:,1), ldxios = lwxios) 
    14461516         ENDIF 
    14471517         ! 
     
    16041674      ll_chk_bot2top   = .TRUE. 
    16051675      ll_chk_top2bot   = .TRUE. 
    1606       dzmin_int  = 0.1_wp   ! Absolute minimum depth in the interior (in meters) 
     1676      dzmin_int  = 1.0_wp   ! Absolute minimum depth in the interior (in meters) 
    16071677      dzmin_surf = 1.0_wp   ! Absolute minimum depth at the surface (in meters) 
    1608       zfrch_stp  = 0.5_wp   ! Maximum fractionnal thickness change in one time step (<= 1.) 
     1678      zfrch_stp  = 5._wp   ! Maximum fractionnal thickness change in one time step (<= 1.) 
    16091679      zfrch_rel  = 0.4_wp   ! Maximum relative thickness change in the vertical (<= 1.) 
    16101680      zfrac_bot  = 0.05_wp  ! Fraction of bottom level allowed to change 
     
    16141684         zvlim   = 0.5_wp   ! max d2h/dh 
    16151685      ll_lapdiff_cond  = .TRUE.  ! Conditionnal Laplacian diffusion of interfaces 
    1616          zhdiff  = 0.01_wp   ! ad. 
     1686         zhdiff  = 0.01_wp  ! ad. 
    16171687         zhlim   = 0.03_wp  ! ad. max lap(z)*e1 
    1618       ll_blpdiff_cond  = .FALSE.  ! Conditionnal Bilaplacian diffusion of interfaces 
     1688      ll_blpdiff_cond  = .TRUE.  ! Conditionnal Bilaplacian diffusion of interfaces 
    16191689         zhdiff2 = 0.2_wp   ! ad. 
    1620 !         zhlim2  = 3.e-11_wp  !  max bilap(z) 
    16211690         zhlim2  = 0.01_wp  ! ad. max bilap(z)*e1**3 
    16221691      ! ---------------------------------------------------------------------------------------  
     
    18191888                  zh_old = zwdw(ji,jj,jk+1) - zwdw(ji,jj,jk) 
    18201889                  zh_min = MIN(zh_0/3._wp, dzmin_int) 
     1890!                  zh_min = MAX(zh_min, zh_min-e3t_a(ji,jj,jk)+e3t_0(ji,jj,jk)) 
    18211891                  !  
    18221892                  ! Set maximum and minimum vertical excursions    
     
    18361906                  ! Ensure minimum layer thickness:                   
    18371907!                  zh_new = MAX((1._wp-zfrch_stp)*zh_bef, zh_new) 
    1838                   zh_new = cush(zh_new, zh_min) 
     1908!                  zh_new = cush(zh_new, zh_min) 
     1909                  zh_new = MAX(zh_new, zh_min)            
    18391910                  ! 
    18401911                  ! Final flux: 
     
    18421913                  ! 
    18431914                  ! Limit thickness change in 1 time step:  
    1844                   ztmp = MIN( ABS(zdiff), zfrch_stp*zh_bef ) 
    1845                   zdiff = SIGN(ztmp, zh_new - zh_old) 
     1915!                  ztmp = MIN( ABS(zdiff), zfrch_stp*zh_bef ) 
     1916!                  zdiff = SIGN(ztmp, zh_new - zh_old) 
    18461917                  zh_new = zdiff + zh_old 
    18471918                  ! 
     
    18641935                  zh_old = zwdw(ji,jj,jk+1) - zwdw(ji,jj,jk) 
    18651936                  zh_min = MIN(zh_0/3._wp, dzmin_int) 
    1866                   ! 
    1867                   zwdw(ji,jj,jk+1) = MAX(zwdw(ji,jj,jk+1), REAL(jk)*dzmin_surf) 
     1937!                  zh_min = MAX(zh_min, zh_min-e3t_a(ji,jj,jk)+e3t_0(ji,jj,jk)) 
     1938                  ! 
     1939!                  zwdw(ji,jj,jk+1) = MAX(zwdw(ji,jj,jk+1), REAL(jk)*dzmin_surf) 
    18681940                  ! 
    18691941                  ! New layer thickness: 
     
    18721944                  ! Ensure minimum layer thickness: 
    18731945!                  zh_new = MAX((1._wp-zfrch_stp)*zh_bef, zh_new) 
    1874                   zh_new = cush(zh_new, zh_min) 
     1946!                  zh_new = cush(zh_new, zh_min) 
     1947                  zh_new = MAX(zh_new, zh_min)    
    18751948                  ! 
    18761949                  ! Final flux: 
     
    20892162      ENDIF 
    20902163 
    2091       IF ( ln_vvl_dbg ) THEN 
    2092          zmin = MINVAL( zwi(:,:,:), mask = tmask(:,:,:) == 1.e0 )  
    2093          IF( lk_mpp )   CALL mpp_min( zmin ) 
    2094          IF( zmin < 0._wp) THEN 
    2095             IF(lwp) CALL ctl_warn('vvl_adv: CFL issue here') 
    2096             IF(lwp) WRITE(numout,*) zmin 
    2097          ENDIF 
    2098       ENDIF 
     2164!      IF ( ln_vvl_dbg ) THEN 
     2165!         zmin = MINVAL( zwi(:,:,:), mask = tmask(:,:,:) == 1.e0 )  
     2166!         IF( lk_mpp )   CALL mpp_min( zmin ) 
     2167!         IF( zmin < 0._wp) THEN 
     2168!            IF(lwp) CALL ctl_warn('vvl_adv: CFL issue here') 
     2169!            IF(lwp) WRITE(numout,*) zmin 
     2170!         ENDIF 
     2171!      ENDIF 
    20992172 
    21002173      ! 3. antidiffusive flux : high order minus low order 
     
    21382211               ! add them to the general tracer trends 
    21392212               pta(ji,jj,jk) = pta(ji,jj,jk) + ztra 
     2213               zwi(ji,jj,jk) = (e3t_b(ji,jj,jk) + z2dtt * pta(ji,jj,jk)* bdytmask(ji,jj) ) * tmask(ji,jj,jk) 
    21402214            END DO 
    21412215         END DO 
    21422216      END DO 
     2217 
     2218      IF ( ln_vvl_dbg ) THEN 
     2219         zmin = MINVAL( zwi(:,:,:), mask = tmask(:,:,:) == 1.e0 )  
     2220         IF( lk_mpp )   CALL mpp_min( zmin ) 
     2221         IF( zmin < 0._wp) THEN 
     2222            IF(lwp) CALL ctl_warn('vvl_adv: CFL issue here') 
     2223            IF(lwp) WRITE(numout,*) zmin 
     2224         ENDIF 
     2225      ENDIF 
    21432226      ! 
    21442227      IF( ln_timing )  CALL timing_stop('dom_vvl_adv_fct') 
  • NEMO/branches/2018/dev_r10057_ENHANCE03_ZTILDE/tests/COMODO-IW/EXPREF/namelist_cfg

    r10128 r10164  
    237237      rn_ahe3_lap = 30.e0 
    238238   ln_vvl_blp = .true.              !  Bilaplacian thickness diffusion 
    239       rn_ahe3_blp = -2.e8 
     239      rn_ahe3_blp = -1.e8 
    240240   rn_rst_e3t    = 30.e0            !  ztilde to zstar restoration timescale [days] 
    241241   rn_lf_cutoff  = 5.0e0            !  cutoff frequency for low-pass filter  [days] 
Note: See TracChangeset for help on using the changeset viewer.