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

Changeset 13002 for NEMO


Ignore:
Timestamp:
2020-06-02T09:03:47+02:00 (4 years ago)
Author:
davestorkey
Message:

UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff : Add two versions of a ramp in addition to the abrupt cutoff on the RR/GS ratio.

Location:
NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff/cfgs/SHARED/field_def_nemo-oce.xml

    r12867 r13002  
    7070        <field id="bn2"          long_name="squared Brunt-Vaisala frequency"                                           unit="s-1"    grid_ref="grid_T_3D" /> 
    7171        <field id="rhop"         long_name="potential density (sigma0)"        standard_name="sea_water_sigma_theta"   unit="kg/m3"  grid_ref="grid_T_3D" /> 
     72        <field id="N_2d"         long_name="Depth-mean N"                                                              unit="m/s" /> 
    7273        <field id="RossRad"      long_name="Rossby radius"                                                             unit="m" /> 
     74        <field id="RR_GS"        long_name="Rossby radius / min(dx,dy)"                                                unit="1" /> 
     75        <field id="GS_RR"        long_name="min(dx,dy) / Rossby radius"                                                unit="1" /> 
    7376 
    7477        <!-- Energy - horizontal divergence --> 
  • NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff/cfgs/SHARED/namelist_ref

    r11715 r13002  
    828828      !                             !   = 30 F(i,j,k)  =ldf_c2d * ldf_c1d 
    829829      !                        !  time invariant coefficients:  aei0 = 1/2  Ue*Le  
    830       rn_Ue        = 0.02           !  lateral diffusive velocity [m/s] (nn_aht_ijk_t= 0, 10, 20, 30) 
    831       rn_Le        = 200.e+3        !  lateral diffusive length   [m]   (nn_aht_ijk_t= 0, 10) 
     830      rn_Ue        = 0.02           !  lateral diffusive velocity [m/s] (nn_aei_ijk_t= 0, 10, 20, 30) 
     831      rn_Le        = 200.e+3        !  lateral diffusive length   [m]   (nn_aei_ijk_t= 0, 10) 
    832832      ! 
    833833      ln_ldfeiv_dia =.false.   ! diagnose eiv stream function and velocities 
     834      nn_ldfeiv_shape = 0           !  shape of bounding coefficient    (nn_aei_ijk_t= 21 only) 
    834835/ 
    835836!----------------------------------------------------------------------- 
  • NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff/src/OCE/LDF/ldftra.F90

    r12861 r13002  
    7272   REAL(wp), PUBLIC ::      rn_Ue               !: lateral diffusive velocity  [m/s] 
    7373   REAL(wp), PUBLIC ::      rn_Le               !: lateral diffusive length    [m] 
     74   INTEGER,  PUBLIC ::   nn_ldfeiv_shape     !: shape of bounding coefficient (Treguier et al formulation only) 
    7475    
    7576   !                                  ! Flag to control the type of lateral diffusive operator 
     
    503504      !! 
    504505      NAMELIST/namtra_eiv/ ln_ldfeiv   , ln_ldfeiv_dia,   &   ! eddy induced velocity (eiv) 
    505          &                 nn_aei_ijk_t, rn_Ue, rn_Le         ! eiv  coefficient 
     506         &                 nn_aei_ijk_t, rn_Ue, rn_Le,    &   ! eiv  coefficient 
     507         &                 nn_ldfeiv_shape 
    506508      !!---------------------------------------------------------------------- 
    507509      ! 
     
    597599            IF(lwp) WRITE(numout,*) '                                       = F( growth rate of baroclinic instability )' 
    598600            IF(lwp) WRITE(numout,*) '           maximum allowed value: aei0 = ', aei0, ' m2/s' 
     601            IF(lwp) WRITE(numout,*) '           shape of bounding coefficient : ',nn_ldfeiv_shape 
    599602            ! 
    600603            l_ldfeiv_time = .TRUE.     ! will be calculated by call to ldf_tra routine in step.F90 
     
    654657      ! 
    655658      INTEGER  ::   ji, jj, jk    ! dummy loop indices 
    656       REAL(wp) ::   zfw, ze3w, zn2, z1_f20, zaht, zaht_min, zzaei    ! local scalars 
    657       REAL(wp), DIMENSION(jpi,jpj) ::   zn, zah, zhw, zRo, zRo_lim, zaeiw   ! 2D workspace 
     659      REAL(wp) ::   zfw, ze3w, zn2, z1_f20, zaht, zaht_min, zzaei, z2_3    ! local scalars 
     660      REAL(wp), DIMENSION(jpi,jpj) ::   zn, zah, zhw, zRo, zRo_lim, zaeiw, zratio   ! 2D workspace 
    658661      !!---------------------------------------------------------------------- 
    659662      ! 
     
    713716         END DO 
    714717      END DO 
     718      IF( iom_use('N_2d') ) CALL iom_put('N_2d',zn(:,:)/ht_0(:,:)) 
    715719      CALL iom_put('RossRad',zRo) 
    716720      !                                         !==  Bound on eiv coeff.  ==! 
    717721      z1_f20 = 1._wp / (  2._wp * omega * sin( rad * 20._wp )  ) 
     722      z2_3 = 2._wp/3._wp 
    718723      DO jj = 2, jpjm1 
    719724         DO ji = fs_2, fs_jpim1   ! vector opt. 
    720725            zzaei = MIN( 1._wp, ABS( ff_t(ji,jj) * z1_f20 ) ) * zaeiw(ji,jj)     ! tropical decrease 
     726            zaeiw(ji,jj) = MIN( zzaei, paei0 ) 
     727         END DO 
     728      END DO 
     729 
     730      !! temporary fudge for diagnostic purposes: 
     731      !! zaeiw(:,:) = paei0  
     732 
     733      SELECT CASE(nn_ldfeiv_shape) 
     734         CASE(1) !! Abrupt cut-off on Rossby radius: 
    721735! JD : modifications here to introduce scaling by local rossby radius of deformation vs local grid scale 
    722736! arbitrary decision that GM is de-activated if local rossy radius larger than 2 times local grid scale 
    723737! based on Hallberg (2013) 
    724             IF ( l_Ro_cutoff .AND. zRo(ji,jj) >= ( 2._wp * MIN( e1t(ji,jj), e2t(ji,jj) ) ) ) THEN 
     738            DO jj = 2, jpjm1 
     739               DO ji = fs_2, fs_jpim1   ! vector opt. 
     740                  IF ( l_Ro_cutoff .AND. zRo(ji,jj) >= ( 2._wp * MIN( e1t(ji,jj), e2t(ji,jj) ) ) ) THEN 
    725741! TODO : use a version of zRo that integrates over a few time steps ? 
    726                 zaeiw(ji,jj) = 0._wp 
    727             ELSE 
    728                 zaeiw(ji,jj) = MIN( zzaei, paei0 ) 
    729             ENDIF 
    730          END DO 
    731       END DO 
     742                      zaeiw(ji,jj) = 0._wp 
     743                  ELSE 
     744                      zaeiw(ji,jj) = MIN( zaeiw(ji,jj), paei0 ) 
     745                  ENDIF 
     746               END DO 
     747            END DO 
     748 
     749         CASE(2) !! Rossby radius ramp type 1: 
     750            DO jj = 2, jpjm1 
     751               DO ji = fs_2, fs_jpim1   ! vector opt. 
     752                  zratio(ji,jj) = zRo(ji,jj)/MIN(e1t(ji,jj),e2t(ji,jj)) 
     753                  zaeiw(ji,jj) = MIN( zaeiw(ji,jj), MAX( 0._wp, MIN( 1._wp, z2_3*(2._wp - zratio(ji,jj)) ) ) * paei0 ) 
     754               END DO 
     755            END DO 
     756            CALL iom_put('RR_GS',zratio) 
     757 
     758         CASE(3) !! Rossby radius ramp type 2: 
     759            DO jj = 2, jpjm1 
     760               DO ji = fs_2, fs_jpim1   ! vector opt. 
     761                  zratio(ji,jj) = MIN(e1t(ji,jj),e2t(ji,jj))/zRo(ji,jj) 
     762                  zaeiw(ji,jj) = MIN( zaeiw(ji,jj), MAX( 0._wp, MIN( 1._wp, z2_3*( zratio(ji,jj) - 0.5_wp ) ) ) * paei0 ) 
     763               END DO 
     764            END DO 
     765            CALL iom_put('GS_RR',zratio) 
     766 
     767      END SELECT 
     768 
    732769      CALL lbc_lnk( 'ldftra', zaeiw(:,:), 'W', 1. )       ! lateral boundary condition 
    733770      !                
Note: See TracChangeset for help on using the changeset viewer.