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 6102 for branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM – NEMO

Ignore:
Timestamp:
2015-12-17T16:48:54+01:00 (8 years ago)
Author:
timgraham
Message:

Reverted changes to TKE and GLS and commented out tke update in agrif user to avoid changing results in v3.6 stable.

Location:
branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM/NEMO
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM/NEMO/NST_SRC/agrif_user.F90

    r5081 r6102  
    316316   ! 
    317317# if defined key_zdftke 
    318    CALL Agrif_Update_tke(0) 
     318!   CALL Agrif_Update_tke(0) 
    319319# endif 
    320320   ! 
  • branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r6092 r6102  
    4343   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avmu , avmv    !: vertical viscosity coef at uw- & vw-pts       [m2/s] 
    4444   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avm  , avt     !: vertical viscosity & diffusivity coef at w-pt [m2/s] 
    45    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k , avm_k  ! not enhanced Kz 
    46    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k, avmv_k ! not enhanced Kz 
    47    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en             !: now turbulent kinetic energy   [m2/s2] 
    48  
     45  
    4946   !!---------------------------------------------------------------------- 
    5047   !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
     
    6259         &     avtb(jpk) , bfrva(jpi,jpj) , avtb_2d(jpi,jpj) ,      & 
    6360         &     tfrua(jpi, jpj), tfrva(jpi, jpj)              ,      & 
    64          &     avmu  (jpi,jpj,jpk), avm   (jpi,jpj,jpk),            & 
    65          &     avmv  (jpi,jpj,jpk), avt   (jpi,jpj,jpk),            & 
    66          &     avt_k (jpi,jpj,jpk), avm_k (jpi,jpj,jpk),            &  
    67          &     avmu_k(jpi,jpj,jpk), avmv_k(jpi,jpj,jpk),            &  
    68          &     en    (jpi,jpj,jpk), STAT = zdf_oce_alloc ) 
     61         &     avmu(jpi,jpj,jpk), avm(jpi,jpj,jpk)           ,      & 
     62         &     avmv(jpi,jpj,jpk), avt(jpi,jpj,jpk)           , STAT = zdf_oce_alloc ) 
    6963         ! 
    7064      IF( zdf_oce_alloc /= 0 )   CALL ctl_warn('zdf_oce_alloc: failed to allocate arrays') 
  • branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r6092 r6102  
    4242   LOGICAL , PUBLIC, PARAMETER ::   lk_zdfgls = .TRUE.   !: TKE vertical mixing flag 
    4343   ! 
     44   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en      !: now turbulent kinetic energy 
    4445   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mxln    !: now mixing length 
    4546   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   zwall   !: wall function 
     47   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k   ! not enhanced Kz 
     48   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avm_k   ! not enhanced Kz 
     49   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k  ! not enhanced Kz 
     50   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmv_k  ! not enhanced Kz 
    4651   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustars2 !: Squared surface velocity scale at T-points 
    4752   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustarb2 !: Squared bottom  velocity scale at T-points 
     
    115120      !!                ***  FUNCTION zdf_gls_alloc  *** 
    116121      !!---------------------------------------------------------------------- 
    117       ALLOCATE( mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) ,     & 
     122      ALLOCATE( en(jpi,jpj,jpk),  mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) ,     & 
     123         &      avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk),                    & 
     124         &      avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk),                    & 
    118125         &      ustars2(jpi,jpj), ustarb2(jpi,jpj)                      , STAT= zdf_gls_alloc ) 
    119126         ! 
     
    322329      !  
    323330      ! One level below 
    324       en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2))/zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 
     331      en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2)) & 
     332          &            / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 
    325333      en(:,:,2) = MAX(en(:,:,2), rn_emin ) 
    326334      z_elem_a(:,:,2) = 0._wp  
     
    343351      z_elem_a(:,:,2) = 0._wp 
    344352      zkar(:,:)       = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:)) )) 
    345       zflxs(:,:)      = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 
     353      zflxs(:,:)      = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) & 
     354           &                      * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 
    346355 
    347356      en(:,:,2) = en(:,:,2) + zflxs(:,:)/fse3w(:,:,2) 
  • branches/2015/dev_agrif_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r6099 r6102  
    5353   USE timing         ! Timing 
    5454   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    55 #if defined key_agrif 
    56    USE agrif_opa_interp 
    57    USE agrif_opa_update 
    58 #endif 
    5955 
    6056   IMPLICIT NONE 
     
    8985   REAL(wp) ::   rhftau_scl = 1.0_wp       ! scale factor applied to HF part of taum  (nn_etau=3) 
    9086 
     87   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en             !: now turbulent kinetic energy   [m2/s2] 
    9188   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   htau           ! depth of tke penetration (nn_htau) 
    9289   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dissl          ! now mixing lenght of dissipation 
    93    REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   apdlr          ! now mixing lenght of dissipation 
     90   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k , avm_k  ! not enhanced Kz 
     91   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k, avmv_k ! not enhanced Kz 
    9492#if defined key_c1d 
    9593   !                                                                        !!** 1D cfg only  **   ('key_c1d') 
     
    102100#  include "vectopt_loop_substitute.h90" 
    103101   !!---------------------------------------------------------------------- 
    104    !! NEMO/OPA 3.6 , NEMO Consortium (2011) 
     102   !! NEMO/OPA 4.0 , NEMO Consortium (2011) 
    105103   !! $Id$ 
    106104   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    117115         &      e_pdl(jpi,jpj,jpk) , e_ric(jpi,jpj,jpk) ,                          & 
    118116#endif 
    119          &      apdlr(jpi,jpj,jpk) , htau  (jpi,jpj)    , dissl(jpi,jpj,jpk) ,     &  
    120          &      STAT= zdf_tke_alloc      ) 
     117         &      en    (jpi,jpj,jpk) , htau  (jpi,jpj)    , dissl(jpi,jpj,jpk) ,     &  
     118         &      avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk),                          & 
     119         &      avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk), STAT= zdf_tke_alloc      ) 
    121120         ! 
    122121      IF( lk_mpp             )   CALL mpp_sum ( zdf_tke_alloc ) 
     
    174173      !!---------------------------------------------------------------------- 
    175174      ! 
    176 #if defined key_agrif  
    177       ! interpolation parent grid => child grid for avm_k ( ex : at west border: update column 1 and 2) 
    178       IF( .NOT.Agrif_Root() )   CALL Agrif_Tke 
    179 #endif 
    180       ! 
    181175      IF( kt /= nit000 ) THEN   ! restore before value to compute tke 
    182176         avt (:,:,:) = avt_k (:,:,:)  
     
    195189      avmv_k(:,:,:) = avmv(:,:,:)  
    196190      ! 
    197 #if defined key_agrif 
    198       ! Update child grid f => parent grid  
    199       IF( .NOT.Agrif_Root() )   CALL Agrif_Update_Tke( kt )      ! children only 
    200 #endif       
    201      !  
    202   END SUBROUTINE zdf_tke 
     191   END SUBROUTINE zdf_tke 
    203192 
    204193 
     
    234223      INTEGER , POINTER, DIMENSION(:,:  ) :: imlc 
    235224      REAL(wp), POINTER, DIMENSION(:,:  ) :: zhlc 
    236       REAL(wp), POINTER, DIMENSION(:,:,:) :: zpelc, zdiag, zd_up, zd_lw, z3du, z3dv 
    237       REAL(wp)                            ::   zri  !   local Richardson number 
     225      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpelc, zdiag, zd_up, zd_lw 
    238226      !!-------------------------------------------------------------------- 
    239227      ! 
     
    242230      CALL wrk_alloc( jpi,jpj, imlc )    ! integer 
    243231      CALL wrk_alloc( jpi,jpj, zhlc )  
    244       CALL wrk_alloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw, z3du, z3dv )  
     232      CALL wrk_alloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw )  
    245233      ! 
    246234      zbbrau = rn_ebb / rau0       ! Local constant initialisation 
     
    344332      ! 
    345333      DO jk = 2, jpkm1           !* Shear production at uw- and vw-points (energy conserving form) 
    346          DO jj = 1, jpjm1 
    347             DO ji = 1, fs_jpim1   ! vector opt. 
    348                z3du(ji,jj,jk) = 0.5 * ( avm(ji,jj,jk  ) + avm(ji+1,jj,jk) )   & 
    349                   &                 * (  un(ji,jj,jk-1) -  un(ji  ,jj,jk) )   & 
    350                   &                 * (  ub(ji,jj,jk-1) -  ub(ji  ,jj,jk) ) * wumask(ji,jj,jk) & 
    351                   &                 / (  fse3uw_n(ji,jj,jk) * fse3uw_b(ji,jj,jk) ) 
    352                z3dv(ji,jj,jk) = 0.5 * ( avm(ji,jj,jk  ) + avm(ji,jj+1,jk) )   & 
    353                   &                 * (  vn(ji,jj,jk-1) -  vn(ji,jj  ,jk) )   & 
    354                   &                 * (  vb(ji,jj,jk-1) -  vb(ji,jj  ,jk) ) * wvmask(ji,jj,jk) & 
    355                   &                 / (  fse3vw_n(ji,jj,jk) * fse3vw_b(ji,jj,jk) ) 
    356             END DO 
    357          END DO 
    358       END DO 
    359       ! 
    360       IF( nn_pdl == 1 ) THEN      !* Prandtl number case: compute apdlr 
    361          ! Note that zesh2 is also computed in the next loop. 
    362          ! We decided to compute it twice to keep code readability and avoid an IF case in the DO loops 
    363          DO jk = 2, jpkm1 
    364             DO jj = 2, jpjm1 
    365                DO ji = fs_2, fs_jpim1   ! vector opt. 
    366                   !                                          ! shear prod. at w-point weightened by mask 
    367                   zesh2  =  ( z3du(ji-1,jj,jk) + z3du(ji,jj,jk) ) / MAX( 1._wp , umask(ji-1,jj,jk) + umask(ji,jj,jk) )   & 
    368                      &    + ( z3dv(ji,jj-1,jk) + z3dv(ji,jj,jk) ) / MAX( 1._wp , vmask(ji,jj-1,jk) + vmask(ji,jj,jk) )     
    369                   !                                          ! local Richardson number 
    370                   zri   = MAX( rn2b(ji,jj,jk), 0._wp ) * avm(ji,jj,jk) / ( zesh2 + rn_bshear ) 
    371                   apdlr(ji,jj,jk) = MAX(  0.1_wp,  ri_cri / MAX( ri_cri , zri )  ) 
    372                    
    373                END DO 
    374             END DO 
    375          END DO 
    376          ! 
    377       ENDIF 
    378          !          
     334         DO jj = 1, jpj                 ! here avmu, avmv used as workspace 
     335            DO ji = 1, jpi 
     336               avmu(ji,jj,jk) = avmu(ji,jj,jk) * ( un(ji,jj,jk-1) - un(ji,jj,jk) )   & 
     337                  &                            * ( ub(ji,jj,jk-1) - ub(ji,jj,jk) )   &  
     338                  &                            / (  fse3uw_n(ji,jj,jk)               & 
     339                  &                              *  fse3uw_b(ji,jj,jk)  ) 
     340               avmv(ji,jj,jk) = avmv(ji,jj,jk) * ( vn(ji,jj,jk-1) - vn(ji,jj,jk) )   & 
     341                  &                            * ( vb(ji,jj,jk-1) - vb(ji,jj,jk) )   & 
     342                  &                            / (  fse3vw_n(ji,jj,jk)               & 
     343                  &                              *  fse3vw_b(ji,jj,jk)  ) 
     344            END DO 
     345         END DO 
     346      END DO 
     347      ! 
    379348      DO jk = 2, jpkm1           !* Matrix and right hand side in en 
    380349         DO jj = 2, jpjm1 
     
    385354               zzd_lw = zcof * ( avm  (ji,jj,jk  ) + avm  (ji,jj,jk-1) )   &  ! lower diagonal 
    386355                  &          / ( fse3t(ji,jj,jk-1) * fse3w(ji,jj,jk  ) ) 
    387                !                                   ! shear prod. at w-point weightened by mask 
    388                zesh2  =  ( z3du(ji-1,jj,jk) + z3du(ji,jj,jk) ) / MAX( 1._wp , umask(ji-1,jj,jk) + umask(ji,jj,jk) )   & 
    389                   &    + ( z3dv(ji,jj-1,jk) + z3dv(ji,jj,jk) ) / MAX( 1._wp , vmask(ji,jj-1,jk) + vmask(ji,jj,jk) )     
    390                ! 
     356                  !                                                           ! shear prod. at w-point weightened by mask 
     357               zesh2  =  ( avmu(ji-1,jj,jk) + avmu(ji,jj,jk) ) / MAX( 1._wp , umask(ji-1,jj,jk) + umask(ji,jj,jk) )   & 
     358                  &    + ( avmv(ji,jj-1,jk) + avmv(ji,jj,jk) ) / MAX( 1._wp , vmask(ji,jj-1,jk) + vmask(ji,jj,jk) )     
     359                  ! 
    391360               zd_up(ji,jj,jk) = zzd_up            ! Matrix (zdiag, zd_up, zd_lw) 
    392361               zd_lw(ji,jj,jk) = zzd_lw 
     
    486455      CALL wrk_dealloc( jpi,jpj, imlc )    ! integer 
    487456      CALL wrk_dealloc( jpi,jpj, zhlc )  
    488       CALL wrk_dealloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw, z3du, z3dv )  
     457      CALL wrk_dealloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw )  
    489458      ! 
    490459      IF( nn_timing == 1 )  CALL timing_stop('tke_tke') 
     
    530499      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    531500      REAL(wp) ::   zrn2, zraug, zcoef, zav     ! local scalars 
    532       REAL(wp) ::   zdku, zri, zsqen     !   -      - 
     501      REAL(wp) ::   zdku, zpdlr, zri, zsqen     !   -      - 
    533502      REAL(wp) ::   zdkv, zemxl, zemlm, zemlp   !   -      - 
    534503      REAL(wp), POINTER, DIMENSION(:,:,:) :: zmpdl, zmxlm, zmxld 
     
    691660            DO jj = 2, jpjm1 
    692661               DO ji = fs_2, fs_jpim1   ! vector opt. 
    693                   avt(ji,jj,jk)   = MAX( apdlr(ji,jj,jk) * avt(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * tmask(ji,jj,jk) 
     662                  zcoef = avm(ji,jj,jk) * 2._wp * fse3w(ji,jj,jk) * fse3w(ji,jj,jk) 
     663                  !                                          ! shear 
     664                  zdku = avmu(ji-1,jj,jk) * ( un(ji-1,jj,jk-1) - un(ji-1,jj,jk) ) * ( ub(ji-1,jj,jk-1) - ub(ji-1,jj,jk) )   & 
     665                    &  + avmu(ji  ,jj,jk) * ( un(ji  ,jj,jk-1) - un(ji  ,jj,jk) ) * ( ub(ji  ,jj,jk-1) - ub(ji  ,jj,jk) ) 
     666                  zdkv = avmv(ji,jj-1,jk) * ( vn(ji,jj-1,jk-1) - vn(ji,jj-1,jk) ) * ( vb(ji,jj-1,jk-1) - vb(ji,jj-1,jk) )   & 
     667                    &  + avmv(ji,jj  ,jk) * ( vn(ji,jj  ,jk-1) - vn(ji,jj  ,jk) ) * ( vb(ji,jj  ,jk-1) - vb(ji,jj  ,jk) ) 
     668                  !                                          ! local Richardson number 
     669                  zri   = MAX( rn2b(ji,jj,jk), 0._wp ) * zcoef / (zdku + zdkv + rn_bshear ) 
     670                  zpdlr = MAX(  0.1_wp,  0.2 / MAX( 0.2 , zri )  ) 
     671!!gm and even better with the use of the "true" ri_crit=0.22222...  (this change the results!) 
     672!!gm              zpdlr = MAX(  0.1_wp,  ri_crit / MAX( ri_crit , zri )  ) 
     673                  avt(ji,jj,jk)   = MAX( zpdlr * avt(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * wmask(ji,jj,jk) 
    694674# if defined key_c1d 
    695                   e_pdl(ji,jj,jk) = apdlr(ji,jj,jk) * tmask(ji,jj,jk)    ! c1d configuration : save masked Prandlt number 
    696                   e_ric(ji,jj,jk) = zri * tmask(ji,jj,jk)                            ! c1d config. : save Ri 
     675                  e_pdl(ji,jj,jk) = zpdlr * wmask(ji,jj,jk)  ! c1d configuration : save masked Prandlt number 
     676                  e_ric(ji,jj,jk) = zri   * wmask(ji,jj,jk)  ! c1d config. : save Ri 
    697677# endif 
    698678              END DO 
Note: See TracChangeset for help on using the changeset viewer.