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

Ignore:
Timestamp:
2014-04-06T17:28:25+02:00 (10 years ago)
Author:
gm
Message:

#1260 : see the associated wiki page for explanation

File:
1 edited

Legend:

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

    r4596 r4616  
    144144      END DO 
    145145      IF( ln_zps ) THEN                           ! partial steps correction at the bottom ocean level 
    146 # if defined key_vectopt_loop 
    147          DO jj = 1, 1 
    148             DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling) 
    149 # else 
    150146         DO jj = 1, jpjm1 
    151147            DO ji = 1, jpim1 
    152 # endif 
    153148               zgru(ji,jj,mbku(ji,jj)) = gru(ji,jj) 
    154149               zgrv(ji,jj,mbkv(ji,jj)) = grv(ji,jj) 
     
    179174            DO ji = fs_2, fs_jpim1   ! vector opt. 
    180175               !                                      ! horizontal and vertical density gradient at u- and v-points 
    181                zau = zgru(ji,jj,jk) / e1u(ji,jj) 
    182                zav = zgrv(ji,jj,jk) / e2v(ji,jj) 
     176               zau = zgru(ji,jj,jk) * r1_e1u(ji,jj) 
     177               zav = zgrv(ji,jj,jk) * r1_e2v(ji,jj) 
    183178               zbu = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji+1,jj  ,jk) ) 
    184179               zbv = 0.5_wp * ( zdzr(ji,jj,jk) + zdzr(ji  ,jj+1,jk) ) 
     
    433428                  zdjt = ( tsb(ji,jj+1,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) )    ! j-gradient of T & S at v-point 
    434429                  zdjs = ( tsb(ji,jj+1,jk,jp_sal) - tsb(ji,jj,jk,jp_sal) ) 
    435                   zdxrho_raw = ( - zalbet(ji+ip,jj   ,jk) * zdit + zbeta0*zdis ) / e1u(ji,jj) 
    436                   zdyrho_raw = ( - zalbet(ji   ,jj+jp,jk) * zdjt + zbeta0*zdjs ) / e2v(ji,jj) 
     430                  zdxrho_raw = ( - zalbet(ji+ip,jj   ,jk) * zdit + zbeta0*zdis ) * r1_e1u(ji,jj) 
     431                  zdyrho_raw = ( - zalbet(ji   ,jj+jp,jk) * zdjt + zbeta0*zdjs ) * r1_e2v(ji,jj) 
    437432                  zdxrho(ji+ip,jj   ,jk,1-ip) = SIGN(  MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw  )   ! keep the sign 
    438433                  zdyrho(ji   ,jj+jp,jk,1-jp) = SIGN(  MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw  ) 
     
    447442                  zdit = gtsu(ji,jj,jp_tem)   ;   zdjt = gtsv(ji,jj,jp_tem)      ! i- & j-gradient of Temperature 
    448443                  zdis = gtsu(ji,jj,jp_sal)   ;   zdjs = gtsv(ji,jj,jp_sal)      ! i- & j-gradient of Salinity 
    449                   zdxrho_raw = ( - zalbet(ji+ip,jj   ,iku) * zdit + zbeta0*zdis ) / e1u(ji,jj) 
    450                   zdyrho_raw = ( - zalbet(ji   ,jj+jp,ikv) * zdjt + zbeta0*zdjs ) / e2v(ji,jj) 
     444                  zdxrho_raw = ( - zalbet(ji+ip,jj   ,iku) * zdit + zbeta0*zdis ) * r1_e1u(ji,jj) 
     445                  zdyrho_raw = ( - zalbet(ji   ,jj+jp,ikv) * zdjt + zbeta0*zdjs ) * r1_e2v(ji,jj) 
    451446                  zdxrho(ji+ip,jj   ,iku,1-ip) = SIGN( MAX( repsln, ABS( zdxrho_raw ) ), zdxrho_raw )   ! keep the sign 
    452447                  zdyrho(ji   ,jj+jp,ikv,1-jp) = SIGN( MAX( repsln, ABS( zdyrho_raw ) ), zdyrho_raw ) 
     
    468463                     zdks = 0._wp 
    469464                  ENDIF 
    470                   zdzrho_raw = ( - zalbet(ji   ,jj   ,jk) * zdkt + zbeta0*zdks ) / fse3w(ji,jj,jk+kp) 
    471                   zdzrho(ji   ,jj   ,jk,  kp) =     - MIN( - repsln,      zdzrho_raw )    ! force zdzrho >= repsln 
     465                  zdzrho_raw = ( - zalbet(ji,jj,jk) * zdkt + zbeta0*zdks ) / fse3w(ji,jj,jk+kp) 
     466                  zdzrho(ji,jj,jk,kp) = - MIN( - repsln , zdzrho_raw )    ! force zdzrho >= repsln 
    472467                 END DO 
    473468            END DO 
     
    507502                     ! Add s-coordinate slope at t-points (do this by *subtracting* gradient of depth) 
    508503                     zti_g_raw = (  zdxrho(ji+ip,jj,jk-kp,1-ip) / zdzrho(ji+ip,jj,jk-kp,kp)      & 
    509                         &          - ( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) / e1u(ji,jj)  ) * umask(ji,jj,jk) 
    510                      ze3_e1    =  fse3w(ji+ip,jj,jk-kp) / e1u(ji,jj)  
     504                        &          - ( fsdept(ji+1,jj,jk-kp) - fsdept(ji,jj,jk-kp) ) * r1_e1u(ji,jj)  ) * umask(ji,jj,jk) 
     505                     ze3_e1    =  fse3w(ji+ip,jj,jk-kp) * r1_e1u(ji,jj)  
    511506                     zti_mlb(ji+ip,jj   ,1-ip,kp) = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1  , ABS( zti_g_raw ) ), zti_g_raw ) 
    512507                  ENDIF 
     
    550545                     ! 
    551546                     ! Must mask contribution to slope for triad jk=1,kp=0 that poke up though ocean surface 
    552                      zti_coord = znot_thru_surface * ( fsdept(ji+1,jj  ,jk) - fsdept(ji,jj,jk) ) / e1u(ji,jj) 
    553                      ztj_coord = znot_thru_surface * ( fsdept(ji  ,jj+1,jk) - fsdept(ji,jj,jk) ) / e2v(ji,jj)                  ! unmasked 
     547                     zti_coord = znot_thru_surface * ( fsdept(ji+1,jj  ,jk) - fsdept(ji,jj,jk) ) * r1_e1u(ji,jj) 
     548                     ztj_coord = znot_thru_surface * ( fsdept(ji  ,jj+1,jk) - fsdept(ji,jj,jk) ) * r1_e2v(ji,jj)     ! unmasked 
    554549                     zti_g_raw = zti_raw - zti_coord      ! ref to geopot surfaces 
    555550                     ztj_g_raw = ztj_raw - ztj_coord 
    556551                     ! additional limit required in bilaplacian case 
    557                      ze3_e1    = fse3w(ji+ip,jj   ,jk+kp) / e1u(ji,jj) 
    558                      ze3_e2    = fse3w(ji   ,jj+jp,jk+kp) / e2v(ji,jj) 
     552                     ze3_e1    = fse3w(ji+ip,jj   ,jk+kp) * r1_e1u(ji,jj) 
     553                     ze3_e2    = fse3w(ji   ,jj+jp,jk+kp) * r1_e2v(ji,jj) 
    559554                     ! NB: hard coded factor 5 (can be a namelist parameter...) 
    560555                     zti_g_lim = SIGN( MIN( rn_slpmax, 5.0_wp * ze3_e1, ABS( zti_g_raw ) ), zti_g_raw ) 
     
    602597#endif 
    603598                     ! 
    604                      zbu  = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji   ,jj   ,jk   ) 
    605                      zbv  = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji   ,jj   ,jk   ) 
    606                      zbti = e1e2t(ji+ip,jj   )      * fse3w(ji+ip,jj   ,jk+kp) 
    607                      zbtj = e1e2t(ji   ,jj+jp)      * fse3w(ji   ,jj+jp,jk+kp) 
     599                     zbu  = e1e2u(ji   ,jj   ) * fse3u(ji   ,jj   ,jk   ) 
     600                     zbv  = e1e2v(ji   ,jj   ) * fse3v(ji   ,jj   ,jk   ) 
     601                     zbti = e1e2t(ji+ip,jj   ) * fse3w(ji+ip,jj   ,jk+kp) 
     602                     zbtj = e1e2t(ji   ,jj+jp) * fse3w(ji   ,jj+jp,jk+kp) 
    608603                     ! 
    609604                     !!gm this may inhibit vectorization on Vect Computers, and even on scalar computers....  ==> to be checked 
     
    673668      !                                            !==   surface mixed layer mask   ! 
    674669      DO jk = 1, jpk                               ! =1 inside the mixed layer, =0 otherwise 
    675 # if defined key_vectopt_loop 
    676          DO jj = 1, 1 
    677             DO ji = 1, jpij                        ! vector opt. (forced unrolling) 
    678 # else 
    679670         DO jj = 1, jpj 
    680671            DO ji = 1, jpi 
    681 # endif 
    682672               ik = nmln(ji,jj) - 1 
    683673               IF( jk <= ik ) THEN   ;   omlmask(ji,jj,jk) = 1._wp 
     
    699689      !----------------------------------------------------------------------- 
    700690      ! 
    701 # if defined key_vectopt_loop 
    702       DO jj = 1, 1 
    703          DO ji = jpi+2, jpij-jpi-1   ! vector opt. (forced unrolling) 
    704 # else 
    705691      DO jj = 2, jpjm1 
    706692         DO ji = 2, jpim1 
    707 # endif 
    708693            !                        !==   Slope at u- & v-points just below the Mixed Layer   ==! 
    709694            ! 
     
    714699            zbv = 0.5_wp * ( p_dzr(ji,jj,ikv) + p_dzr(ji  ,jj+1,ikv) ) 
    715700            !                        !- horizontal density gradient at u- & v-points 
    716             zau = p_gru(ji,jj,iku) / e1u(ji,jj) 
    717             zav = p_grv(ji,jj,ikv) / e2v(ji,jj) 
     701            zau = p_gru(ji,jj,iku) * r1_e1u(ji,jj) 
     702            zav = p_grv(ji,jj,ikv) * r1_e2v(ji,jj) 
    718703            !                        !- bound the slopes: abs(zw.)<= 1/100 and zb..<0 
    719704            !                           kxz max= ah slope max =< e1 e3 /(pi**2 2 dt) 
     
    816801!               DO jj = 2, jpjm1 
    817802!                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    818 !                     uslp (ji,jj,jk) = -1./e1u(ji,jj) * ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * umask(ji,jj,jk) 
    819 !                     vslp (ji,jj,jk) = -1./e2v(ji,jj) * ( fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * vmask(ji,jj,jk) 
    820 !                     wslpi(ji,jj,jk) = -1./e1t(ji,jj) * ( fsdepw(ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * tmask(ji,jj,jk) * 0.5 
    821 !                     wslpj(ji,jj,jk) = -1./e2t(ji,jj) * ( fsdepw(ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * tmask(ji,jj,jk) * 0.5 
     803!                     uslp (ji,jj,jk) = - ( fsdept(ji+1,jj,jk) - fsdept(ji ,jj ,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 
     804!                     vslp (ji,jj,jk) = - ( fsdept(ji,jj+1,jk) - fsdept(ji ,jj ,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk) 
     805!                     wslpi(ji,jj,jk) = - ( fsdepw(ji+1,jj,jk) - fsdepw(ji-1,jj,jk) ) * r1_e1t(ji,jj) * tmask(ji,jj,jk) * 0.5 
     806!                     wslpj(ji,jj,jk) = - ( fsdepw(ji,jj+1,jk) - fsdepw(ji,jj-1,jk) ) * r1_e2t(ji,jj) * tmask(ji,jj,jk) * 0.5 
    822807!                  END DO 
    823808!               END DO 
Note: See TracChangeset for help on using the changeset viewer.