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 5075 for branches/2014/dev_r4650_UKMO7_STARTHOUR/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90 – NEMO

Ignore:
Timestamp:
2015-02-11T11:50:34+01:00 (9 years ago)
Author:
timgraham
Message:

Upgraded branch to current head of trunk (r5072) so it can be used with the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4650_UKMO7_STARTHOUR/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90

    r4488 r5075  
    2828   USE zdfmxl         ! mixed layer depth 
    2929   USE eosbn2         ! equation of states 
     30   ! 
     31   USE in_out_manager ! I/O manager 
    3032   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    31    USE in_out_manager ! I/O manager 
    3233   USE prtctl         ! Print control 
    3334   USE wrk_nemo       ! work arrays 
     
    108109      REAL(wp) ::   zcj, zfj, zav, zbv, zaj, zbj   !   -      - 
    109110      REAL(wp) ::   zck, zfk,      zbw             !   -      - 
     111      REAL(wp) ::   zdepv, zdepu         !   -      - 
    110112      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zwz, zww 
    111113      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zdzr 
    112114      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zgru, zgrv 
     115      REAL(wp), POINTER, DIMENSION(:,:  ) :: zhmlpu, zhmlpv 
    113116      !!---------------------------------------------------------------------- 
    114117      ! 
     
    116119      ! 
    117120      CALL wrk_alloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 
     121      CALL wrk_alloc( jpi,jpj, zhmlpu, zhmlpv ) 
    118122 
    119123      IF ( ln_traldf_iso .OR. ln_dynldf_iso ) THEN  
     
    136140         END DO 
    137141         IF( ln_zps ) THEN                           ! partial steps correction at the bottom ocean level 
    138 # if defined key_vectopt_loop 
    139             DO jj = 1, 1 
    140                DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    141 # else 
    142142            DO jj = 1, jpjm1 
    143143               DO ji = 1, jpim1 
    144 # endif 
    145                   zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 
    146                   zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 
     144! IF should be useless check zpshde (PM) 
     145               IF ( mbku(ji,jj) > 1 ) zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 
     146               IF ( mbkv(ji,jj) > 1 ) zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 
     147               IF ( miku(ji,jj) > 1 ) zgru(ji,jj,miku(ji,jj)) = grui(ji,jj)  
     148               IF ( mikv(ji,jj) > 1 ) zgrv(ji,jj,mikv(ji,jj)) = grvi(ji,jj) 
    147149               END DO 
    148150            END DO 
     
    150152         ! 
    151153         zdzr(:,:,1) = 0._wp        !==   Local vertical density gradient at T-point   == !   (evaluated from N^2) 
    152          DO jk = 2, jpkm1 
     154         DO jk = 1, jpkm1 
    153155            !                                ! zdzr = d/dz(prd)= - ( prd ) / grav * mk(pn2) -- at t point 
    154156            !                                !   trick: tmask(ik  )  = 0   =>   all pn2   = 0   =>   zdzr = 0 
     
    159161               &                 * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 
    160162         END DO 
     163         ! surface initialisation  
     164         DO jj = 1, jpjm1 
     165            DO ji = 1, jpim1 
     166              zdzr(ji,jj,1:mikt(ji,jj)) = 0._wp 
     167            END DO 
     168         END DO 
    161169         ! 
    162170         !                          !==   Slopes just below the mixed layer   ==! 
     
    167175         ! ===========================      | vslp = d/dj( prd ) / d/dz( prd ) 
    168176         ! 
     177         DO jj = 2, jpjm1 
     178            DO ji = fs_2, fs_jpim1   ! vector opt. 
     179               IF (miku(ji,jj) .GT. miku(ji+1,jj)) zhmlpu(ji,jj) = hmlpt(ji  ,jj) 
     180               IF (miku(ji,jj) .LT. miku(ji+1,jj)) zhmlpu(ji,jj) = hmlpt(ji+1,jj) 
     181               IF (miku(ji,jj) .EQ. miku(ji+1,jj)) zhmlpu(ji,jj) = MAX(hmlpt(ji  ,jj), hmlpt(ji+1,jj)) 
     182               IF (mikv(ji,jj) .GT. miku(ji,jj+1)) zhmlpv(ji,jj) = hmlpt(ji  ,jj) 
     183               IF (mikv(ji,jj) .LT. miku(ji,jj+1)) zhmlpv(ji,jj) = hmlpt(ji,jj+1) 
     184               IF (mikv(ji,jj) .EQ. miku(ji,jj+1)) zhmlpv(ji,jj) = MAX(hmlpt(ji,jj), hmlpt(ji,jj+1)) 
     185            ENDDO 
     186         ENDDO 
    169187         DO jk = 2, jpkm1                            !* Slopes at u and v points 
    170188            DO jj = 2, jpjm1 
     
    180198                  zbv = MIN(  zbv, -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,jk)* ABS( zav )  ) 
    181199                  !                                      ! uslp and vslp output in zwz and zww, resp. 
    182                   zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) ) 
    183                   zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) ) 
    184                   zwz(ji,jj,jk) = ( ( 1. - zfi) * zau / ( zbu - zeps )                                              & 
    185                      &                   + zfi  * uslpml(ji,jj)                                                     & 
    186                      &                          * 0.5_wp * ( fsdept(ji+1,jj,jk)+fsdept(ji,jj,jk)-fse3u(ji,jj,1) )   & 
    187                      &                          / MAX( hmlpt(ji,jj), hmlpt(ji+1,jj), 5._wp ) ) * umask(ji,jj,jk) 
    188                   zww(ji,jj,jk) = ( ( 1. - zfj) * zav / ( zbv - zeps )                                              & 
    189                      &                   + zfj  * vslpml(ji,jj)                                                     & 
    190                      &                          * 0.5_wp * ( fsdept(ji,jj+1,jk)+fsdept(ji,jj,jk)-fse3v(ji,jj,1) )   & 
    191                      &                          / MAX( hmlpt(ji,jj), hmlpt(ji,jj+1), 5. ) ) * vmask(ji,jj,jk) 
     200                  zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) )  
     201                  zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) )  
     202                  ! thickness of water column between surface and level k at u/v point 
     203                  zdepu = 0.5_wp * (( fsdept(ji,jj,jk) + fsdept(ji+1,jj  ,jk) )                   & 
     204                             - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj  ) )  & 
     205                             - fse3u(ji,jj,miku(ji,jj))                                         ) 
     206                  zdepv = 0.5_wp * (( fsdept(ji,jj,jk) + fsdept(ji  ,jj+1,jk) )                   & 
     207                             - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) ) & 
     208                             - fse3v(ji,jj,mikv(ji,jj))                                         ) 
     209                  zwz(ji,jj,jk) = ( 1. - zfi) * zau / ( zbu - zeps )                                              & 
     210                     &                 + zfi  * uslpml(ji,jj)                                                     & 
     211                     &                        * zdepu / MAX( zhmlpu(ji,jj), 5._wp ) 
     212                  zwz(ji,jj,jk) = zwz(ji,jj,jk) * umask(ji,jj,jk) * umask(ji,jj,jk-1) 
     213                  zww(ji,jj,jk) = ( 1. - zfj) * zav / ( zbv - zeps )                                              & 
     214                     &                 + zfj  * vslpml(ji,jj)                                                     & 
     215                     &                        * zdepv / MAX( zhmlpv(ji,jj), 5._wp )  
     216                  zww(ji,jj,jk) = zww(ji,jj,jk) * vmask(ji,jj,jk) * vmask(ji,jj,jk-1) 
     217                   
     218                  
    192219!!gm  modif to suppress omlmask.... (as in Griffies case) 
    193220!                  !                                         ! jk must be >= ML level for zf=1. otherwise  zf=0. 
     
    238265               DO ji = fs_2, fs_jpim1   ! vector opt. 
    239266                  uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk  ) ) * 0.5_wp   & 
    240                      &                            * ( umask(ji,jj  ,jk) + umask(ji,jj  ,jk+1) ) * 0.5_wp 
     267                     &                            * ( umask(ji,jj  ,jk) + umask(ji,jj  ,jk+1) ) * 0.5_wp   & 
     268                     &                            *   umask(ji,jj,jk-1) !* umask(ji,jj,jk) * umask(ji,jj,jk+1) 
    241269                  vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk  ) ) * 0.5_wp   & 
    242                      &                            * ( vmask(ji  ,jj,jk) + vmask(ji  ,jj,jk+1) ) * 0.5_wp 
     270                     &                            * ( vmask(ji  ,jj,jk) + vmask(ji  ,jj,jk+1) ) * 0.5_wp   & 
     271                     &                            *   vmask(ji,jj,jk-1) !* vmask(ji,jj,jk) * vmask(ji,jj,jk+1) 
    243272               END DO 
    244273            END DO 
     
    253282               DO ji = fs_2, fs_jpim1   ! vector opt. 
    254283                  !                                  !* Local vertical density gradient evaluated from N^2 
    255                   zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) 
     284                  zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
    256285                  !                                  !* Slopes at w point 
    257286                  !                                        ! i- & j-gradient of density at w-points 
     
    261290                     &      + vmask(ji,jj-1,jk-1) + vmask(ji,jj,jk  ) , zeps  ) *  e2t(ji,jj) 
    262291                  zai =    (  zgru (ji-1,jj,jk  ) + zgru (ji,jj,jk-1)           & 
    263                      &      + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk  )   ) / zci * tmask (ji,jj,jk) 
     292                     &      + zgru (ji-1,jj,jk-1) + zgru (ji,jj,jk  )   ) / zci 
    264293                  zaj =    (  zgrv (ji,jj-1,jk  ) + zgrv (ji,jj,jk-1)           & 
    265                      &      + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk  )   ) / zcj * tmask (ji,jj,jk) 
     294                     &      + zgrv (ji,jj-1,jk-1) + zgrv (ji,jj,jk  )   ) / zcj 
    266295                  !                                        ! bound the slopes: abs(zw.)<= 1/100 and zb..<0. 
    267296                  !                                        ! + kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 
     
    269298                  zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zaj )  ) 
    270299                  !                                        ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 
    271                   zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) )   ! zfk=1 in the ML otherwise zfk=0 
    272                   zck = fsdepw(ji,jj,jk) / MAX( hmlp(ji,jj), 10._wp ) 
    273                   zwz(ji,jj,jk) = (  zai / ( zbi - zeps ) * ( 1._wp - zfk ) + zck * wslpiml(ji,jj) * zfk  ) * tmask(ji,jj,jk) 
    274                   zww(ji,jj,jk) = (  zaj / ( zbj - zeps ) * ( 1._wp - zfk ) + zck * wslpjml(ji,jj) * zfk  ) * tmask(ji,jj,jk) 
     300                  zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) )    ! zfk=1 in the ML otherwise zfk=0 
     301                  zck = ( fsdepw(ji,jj,jk) - fsdepw(ji,jj,mikt(ji,jj) ) ) / MAX( hmlp(ji,jj), 10._wp ) 
     302                  zwz(ji,jj,jk) = (  zai / ( zbi - zeps ) * ( 1._wp - zfk ) & 
     303                     &            + zck * wslpiml(ji,jj) * zfk  ) * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
     304                  zww(ji,jj,jk) = (  zaj / ( zbj - zeps ) * ( 1._wp - zfk ) & 
     305                     &            + zck * wslpjml(ji,jj) * zfk  ) * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
    275306 
    276307!!gm  modif to suppress omlmask....  (as in Griffies operator) 
     
    325356                  zck =   ( umask(ji,jj,jk) + umask(ji-1,jj,jk) )   & 
    326357                     &  * ( vmask(ji,jj,jk) + vmask(ji,jj-1,jk) ) * 0.25 
    327                   wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck 
    328                   wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck 
     358                  wslpi(ji,jj,jk) = wslpi(ji,jj,jk) * zck * tmask(ji,jj,jk-1) * tmask(ji,jj,jk) 
     359                  wslpj(ji,jj,jk) = wslpj(ji,jj,jk) * zck * tmask(ji,jj,jk-1) * tmask(ji,jj,jk) 
    329360               END DO 
    330361            END DO 
     
    377408         ! set the slope of diffusion to the slope of s-surfaces  
    378409         !      ( c a u t i o n : minus sign as fsdep has positive value )  
    379          DO jk = 1, jpk  
     410         DO jj = 2, jpjm1  
     411            DO ji = fs_2, fs_jpim1   ! vector opt.  
     412               uslp(ji,jj,1) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,1) - fsdept_b(ji ,jj ,1) )  * umask(ji,jj,1)  
     413               vslp(ji,jj,1) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,1) - fsdept_b(ji ,jj ,1) )  * vmask(ji,jj,1)  
     414               wslpi(ji,jj,1) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,1) - fsdepw_b(ji-1,jj,1) ) * tmask(ji,jj,1) * 0.5  
     415               wslpj(ji,jj,1) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,1) - fsdepw_b(ji,jj-1,1) ) * tmask(ji,jj,1) * 0.5  
     416            END DO  
     417         END DO  
     418 
     419         DO jk = 2, jpk  
    380420            DO jj = 2, jpjm1  
    381421               DO ji = fs_2, fs_jpim1   ! vector opt.  
    382422                  uslp(ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept_b(ji+1,jj,jk) - fsdept_b(ji ,jj ,jk) ) * umask(ji,jj,jk)  
    383423                  vslp(ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept_b(ji,jj+1,jk) - fsdept_b(ji ,jj ,jk) ) * vmask(ji,jj,jk)  
    384                   wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5  
    385                   wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5  
     424                  wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw_b(ji+1,jj,jk) - fsdepw_b(ji-1,jj,jk) ) & 
     425                    &                              * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) * 0.5  
     426                  wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw_b(ji,jj+1,jk) - fsdepw_b(ji,jj-1,jk) ) & 
     427                    &                              * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) * 0.5  
    386428               END DO  
    387429            END DO  
     
    405447       
    406448      CALL wrk_dealloc( jpi,jpj,jpk, zwz, zww, zdzr, zgru, zgrv ) 
     449      CALL wrk_dealloc( jpi,jpj,     zhmlpu, zhmlpv) 
    407450      ! 
    408451      IF( nn_timing == 1 )  CALL timing_stop('ldf_slp') 
     
    435478      REAL(wp) ::   zdyrho_raw, ztj_coord, ztj_raw, ztj_lim, ztj_g_raw, ztj_g_lim 
    436479      REAL(wp) ::   zdzrho_raw 
    437       REAL(wp) ::   zbeta0 
    438480      REAL(wp), POINTER, DIMENSION(:,:)     ::   z1_mlbw 
    439       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalbet 
    440481      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   zdxrho , zdyrho, zdzrho     ! Horizontal and vertical density gradients 
    441482      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   zti_mlb, ztj_mlb            ! for Griffies operator only 
     
    445486      ! 
    446487      CALL wrk_alloc( jpi,jpj, z1_mlbw ) 
    447       CALL wrk_alloc( jpi,jpj,jpk, zalbet ) 
    448488      CALL wrk_alloc( jpi,jpj,jpk,2, zdxrho , zdyrho, zdzrho,              klstart = 0  ) 
    449489      CALL wrk_alloc( jpi,jpj,  2,2, zti_mlb, ztj_mlb,        kkstart = 0, klstart = 0  ) 
     
    452492      !  Some preliminary calculation  ! 
    453493      !--------------------------------! 
    454       ! 
    455       CALL eos_alpbet( tsb, zalbet, zbeta0 )  !==  before local thermal/haline expension ratio at T-points  ==! 
    456494      ! 
    457495      DO jl = 0, 1                            !==  unmasked before density i- j-, k-gradients  ==! 
     
    465503                  zdjt = ( tsb(ji,jj+1,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) )    ! j-gradient of T & S at v-point 
    466504                  zdjs = ( tsb(ji,jj+1,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) 
    467                   zdxrho_raw = ( - zalbet(ji+ip,jj   ,jk) * zdit + zbeta0*zdis ) / e1u(ji,jj) 
    468                   zdyrho_raw = ( - zalbet(ji   ,jj+jp,jk) * zdjt + zbeta0*zdjs ) / e2v(ji,jj) 
    469                   zdxrho(ji+ip,jj   ,jk,1-ip) = SIGN( MAX(   repsln, ABS( zdxrho_raw ) ), zdxrho_raw )   ! keep the sign 
     505                  zdxrho_raw = ( - rab_b(ji+ip,jj   ,jk,jp_tem) * zdit + rab_b(ji+ip,jj   ,jk,jp_sal) * zdis ) / e1u(ji,jj) 
     506                  zdyrho_raw = ( - rab_b(ji   ,jj+jp,jk,jp_tem) * zdjt + rab_b(ji   ,jj+jp,jk,jp_sal) * zdjs ) / e2v(ji,jj) 
     507                  zdxrho(ji+ip,jj   ,jk,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw )   ! keep the sign 
    470508                  zdyrho(ji   ,jj+jp,jk,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 
    471509               END DO 
     
    473511         END DO 
    474512         ! 
    475          IF( ln_zps.and.l_grad_zps ) THEN     ! partial steps: correction of i- & j-grad on bottom 
    476 # if defined key_vectopt_loop 
    477             DO jj = 1, 1 
    478                DO ji = 1, jpij-jpi            ! vector opt. (forced unrolling) 
    479 # else 
     513         IF( ln_zps .AND. l_grad_zps ) THEN     ! partial steps: correction of i- & j-grad on bottom 
    480514            DO jj = 1, jpjm1 
    481515               DO ji = 1, jpim1 
    482 # endif 
    483516                  iku  = mbku(ji,jj)          ;   ikv  = mbkv(ji,jj)             ! last ocean level (u- & v-points) 
    484517                  zdit = gtsu(ji,jj,jp_tem)   ;   zdjt = gtsv(ji,jj,jp_tem)      ! i- & j-gradient of Temperature 
    485518                  zdis = gtsu(ji,jj,jp_sal)   ;   zdjs = gtsv(ji,jj,jp_sal)      ! i- & j-gradient of Salinity 
    486                   zdxrho_raw = ( - zalbet(ji+ip,jj   ,iku) * zdit + zbeta0*zdis ) / e1u(ji,jj) 
    487                   zdyrho_raw = ( - zalbet(ji   ,jj+jp,ikv) * zdjt + zbeta0*zdjs ) / e2v(ji,jj) 
     519                  zdxrho_raw = ( - rab_b(ji+ip,jj   ,iku,jp_tem) * zdit + rab_b(ji+ip,jj   ,iku,jp_sal) * zdis ) / e1u(ji,jj) 
     520                  zdyrho_raw = ( - rab_b(ji   ,jj+jp,ikv,jp_tem) * zdjt + rab_b(ji   ,jj+jp,ikv,jp_sal) * zdjs ) / e2v(ji,jj) 
    488521                  zdxrho(ji+ip,jj   ,iku,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw )   ! keep the sign 
    489522                  zdyrho(ji   ,jj+jp,ikv,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 
     
    505538                     zdks = 0._wp 
    506539                  ENDIF 
    507                   zdzrho_raw = ( - zalbet(ji   ,jj   ,jk) * zdkt + zbeta0*zdks ) / fse3w(ji,jj,jk+kp) 
    508                   zdzrho(ji   ,jj   ,jk,  kp) =     - MIN( - repsln,      zdzrho_raw )    ! force zdzrho >= repsln 
     540                  zdzrho_raw = ( - rab_b(ji,jj,jk,jp_tem) * zdkt + rab_b(ji,jj,jk,jp_sal) * zdks ) / fse3w(ji,jj,jk+kp) 
     541                  zdzrho(ji,jj,jk,kp) = - MIN( - repsln, zdzrho_raw )    ! force zdzrho >= repsln 
    509542                 END DO 
    510543            END DO 
     
    650683      ! 
    651684      CALL wrk_dealloc( jpi,jpj, z1_mlbw ) 
    652       CALL wrk_dealloc( jpi,jpj,jpk, zalbet ) 
    653685      CALL wrk_dealloc( jpi,jpj,jpk,2, zdxrho , zdyrho, zdzrho,              klstart = 0  ) 
    654686      CALL wrk_dealloc( jpi,jpj,  2,2, zti_mlb, ztj_mlb,        kkstart = 0, klstart = 0  ) 
     
    701733      !                                            !==   surface mixed layer mask   ! 
    702734      DO jk = 1, jpk                               ! =1 inside the mixed layer, =0 otherwise 
    703 # if defined key_vectopt_loop 
    704          DO jj = 1, 1 
    705             DO ji = 1, jpij                        ! vector opt. (forced unrolling) 
    706 # else 
    707735         DO jj = 1, jpj 
    708736            DO ji = 1, jpi 
    709 # endif 
    710737               ik = nmln(ji,jj) - 1 
    711                IF( jk <= ik ) THEN   ;   omlmask(ji,jj,jk) = 1._wp 
     738               IF( jk <= ik .AND. jk >= mikt(ji,jj) ) THEN   ;   omlmask(ji,jj,jk) = 1._wp 
    712739               ELSE                  ;   omlmask(ji,jj,jk) = 0._wp 
    713740               ENDIF 
     
    727754      !----------------------------------------------------------------------- 
    728755      ! 
    729 # if defined key_vectopt_loop 
    730       DO jj = 1, 1 
    731          DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    732 # else 
    733756      DO jj = 2, jpjm1 
    734757         DO ji = 2, jpim1 
    735 # endif 
    736758            !                        !==   Slope at u- & v-points just below the Mixed Layer   ==! 
    737759            ! 
    738760            !                        !- vertical density gradient for u- and v-slopes (from dzr at T-point) 
    739             iku = MIN(  MAX( 1, nmln(ji,jj) , nmln(ji+1,jj) ) , jpkm1  )   ! ML (MAX of T-pts, bound by jpkm1) 
    740             ikv = MIN(  MAX( 1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1  )   ! 
     761            iku = MIN(  MAX( miku(ji,jj)+1, nmln(ji,jj) , nmln(ji+1,jj) ) , jpkm1  )   ! ML (MAX of T-pts, bound by jpkm1) 
     762            ikv = MIN(  MAX( mikv(ji,jj)+1, nmln(ji,jj) , nmln(ji,jj+1) ) , jpkm1  )   ! 
    741763            zbu = 0.5_wp * ( p_dzr(ji,jj,iku) + p_dzr(ji+1,jj  ,iku) ) 
    742764            zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji  ,jj+1,ikv) ) 
     
    825847         wslpj(:,:,:) = 0._wp   ;   wslpjml(:,:) = 0._wp 
    826848 
    827          IF( ln_traldf_hor .OR. ln_dynldf_hor ) THEN 
     849         IF(ln_sco .AND.  (ln_traldf_hor .OR. ln_dynldf_hor )) THEN 
    828850            IF(lwp)   WRITE(numout,*) '          Horizontal mixing in s-coordinate: slope = slope of s-surfaces' 
    829851 
Note: See TracChangeset for help on using the changeset viewer.