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 461 for trunk/NEMO/OPA_SRC/LDF/ldfslp.F90 – NEMO

Ignore:
Timestamp:
2006-05-10T19:15:54+02:00 (18 years ago)
Author:
opalod
Message:

nemo_v1_update_052:RB: update lateral diffusion computation following the reorganization of both dynamics and tracers

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMO/OPA_SRC/LDF/ldfslp.F90

    r258 r461  
    4949   !!---------------------------------------------------------------------- 
    5050   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    51    !! $Header$  
    52    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    5351   !!---------------------------------------------------------------------- 
    5452 
     
    7169      !!      of 10cm/s) 
    7270      !!        A horizontal shapiro filter is applied to the slopes 
    73       !!        'key_s_coord' defined: add to the previously computed slopes 
     71      !!        ln_sco=T, s-coordinate, add to the previously computed slopes 
    7472      !!      the slope of the model level surface. 
    7573      !!        macro-tasked on horizontal slab (jk-loop)  (2, jpk-1) 
    7674      !!      [slopes already set to zero at level 1, and to zero or the ocean 
    77       !!      bottom slope ('key_s_coord' defined) at level jpk in inildf] 
     75      !!      bottom slope (ln_sco=T) at level jpk in inildf] 
    7876      !! 
    7977      !! ** Action : - uslp, wslpi, and vslp, wslpj, the i- and  j-slopes  
     
    8583      !!   8.1  !  99-10  (A. Jouzeau)  NEW profile 
    8684      !!   8.5  !  99-10  (G. Madec)  Free form, F90 
     85      !!   9.0  !  05-10  (A. Beckmann)  correction for s-coordinates 
    8786      !!---------------------------------------------------------------------- 
    8887      !! * Modules used 
     
    9998      !! * Local declarations 
    10099      INTEGER  ::   ji, jj, jk              ! dummy loop indices 
    101       INTEGER  ::   ii0, ii1, ij0, ij1      ! temporary integer 
    102 #if defined key_partial_steps 
    103       INTEGER  ::   iku, ikv  ! temporary integers 
    104 #endif 
     100      INTEGER  ::   ii0, ii1, ij0, ij1,  &  ! temporary integer 
     101         &          iku, ikv                !    "          " 
    105102      REAL(wp) ::   & 
    106          zeps, zmg, zm05g, zcoef1, zcoef2,   &  ! temporary scalars 
    107          zau, zbu, zav, zbv,                 & 
    108          zai, zbi, zaj, zbj,                & 
    109          zcofu, zcofv, zcofw,                & 
    110          z1u, z1v, z1wu, z1wv,               & 
     103         zeps, zmg, zm05g,               &  ! temporary scalars 
     104         zcoef1, zcoef2, zcoef3,         &  ! 
     105         zau, zbu, zav, zbv,             & 
     106         zai, zbi, zaj, zbj,             & 
     107         zcofu, zcofv, zcofw,            & 
     108         z1u, z1v, z1wu, z1wv,           & 
    111109         zalpha 
    112110      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zww 
     
    140138      END DO 
    141139 
    142 #if defined key_partial_steps 
    143       ! partial steps correction at the bottom ocean level (zps_hde routine) 
    144 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
    145       jj = 1 
    146       DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
     140      IF( ln_zps ) THEN      ! partial steps correction at the bottom ocean level (zps_hde routine) 
     141# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
     142         jj = 1 
     143         DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    147144# else 
    148       DO jj = 1, jpjm1 
    149          DO ji = 1, jpim1 
    150 # endif 
    151             ! last ocean level 
    152             iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj) ) - 1 
    153             ikv = MIN ( mbathy(ji,jj), mbathy(ji,jj+1) ) - 1 
    154             zgru(ji,jj,iku) = gru(ji,jj)  
    155             zgrv(ji,jj,ikv) = grv(ji,jj)                
    156 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
    157          END DO 
    158 # endif 
    159       END DO 
    160 #endif 
     145         DO jj = 1, jpjm1 
     146            DO ji = 1, jpim1 
     147# endif 
     148               ! last ocean level 
     149               iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj) ) - 1 
     150               ikv = MIN ( mbathy(ji,jj), mbathy(ji,jj+1) ) - 1 
     151               zgru(ji,jj,iku) = gru(ji,jj)  
     152               zgrv(ji,jj,ikv) = grv(ji,jj)                
     153# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
     154            END DO 
     155# endif 
     156         END DO 
     157      ENDIF 
    161158 
    162159      ! Slopes of isopycnal surfaces just below the mixed layer 
     
    205202               ! uslp and vslp output in zwz and zww, resp. 
    206203               zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 
    207 #if defined key_s_coord  
    208204               zwz (ji,jj,jk) = ( zau / ( zbu - zeps ) * ( 1. - zalpha)   & 
    209                   &        + zalpha * uslpml(ji,jj)   & 
    210                   &        * ( fsdepu(ji,jj,jk) - .5*fse3u(ji,jj,1) )   & 
    211                   &        / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5. ) )   & 
    212                   &        * umask(ji,jj,jk) 
     205                  &           + zalpha * uslpml(ji,jj)                    & 
     206                  &                    * 0.5 * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) )   & 
     207                  &                          / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5. ) ) * umask(ji,jj,jk) 
    213208               zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 
    214209               zww (ji,jj,jk) = ( zav / ( zbv - zeps ) * ( 1. - zalpha)   & 
    215                   &        + zalpha * vslpml(ji,jj)   & 
    216                   &        * ( fsdepv(ji,jj,jk) - .5*fse3v(ji,jj,1) )   & 
    217                   &         / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) )   & 
    218                   &        * vmask(ji,jj,jk) 
    219 #else 
    220                ! z-coord and partial steps slope computed in the same way 
    221                zwz (ji,jj,jk) = ( zau / ( zbu - zeps ) * ( 1. - zalpha)    & 
    222                   &        + zalpha * uslpml(ji,jj)    & 
    223                   &        * ( fsdept(ji,jj,jk) - .5*fse3u(ji,jj,1))    & 
    224                   &        / MAX (hmlpt(ji,jj),hmlpt(ji+1,jj),5.) )    & 
    225                   &        * umask (ji,jj,jk) 
    226                zalpha = MAX(omlmask(ji,jj,jk),omlmask(ji,jj+1,jk)) 
    227                zww (ji,jj,jk) = ( zav / ( zbv - zeps ) * ( 1. - zalpha)    & 
    228                   &        + zalpha * vslpml(ji,jj)    & 
    229                   &        * ( fsdept(ji,jj,jk) - .5*fse3v(ji,jj,1))    & 
    230                   &        / MAX(hmlpt(ji,jj),hmlpt(ji,jj+1),5.) )    & 
    231                   &        * vmask (ji,jj,jk) 
    232 #endif 
     210                  &           + zalpha * vslpml(ji,jj)                    & 
     211                  &                    * 0.5 * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) )   & 
     212                  &                          / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 
    233213            END DO 
    234214         END DO 
     
    294274            END DO 
    295275         END DO 
    296  
    297  
    298          IF( lk_sco ) THEN 
    299             ! Add the slope of level surfaces 
    300             ! ----------------------------------- 
    301             ! 'key_s_coord' defined but not 'key_traldfiso' the computation is done 
    302             ! in inildf, ldfslp never called 
    303             ! 'key_s_coord' and 'key_traldfiso' defined, the slope of level surfaces 
    304             ! is added to the slope of isopycnal surfaces. 
    305             ! c a u t i o n : minus sign as fsdep has positive value  
    306           
    307             DO jj = 2, jpjm1 
    308                DO ji = fs_2, fs_jpim1   ! vector opt. 
    309                   uslp(ji,jj,jk) = uslp(ji,jj,jk) - 1. / e1u(ji,jj)   & 
    310                      &           * ( fsdept(ji+1,jj,jk) - fsdept(ji,jj,jk) ) 
    311                   vslp(ji,jj,jk) = vslp(ji,jj,jk) - 1. / e2v(ji,jj)   & 
    312                      &           * ( fsdept(ji,jj+1,jk) - fsdept(ji,jj,jk) ) 
    313                END DO 
    314             END DO 
    315          ENDIF 
    316276 
    317277 
     
    356316               zbj = MIN( zwy (ji,jj,jk), -100.*ABS(zaj), -7.e+3/fse3w(ji,jj,jk)*ABS(zaj) ) 
    357317               ! wslpi and wslpj output in zwz and zww, resp. 
    358                zalpha = MAX(omlmask(ji,jj,jk),omlmask(ji,jj,jk-1)) 
    359                zwz(ji,jj,jk) = ( zai / ( zbi - zeps) * ( 1. - zalpha )   & 
    360                   &            + zalpha * wslpiml(ji,jj)   & 
    361                   &            * fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj),10. ) )   & 
    362                   &            * tmask (ji,jj,jk) 
    363                zww(ji,jj,jk) = ( zaj / ( zbj - zeps) * ( 1. - zalpha )   & 
    364                   &            + zalpha * wslpjml(ji,jj)   & 
    365                   &            * fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj),10. ) )   & 
    366                   &            * tmask (ji,jj,jk) 
     318               zalpha = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) ) 
     319               zcoef3 = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10. ) 
     320               zwz(ji,jj,jk) = (     zai / ( zbi - zeps)  * ( 1. - zalpha )   & 
     321                  &             + zcoef3 * wslpiml(ji,jj) *        zalpha   ) * tmask (ji,jj,jk) 
     322               zww(ji,jj,jk) = (     zaj / ( zbj - zeps)  * ( 1. - zalpha )   & 
     323                  &             + zcoef3 * wslpjml(ji,jj) *        zalpha   ) * tmask (ji,jj,jk) 
    367324            END DO 
    368325         END DO 
     
    426383         END DO 
    427384          
    428          IF( lk_sco ) THEN 
    429           
    430             ! Slope of level surfaces 
    431             ! ----------------------- 
    432             ! 'key_s_coord' defined but not 'key_traldfiso' the computation is done 
    433             ! in inildf, ldfslp never called 
    434             ! 'key_s_coord' and 'key_traldfiso' defined, the slope of level surfaces 
    435             ! is added to the slope of isopycnal surfaces. 
    436           
    437             DO jj = 2, jpjm1 
    438                DO ji = fs_2, fs_jpim1   ! vector opt. 
    439                   wslpi(ji,jj,jk) = wslpi(ji,jj,jk) - 1. / e1t(ji,jj)   & 
    440                      &                                   * ( fsdepuw(ji+1,jj,jk) - fsdepuw(ji,jj,jk) ) 
    441                   wslpj(ji,jj,jk) = wslpj(ji,jj,jk) - 1. / e2t(ji,jj)   & 
    442                      &                                   * ( fsdepvw(ji,jj+1,jk) - fsdepvw(ji,jj,jk) ) 
    443                END DO 
    444             END DO 
    445          ENDIF 
    446385          
    447386         ! III. Specific grid points 
     
    476415      ! III Lateral boundary conditions on all slopes (uslp , vslp,  
    477416      ! -------------------------------                wslpi, wslpj ) 
    478       CALL lbc_lnk( uslp , 'U', -1. ) 
    479       CALL lbc_lnk( vslp , 'V', -1. ) 
    480       CALL lbc_lnk( wslpi, 'W', -1. ) 
    481       CALL lbc_lnk( wslpj, 'W', -1. ) 
     417      CALL lbc_lnk( uslp , 'U', -1. )      ;      CALL lbc_lnk( vslp , 'V', -1. ) 
     418      CALL lbc_lnk( wslpi, 'W', -1. )      ;      CALL lbc_lnk( wslpj, 'W', -1. ) 
    482419 
    483420      IF(ln_ctl) THEN 
     
    553490      ! mask for mixed layer 
    554491      DO jk = 1, jpk 
    555 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     492# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    556493         jj = 1 
    557494         DO ji = 1, jpij   ! vector opt. (forced unrolling) 
     
    567504                  omlmask(ji,jj,jk) = 0.e0 
    568505               ENDIF 
    569 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     506# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    570507            END DO 
    571508# endif 
     
    585522      zwy(:,jpj) = 0.e0 
    586523      zwy(jpi,:) = 0.e0 
    587 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     524# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    588525      jj = 1 
    589526      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
     
    598535               &             * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) )   & 
    599536               &             / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 
    600 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     537# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    601538         END DO 
    602539# endif 
     
    606543 
    607544      ! Slope at u points 
    608 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     545# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    609546      jj = 1 
    610547      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
     
    623560            ! uslpml 
    624561            uslpml (ji,jj) = zau / ( zbu - zeps ) * umask (ji,jj,ik) 
    625 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     562# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    626563         END DO 
    627564# endif 
     
    635572      zwy ( :, jpj) = 0.e0 
    636573      zwy ( jpi, :) = 0.e0 
    637 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     574# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    638575      jj = 1 
    639576      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
     
    647584               &             * ( pn2(ji,jj,ik) + pn2(ji,jj,ik+1) )   & 
    648585               &             / MAX( tmask(ji,jj,ik) + tmask (ji,jj,ik+1), 1. ) 
    649 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     586# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    650587         END DO 
    651588# endif 
     
    656593 
    657594      ! Slope at v points 
    658 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     595# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    659596      jj = 1 
    660597      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
     
    673610            ! vslpml 
    674611            vslpml (ji,jj) = zav / ( zbv - zeps ) * vmask (ji,jj,ik) 
    675 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     612# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    676613         END DO 
    677614# endif 
     
    687624      ! Local vertical density gradient evaluated from N^2 
    688625      ! zwy = d/dz(prd)= - mk ( prd ) / grav * pn2 -- at w point 
    689 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     626# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    690627      jj = 1 
    691628      DO ji = 1, jpij   ! vector opt. (forced unrolling) 
     
    699636            zwy (ji,jj) = zm05g * pn2 (ji,jj,ik) *     & 
    700637               &             ( prd (ji,jj,ik) + prd (ji,jj,ikm1) + 2. ) 
    701 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     638# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    702639         END DO 
    703640# endif 
     
    705642 
    706643      ! Slope at w point 
    707 # if defined key_vectopt_loop   &&   ! defined key_autotasking 
     644# if defined key_vectopt_loop   &&   ! defined key_mpp_omp 
    708645      jj = 1 
    709646      DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
     
    735672            wslpiml (ji,jj) = zai / ( zbi - zeps) * tmask (ji,jj,ik) 
    736673            wslpjml (ji,jj) = zaj / ( zbj - zeps) * tmask (ji,jj,ik) 
    737 # if ! defined key_vectopt_loop   ||   defined key_autotasking 
     674# if ! defined key_vectopt_loop   ||   defined key_mpp_omp 
    738675         END DO 
    739676# endif 
     
    787724 
    788725      IF( ln_traldf_hor .OR. ln_dynldf_hor ) THEN 
     726         IF(lwp) THEN 
     727            WRITE(numout,*) '          Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 
     728         ENDIF 
    789729 
    790730         ! geopotential diffusion in s-coordinates on tracers and/or momentum 
     
    797737            DO jj = 2, jpjm1 
    798738               DO ji = fs_2, fs_jpim1   ! vector opt. 
    799                   uslp (ji,jj,jk) = -1. / e1u(ji,jj) * umask(ji,jj,jk)   & 
    800                      &                               * ( fsdept(ji+1,jj,jk) - fsdept(ji,jj,jk) ) 
    801                   vslp (ji,jj,jk) = -1. / e2v(ji,jj) * vmask(ji,jj,jk)   & 
    802                      &                               * ( fsdept(ji,jj+1,jk) - fsdept(ji,jj,jk) ) 
    803                   wslpi(ji,jj,jk) = -1. / e1t(ji,jj) * tmask(ji,jj,jk)   & 
    804                      &                               * ( fsdepuw(ji+1,jj,jk) - fsdepuw(ji,jj,jk) ) 
    805                   wslpj(ji,jj,jk) = -1. / e2t(ji,jj) * tmask(ji,jj,jk)   & 
    806                      &                               * ( fsdepvw(ji,jj+1,jk) - fsdepvw(ji,jj,jk) ) 
     739                  uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk) 
     740                  vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk) 
     741                  wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw(ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 
     742                  wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw(ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 
    807743               END DO 
    808744            END DO 
     
    810746 
    811747         ! Lateral boundary conditions on the slopes 
    812          CALL lbc_lnk( uslp , 'U', -1. ) 
    813          CALL lbc_lnk( vslp , 'V', -1. ) 
    814          CALL lbc_lnk( wslpi, 'W', -1. ) 
    815          CALL lbc_lnk( wslpj, 'W', -1. ) 
     748         CALL lbc_lnk( uslp , 'U', -1. )      ;      CALL lbc_lnk( vslp , 'V', -1. ) 
     749         CALL lbc_lnk( wslpi, 'W', -1. )      ;      CALL lbc_lnk( wslpj, 'W', -1. ) 
    816750      ENDIF 
    817751 
Note: See TracChangeset for help on using the changeset viewer.