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 5737 for branches/2015/dev_r5721_CNRS9_NOC3_LDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90 – NEMO

Ignore:
Timestamp:
2015-09-13T09:42:41+02:00 (9 years ago)
Author:
gm
Message:

#1593: LDF-ADV, step I: Phasing of horizontal scale factors correct 2

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5721_CNRS9_NOC3_LDF/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90

    r5029 r5737  
    213213      !! References : Sadourny, r., 1975, j. atmos. sciences, 32, 680-689. 
    214214      !!---------------------------------------------------------------------- 
    215       ! 
    216215      INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    217216      INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
     
    237236      zfact2 = 0.5 * 0.5      ! Local constant initialization 
    238237 
    239 !CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
    240238      !                                                ! =============== 
    241239      DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    252250                  zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    253251                       &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    254                        &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) 
     252                       &     * 0.5 * r1_e1e2f(ji,jj) 
    255253               END DO 
    256254            END DO 
     
    262260                       &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    263261                       &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    264                        &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                               & 
     262                       &       * 0.5 * r1_e1e2f(ji,jj)                                              & 
    265263                       &       ) 
    266264               END DO 
     
    285283               zx1 = zwx(ji-1,jj) + zwx(ji-1,jj+1) 
    286284               zx2 = zwx(ji  ,jj) + zwx(ji  ,jj+1) 
    287                pua(ji,jj,jk) = pua(ji,jj,jk) + zfact2 / e1u(ji,jj) * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
    288                pva(ji,jj,jk) = pva(ji,jj,jk) - zfact2 / e2v(ji,jj) * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 )  
     285               pua(ji,jj,jk) = pua(ji,jj,jk) + zfact2 * r1_e1u(ji,jj) * ( zwz(ji  ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 
     286               pva(ji,jj,jk) = pva(ji,jj,jk) - zfact2 * r1_e2v(ji,jj) * ( zwz(ji-1,jj  ) * zx1 + zwz(ji,jj) * zx2 )  
    289287            END DO   
    290288         END DO   
     
    365363                     zww(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    366364                        &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    367                         &       * 0.5 / ( e1f(ji,jj) * e2f (ji,jj) * fse3f(ji,jj,jk) ) 
     365                        &       * 0.5 / ( e1e2f (ji,jj) * fse3f(ji,jj,jk) ) 
    368366                  END DO 
    369367               END DO 
     
    380378                     zww(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    381379                        &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    382                         &       * 0.5 / ( e1f(ji,jj) * e2f (ji,jj) ) 
     380                        &       * 0.5 * r1_e1e2f(ji,jj) 
    383381                  END DO 
    384382               END DO 
     
    393391         DO jj = 2, jpjm1 
    394392            DO ji = fs_2, fs_jpim1   ! vector opt. 
    395                zy1 = ( zwy(ji,jj-1) + zwy(ji+1,jj-1) ) / e1u(ji,jj) 
    396                zy2 = ( zwy(ji,jj  ) + zwy(ji+1,jj  ) ) / e1u(ji,jj) 
    397                zx1 = ( zwx(ji-1,jj) + zwx(ji-1,jj+1) ) / e2v(ji,jj) 
    398                zx2 = ( zwx(ji  ,jj) + zwx(ji  ,jj+1) ) / e2v(ji,jj) 
     393               zy1 = ( zwy(ji,jj-1) + zwy(ji+1,jj-1) ) * r1_e1u(ji,jj) 
     394               zy2 = ( zwy(ji,jj  ) + zwy(ji+1,jj  ) ) * r1_e1u(ji,jj) 
     395               zx1 = ( zwx(ji-1,jj) + zwx(ji-1,jj+1) ) * r1_e2v(ji,jj) 
     396               zx2 = ( zwx(ji  ,jj) + zwx(ji  ,jj+1) ) * r1_e2v(ji,jj) 
    399397               ! enstrophy conserving formulation for relative vorticity term 
    400398               zua = zfact1 * ( zww(ji  ,jj-1) + zww(ji,jj) ) * ( zy1 + zy2 ) 
     
    481479                  zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    482480                       &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    483                        &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) 
     481                       &     * 0.5 * r1_e1e2f(ji,jj) 
    484482               END DO 
    485483            END DO 
     
    491489                       &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    492490                       &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    493                        &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
     491                       &       * 0.5 * r1_e1e2f(ji,jj)                                                & 
    494492                       &       ) 
    495493               END DO 
     
    497495         END SELECT 
    498496         ! 
    499          IF( ln_sco ) THEN 
    500             DO jj = 1, jpj                      ! caution: don't use (:,:) for this loop  
    501                DO ji = 1, jpi                   ! it causes optimization problems on NEC in auto-tasking 
    502                   zwz(ji,jj) = zwz(ji,jj) / fse3f(ji,jj,jk) 
    503                   zwx(ji,jj) = e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk) 
    504                   zwy(ji,jj) = e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk) 
    505                END DO 
    506             END DO 
     497         IF( ln_sco ) THEN                   !==  horizontal fluxes  ==! 
     498            zwz(:,:) = zwz(:,:) / fse3f(:,:,jk) 
     499            zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk) 
     500            zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
    507501         ELSE 
    508             DO jj = 1, jpj                      ! caution: don't use (:,:) for this loop  
    509                DO ji = 1, jpi                   ! it causes optimization problems on NEC in auto-tasking 
    510                   zwx(ji,jj) = e2u(ji,jj) * un(ji,jj,jk) 
    511                   zwy(ji,jj) = e1v(ji,jj) * vn(ji,jj,jk) 
    512                END DO 
    513             END DO 
    514          ENDIF 
    515          ! 
    516          ! Compute and add the vorticity term trend 
    517          ! ---------------------------------------- 
     502            zwx(:,:) = e2u(:,:) * un(:,:,jk) 
     503            zwy(:,:) = e1v(:,:) * vn(:,:,jk) 
     504         ENDIF 
     505         !                                   !==  compute and add the vorticity term trend  =! 
    518506         DO jj = 2, jpjm1 
    519507            DO ji = fs_2, fs_jpim1   ! vector opt. 
    520                zuav = zfact1 / e1u(ji,jj) * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1)   & 
    521                   &                         + zwy(ji  ,jj  ) + zwy(ji+1,jj  ) ) 
    522                zvau =-zfact1 / e2v(ji,jj) * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1)   & 
    523                   &                         + zwx(ji  ,jj  ) + zwx(ji  ,jj+1) ) 
     508               zuav = zfact1 * r1_e1u(ji,jj) * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1)   & 
     509                  &                            + zwy(ji  ,jj  ) + zwy(ji+1,jj  ) ) 
     510               zvau =-zfact1 * r1_e2v(ji,jj) * ( zwx(ji-1,jj  ) + zwx(ji-1,jj+1)   & 
     511                  &                            + zwx(ji  ,jj  ) + zwx(ji  ,jj+1) ) 
    524512               pua(ji,jj,jk) = pua(ji,jj,jk) + zuav * ( zwz(ji  ,jj-1) + zwz(ji,jj) ) 
    525513               pva(ji,jj,jk) = pva(ji,jj,jk) + zvau * ( zwz(ji-1,jj  ) + zwz(ji,jj) ) 
     
    553541      !! References : Arakawa and Lamb 1980, Mon. Wea. Rev., 109, 18-36 
    554542      !!---------------------------------------------------------------------- 
    555       ! 
    556543      INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    557544      INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
     
    627614      zfac12 = 1._wp / 12._wp    ! Local constant initialization 
    628615 
    629        
    630 !CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz, ztnw, ztne, ztsw, ztse ) 
    631616      !                                                ! =============== 
    632617      DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    645630                  zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    646631                       &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    647                        &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) * ze3f(ji,jj,jk) 
     632                       &     * 0.5 * r1_e1e2f(ji,jj) * ze3f(ji,jj,jk) 
    648633               END DO 
    649634            END DO 
     
    657642                       &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    658643                       &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    659                        &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
    660                        &       ) * ze3f(ji,jj,jk) 
     644                       &       * 0.5 * r1_e1e2f(ji,jj)   ) * ze3f(ji,jj,jk) 
    661645               END DO 
    662646            END DO 
    663647            CALL lbc_lnk( zwz, 'F', 1. ) 
    664648         END SELECT 
    665  
     649         ! 
     650         !                                   !==  horizontal fluxes  ==! 
    666651         zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk) 
    667652         zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
    668653 
    669          ! Compute and add the vorticity term trend 
    670          ! ---------------------------------------- 
     654         !                                   !==  compute and add the vorticity term trend  =! 
    671655         jj = 2 
    672656         ztne(1,:) = 0   ;   ztnw(1,:) = 0   ;   ztse(1,:) = 0   ;   ztsw(1,:) = 0 
    673          DO ji = 2, jpi    
     657         DO ji = 2, jpi          ! split in 2 parts due to vector opt. 
    674658               ztne(ji,jj) = zwz(ji-1,jj  ) + zwz(ji  ,jj  ) + zwz(ji  ,jj-1) 
    675659               ztnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj  ) + zwz(ji  ,jj  ) 
     
    687671         DO jj = 2, jpjm1 
    688672            DO ji = fs_2, fs_jpim1   ! vector opt. 
    689                zua = + zfac12 / e1u(ji,jj) * (  ztne(ji,jj  ) * zwy(ji  ,jj  ) + ztnw(ji+1,jj) * zwy(ji+1,jj  )   & 
     673               zua = + zfac12 * r1_e1u(ji,jj) * (  ztne(ji,jj  ) * zwy(ji  ,jj  ) + ztnw(ji+1,jj) * zwy(ji+1,jj  )   & 
    690674                  &                           + ztse(ji,jj  ) * zwy(ji  ,jj-1) + ztsw(ji+1,jj) * zwy(ji+1,jj-1) ) 
    691                zva = - zfac12 / e2v(ji,jj) * (  ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji  ,jj+1)   & 
     675               zva = - zfac12 * r1_e2v(ji,jj) * (  ztsw(ji,jj+1) * zwx(ji-1,jj+1) + ztse(ji,jj+1) * zwx(ji  ,jj+1)   & 
    692676                  &                           + ztnw(ji,jj  ) * zwx(ji-1,jj  ) + ztne(ji,jj  ) * zwx(ji  ,jj  ) ) 
    693677               pua(ji,jj,jk) = pua(ji,jj,jk) + zua 
Note: See TracChangeset for help on using the changeset viewer.