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 13239 for NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff/src/OCE/LDF/ldftra.F90 – NEMO

Ignore:
Timestamp:
2020-07-03T12:30:58+02:00 (4 years ago)
Author:
davestorkey
Message:

UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff : Yet another shape option, and new slope diagnostic.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0.1_GM_rossby_radius_cutoff/src/OCE/LDF/ldftra.F90

    r13161 r13239  
    647647      REAL(wp) ::   zfw, ze3w, zn2, z1_f20, zaht, zaht_min, zzaei, z2_3    ! local scalars 
    648648      REAL(wp), DIMENSION(jpi,jpj) ::   zn, zah, zhw, zRo, zRo_lim, zTclinic_recip, zaeiw, zratio   ! 2D workspace 
     649      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmodslp ! 3D workspace 
    649650      !!---------------------------------------------------------------------- 
    650651      ! 
    651652      zn (:,:) = 0._wp        ! Local initialization 
     653      zmodslp(:,:,:) = 0._wp 
    652654      zhw(:,:) = 5._wp 
    653655      zah(:,:) = 0._wp 
     
    689691                  ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w)) 
    690692                  ze3w = e3w_n(ji,jj,jk) * tmask(ji,jj,jk) 
    691                   zah(ji,jj) = zah(ji,jj) + zn2 * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
    692                      &                            + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) * ze3w 
     693                  zmodslp(ji,jj,jk) =  wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
     694                     &               + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) 
     695                  zah(ji,jj) = zah(ji,jj) + zn2 * zmodslp(ji,jj,jk) * ze3w 
    693696                  zhw(ji,jj) = zhw(ji,jj) + ze3w 
    694697               END DO 
     
    709712      END DO 
    710713      IF( iom_use('N_2d') ) CALL iom_put('N_2d',zn(:,:)/zhw(:,:)) 
     714      IF( iom_use('modslp') ) CALL iom_put('modslp',SQRT(zmodslp(:,:,:)) ) 
    711715      CALL iom_put('RossRad',zRo) 
    712716      CALL iom_put('RossRadlim',zRo_lim) 
     
    715719      z1_f20 = 1._wp / (  2._wp * omega * sin( rad * 20._wp )  ) 
    716720      z2_3 = 2._wp/3._wp 
    717       DO jj = 2, jpjm1 
    718          DO ji = fs_2, fs_jpim1   ! vector opt. 
    719             zzaei = MIN( 1._wp, ABS( ff_t(ji,jj) * z1_f20 ) ) * zaeiw(ji,jj)     ! tropical decrease 
    720             zaeiw(ji,jj) = MIN( zzaei, paei0 ) 
    721          END DO 
    722       END DO 
    723721 
    724722      SELECT CASE(nn_ldfeiv_shape) 
     723         CASE(0) !! Standard shape applied - decrease in tropics and cap.  
     724            DO jj = 2, jpjm1 
     725               DO ji = fs_2, fs_jpim1   ! vector opt. 
     726                  zzaei = MIN( 1._wp, ABS( ff_t(ji,jj) * z1_f20 ) ) * zaeiw(ji,jj)     ! tropical decrease 
     727                  zaeiw(ji,jj) = MIN( zzaei, paei0 ) 
     728               END DO 
     729            END DO 
     730 
    725731         CASE(1) !! Abrupt cut-off on Rossby radius: 
    726732! JD : modifications here to introduce scaling by local rossby radius of deformation vs local grid scale 
     
    761767               END DO 
    762768            END DO 
     769 
     770         CASE(5) !! Rossby radius ramp type 1 applied to Treguier et al coefficient rather than cap: 
     771                 !! Note the ramp is RR/GS=[2.0,1.0] (not [2.0,0.5] as for cases 2,3) and we ramp up  
     772                 !! to 5% of the Treguier et al coefficient, aiming for peak values of around 100m2/s 
     773                 !! at high latitudes rather than 2000m2/s which is what you get in eORCA025 with an  
     774                 !! uncapped coefficient. 
     775            DO jj = 2, jpjm1 
     776               DO ji = fs_2, fs_jpim1   ! vector opt. 
     777                  zratio(ji,jj) = zRo(ji,jj)/MIN(e1t(ji,jj),e2t(ji,jj)) 
     778                  zaeiw(ji,jj) = MAX( 0._wp, MIN( 1._wp, 2._wp - zratio(ji,jj) ) ) * 0.05 * zaeiw(ji,jj) 
     779                  zaeiw(ji,jj) = MIN( zaeiw(ji,jj), paei0 ) 
     780               END DO 
     781            END DO 
     782            CALL iom_put('RR_GS',zratio) 
     783 
     784         CASE DEFAULT 
     785               CALL ctl_stop('ldf_eiv: Unrecognised option for nn_ldfeiv_shape.')          
    763786 
    764787      END SELECT 
Note: See TracChangeset for help on using the changeset viewer.