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 3077 for branches/2011/dev_NOC_2011_MERGE – NEMO

Ignore:
Timestamp:
2011-11-10T18:23:33+01:00 (12 years ago)
Author:
acc
Message:

Branch dev_NOC_2011_MERGE. Minor correction to ldfslp (Griffies)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2011/dev_NOC_2011_MERGE/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90

    r2980 r3077  
    410410      INTEGER  ::   iku, ikv                    ! local integer 
    411411      REAL(wp) ::   zfacti, zfactj              ! local scalars 
     412      REAL(wp) ::   znot_thru_surface           ! local scalars 
    412413      REAL(wp) ::   zdit, zdis, zdjt, zdjs, zdkt, zdks, zbu, zbv, zbti, zbtj 
    413       REAL(wp) ::   zdxrho_raw, zti_coord, zti_raw, zti_lim, zti_lim2, zti_g_raw, zti_g_lim 
    414       REAL(wp) ::   zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_lim2, ztj_g_raw, ztj_g_lim 
     414      REAL(wp) ::   zdxrho_raw, zti_coord, zti_raw, zti_lim, zti_g_raw, zti_g_lim 
     415      REAL(wp) ::   zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_g_raw, ztj_g_lim 
    415416      REAL(wp) ::   zdzrho_raw 
    416417      REAL(wp) ::   zbeta0 
     
    531532            ip = jl   ;   jp = jl             ! i- and j-indices of triads (i-k and j-k planes) 
    532533            DO jk = 1, jpkm1 
     534               ! Must mask contribution to slope from dz/dx at constant s for triads jk=1,kp=0 that poke up though ocean surface 
     535               znot_thru_surface = REAL( 1-1/(jk+kp), wp )  !jk+kp=1,=0.; otherwise=1.0 
    533536               DO jj = 1, jpjm1 
    534537                  DO ji = 1, fs_jpim1         ! vector opt. 
     
    543546                     zti_raw   = zdxrho(ji+ip,jj   ,jk,1-ip) / zdzrho(ji+ip,jj   ,jk,kp)                   ! unmasked 
    544547                     ztj_raw   = zdyrho(ji   ,jj+jp,jk,1-jp) / zdzrho(ji   ,jj+jp,jk,kp) 
    545                      zti_coord = ( fsdept(ji+1,jj  ,jk) - fsdept(ji,jj,jk) ) / e1u(ji,jj) 
    546                      ztj_coord = ( fsdept(ji  ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj)                  ! unmasked 
     548 
     549                     ! Must mask contribution to slope for triad jk=1,kp=0 that poke up though ocean surface 
     550                     zti_coord = znot_thru_surface * ( fsdept(ji+1,jj  ,jk) - fsdept(ji,jj,jk) ) / e1u(ji,jj) 
     551                     ztj_coord = znot_thru_surface * ( fsdept(ji  ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj)                  ! unmasked 
    547552                     zti_g_raw = zti_raw - zti_coord      ! ref to geopot surfaces 
    548553                     ztj_g_raw = ztj_raw - ztj_coord 
     
    550555                     ztj_g_lim = SIGN( MIN( rn_slpmax, ABS( ztj_g_raw ) ), ztj_g_raw ) 
    551556                     ! 
    552                      ! Below  ML use limited zti_g as is 
    553                      ! Inside ML replace by linearly reducing sx_mlb towards surface 
     557                     ! Below  ML use limited zti_g as is & mask 
     558                     ! Inside ML replace by linearly reducing sx_mlb towards surface & mask 
    554559                     ! 
    555560                     zfacti = REAL( 1 - 1/(1 + (jk+kp-1)/nmln(ji+ip,jj)), wp )  ! k index of uppermost point(s) of triad is jk+kp-1 
    556561                     zfactj = REAL( 1 - 1/(1 + (jk+kp-1)/nmln(ji,jj+jp)), wp )  ! must be .ge. nmln(ji,jj) for zfact=1 
    557562                     !                                                          !                   otherwise  zfact=0 
    558                      zti_g_lim =          zfacti   * zti_g_lim                       & 
     563                     zti_g_lim =          ( zfacti   * zti_g_lim                       & 
    559564                        &      + ( 1._wp - zfacti ) * zti_mlb(ji+ip,jj,1-ip,kp)   & 
    560                         &                           * fsdepw(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) 
    561                      ztj_g_lim =          zfactj   * ztj_g_lim                       & 
     565                        &                           * fsdepw(ji+ip,jj,jk+kp) * z1_mlbw(ji+ip,jj) ) * umask(ji,jj,jk+kp) 
     566                     ztj_g_lim =          ( zfactj   * ztj_g_lim                       & 
    562567                        &      + ( 1._wp - zfactj ) * ztj_mlb(ji,jj+jp,1-jp,kp)   & 
    563                         &                           * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) 
    564                      ! 
    565                      triadi_g(ji+ip,jj   ,jk,1-ip,kp) = zti_g_lim * umask(ji,jj,jk+kp)                     ! masked 
    566                      triadj_g(ji   ,jj+jp,jk,1-jp,kp) = ztj_g_lim * vmask(ji,jj,jk+kp) 
     568                        &                           * fsdepw(ji,jj+jp,jk+kp) * z1_mlbw(ji,jj+jp) ) * vmask(ji,jj,jk+kp) 
     569                     ! 
     570                     triadi_g(ji+ip,jj   ,jk,1-ip,kp) = zti_g_lim 
     571                     triadj_g(ji   ,jj+jp,jk,1-jp,kp) = ztj_g_lim 
    567572                     ! 
    568573                     ! Get coefficients of isoneutral diffusion tensor 
     
    575580                     zti_lim  = ( zti_g_lim + zti_coord ) * umask(ji,jj,jk+kp)    ! remove coordinate slope => relative to coordinate surfaces 
    576581                     ztj_lim  = ( ztj_g_lim + ztj_coord ) * vmask(ji,jj,jk+kp) 
    577                      zti_lim2 = zti_lim * zti_lim      ! square of limited slopes                          ! masked <<== 
    578                      ztj_lim2 = ztj_lim * ztj_lim 
    579582                     ! 
    580583                     IF( ln_triad_iso ) THEN 
    581                         zti_raw = zti_lim2 / zti_raw 
    582                         ztj_raw = ztj_lim2 / ztj_raw 
     584                        zti_raw = ( zti_lim*zti_lim ) / zti_raw 
     585                        ztj_raw = ( ztj_lim*ztj_lim ) / ztj_raw 
    583586                        zti_raw = SIGN( MIN( ABS(zti_lim), ABS( zti_raw ) ), zti_raw ) 
    584587                        ztj_raw = SIGN( MIN( ABS(ztj_lim), ABS( ztj_raw ) ), ztj_raw ) 
     
    597600                     ! 
    598601                     !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers....  ==> to be checked 
    599                      ! agn may need to change to using ztj_g_lim**2, as wslp2 just used for eddy growth rate, needs *real* slope 
    600                      wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * zti_lim2      ! masked 
    601                      wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ztj_lim2 
     602                     wslp2 (ji+ip,jj,jk+kp) = wslp2(ji+ip,jj,jk+kp) + 0.25_wp * zbu / zbti * ( zti_g_lim * zti_g_lim )      ! masked 
     603                     wslp2 (ji,jj+jp,jk+kp) = wslp2(ji,jj+jp,jk+kp) + 0.25_wp * zbv / zbtj * ( ztj_g_lim * ztj_g_lim ) 
    602604                  END DO 
    603605               END DO 
Note: See TracChangeset for help on using the changeset viewer.