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.
Diff from branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90@8627 to trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90@4990 – NEMO

Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilap.F90

    r8627 r4990  
    1818   USE oce             ! ocean dynamics and tracers 
    1919   USE dom_oce         ! ocean space and time domain 
    20    USE phycst          ! physical constants 
    2120   USE ldfdyn_oce      ! ocean dynamics: lateral physics 
    2221   ! 
    2322   USE in_out_manager  ! I/O manager 
    24    USE iom             ! I/O library 
    2523   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2624   USE wrk_nemo        ! Memory Allocation 
     
    7775      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    7876      ! 
    79       INTEGER  ::   ji, jj, jk                       ! dummy loop indices 
    80       REAL(wp) ::   zua, zva, zbt, ze2u, ze2v, zzz   ! local scalar 
    81       REAL(wp), POINTER, DIMENSION(:,:  ) ::   zcu, zcv 
    82       REAL(wp), POINTER, DIMENSION(:,:,:) ::   zuf, zut, zlu, zlv 
    83       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   z2d   ! 2D workspace  
     77      INTEGER  ::   ji, jj, jk                  ! dummy loop indices 
     78      REAL(wp) ::   zua, zva, zbt, ze2u, ze2v   ! temporary scalar 
     79      REAL(wp), POINTER, DIMENSION(:,:  ) :: zcu, zcv 
     80      REAL(wp), POINTER, DIMENSION(:,:,:) :: zuf, zut, zlu, zlv 
    8481      !!---------------------------------------------------------------------- 
    8582      ! 
     
    115112            DO jj = 2, jpjm1 
    116113               DO ji = fs_2, fs_jpim1   ! vector opt. 
    117                   zlu(ji,jj,jk) = - (   zuf(ji  ,jj,jk) -  zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   & 
    118                      &            + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) /  e1u(ji,jj) 
    119     
    120                   zlv(ji,jj,jk) = + (   zuf(ji,jj  ,jk) -  zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   & 
    121                      &            + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) /  e2v(ji,jj) 
     114                  zlu(ji,jj,jk) = - ( zuf(ji,jj,jk) - zuf(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   & 
     115                     &         + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj,jk) ) / e1u(ji,jj) 
     116    
     117                  zlv(ji,jj,jk) = + ( zuf(ji,jj,jk) - zuf(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   & 
     118                     &         + ( hdivb(ji,jj+1,jk) - hdivb(ji,jj,jk) ) / e2v(ji,jj) 
    122119               END DO 
    123120            END DO 
     
    126123               DO ji = fs_2, fs_jpim1   ! vector opt. 
    127124                  zlu(ji,jj,jk) = - ( rotb (ji  ,jj,jk) - rotb (ji,jj-1,jk) ) / e2u(ji,jj)   & 
    128                      &            + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) / e1u(ji,jj) 
     125                     &         + ( hdivb(ji+1,jj,jk) - hdivb(ji,jj  ,jk) ) / e1u(ji,jj) 
    129126    
    130127                  zlv(ji,jj,jk) = + ( rotb (ji,jj  ,jk) - rotb (ji-1,jj,jk) ) / e1v(ji,jj)   & 
    131                      &            + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) / e2v(ji,jj) 
     128                     &         + ( hdivb(ji,jj+1,jk) - hdivb(ji  ,jj,jk) ) / e2v(ji,jj) 
    132129               END DO   
    133130            END DO   
     
    136133      CALL lbc_lnk( zlu, 'U', -1. )   ;   CALL lbc_lnk( zlv, 'V', -1. )   ! Boundary conditions 
    137134 
    138       IF( iom_use('dispkexyfo') ) THEN   ! ocean kinetic energy dissipation per unit area 
    139          !                               ! due to xy friction (xy=lateral)  
    140          !   see NEMO_book appendix C, §C.7.2 (N.B. here averaged at t-points) 
    141          !   local dissipation of KE at t-point due to bilaplacian operator is given by : 
    142          !      + ahmu mi( zlu**2 ) + mj( ahmv zlv**2 ) 
    143          !   Note that a sign + is used as in v3.6 ahm is negative for bilaplacian viscosity 
    144          ! 
    145          ! NB: ahm is negative when bilaplacian is used 
    146          ALLOCATE( z2d(jpi,jpj) ) 
    147          z2d(:,:) = 0._wp 
    148          DO jk = 1, jpkm1 
    149             DO jj = 2, jpjm1 
    150                DO ji = 2, jpim1 
    151                   z2d(ji,jj) = z2d(ji,jj)                                                                  & 
    152                      &  +  (  fsahmu(ji,jj,jk)*zlu(ji,jj,jk)**2 + fsahmu(ji-1,jj,jk)*zlu(ji-1,jj,jk)**2    & 
    153                      &      + fsahmv(ji,jj,jk)*zlv(ji,jj,jk)**2 + fsahmv(ji,jj-1,jk)*zlv(ji,jj-1,jk)**2  ) * tmask(ji,jj,jk) 
    154                END DO 
    155             END DO 
    156          END DO 
    157          zzz = 0.5_wp * rau0 
    158          z2d(:,:) = zzz * z2d(:,:)  
    159          CALL lbc_lnk( z2d,'T', 1. ) 
    160          CALL iom_put( 'dispkexyfo', z2d ) 
    161          DEALLOCATE( z2d ) 
    162       ENDIF 
    163        
    164     
    165       ! Third derivative 
    166       ! ---------------- 
    167       ! 
     135          
    168136      DO jk = 1, jpkm1 
    169          ! 
     137    
     138         ! Third derivative 
     139         ! ---------------- 
     140          
    170141         ! Multiply by the eddy viscosity coef. (at u- and v-points) 
    171          zlu(:,:,jk) = zlu(:,:,jk) * fsahmu(:,:,jk)  
    172          zlv(:,:,jk) = zlv(:,:,jk) * fsahmv(:,:,jk) 
    173          !          
     142         zlu(:,:,jk) = zlu(:,:,jk) * ( fsahmu(:,:,jk) * (1-nkahm_smag) + nkahm_smag) 
     143 
     144         zlv(:,:,jk) = zlv(:,:,jk) * ( fsahmv(:,:,jk) * (1-nkahm_smag) + nkahm_smag) 
     145          
    174146         ! Contravariant "laplacian" 
    175147         zcu(:,:) = e1u(:,:) * zlu(:,:,jk) 
     
    180152            DO ji = 1, fs_jpim1   ! vector opt. 
    181153               zuf(ji,jj,jk) = fmask(ji,jj,jk) * (  zcv(ji+1,jj  ) - zcv(ji,jj)      & 
    182                   &                               - zcu(ji  ,jj+1) + zcu(ji,jj)  )   & 
    183                   &       * fse3f(ji,jj,jk) * r1_e12f(ji,jj) 
     154                  &                            - zcu(ji  ,jj+1) + zcu(ji,jj)  )   & 
     155                  &       * fse3f(ji,jj,jk) / ( e1f(ji,jj)*e2f(ji,jj) ) 
    184156            END DO   
    185157         END DO   
     
    203175      END DO 
    204176 
     177 
    205178      ! boundary conditions on the laplacian curl and div (zuf,zut) 
    206179!!bug gm no need to do this 2 following lbc... 
     
    208181      CALL lbc_lnk( zut, 'T', 1. ) 
    209182 
    210       DO jk = 1, jpkm1               ! Bilaplacian 
     183      DO jk = 1, jpkm1       
     184    
     185         ! Bilaplacian 
     186         ! ----------- 
     187 
    211188         DO jj = 2, jpjm1 
    212189            DO ji = fs_2, fs_jpim1   ! vector opt. 
     
    220197                  &  + ( zut(ji,jj+1,jk) - zut(ji  ,jj,jk) ) / e2v(ji,jj) 
    221198               ! add it to the general momentum trends 
    222                ua(ji,jj,jk) = ua(ji,jj,jk) + zua 
    223                va(ji,jj,jk) = va(ji,jj,jk) + zva 
     199               ua(ji,jj,jk) = ua(ji,jj,jk) + zua * ( fsahmu(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag )) 
     200               va(ji,jj,jk) = va(ji,jj,jk) + zva * ( fsahmv(ji,jj,jk)*nkahm_smag +(1 -nkahm_smag )) 
    224201            END DO 
    225202         END DO 
     203 
    226204         !                                             ! =============== 
    227205      END DO                                           !   End of slab 
Note: See TracChangeset for help on using the changeset viewer.