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 12590 for NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/TRA/traldf_iso.F90 – NEMO

Ignore:
Timestamp:
2020-03-23T22:16:19+01:00 (4 years ago)
Author:
techene
Message:

all: add e3 substitute, OCE/DOM/domzgr_substitute.h90: correct a bug for e3f

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/TRA/traldf_iso.F90

    r12377 r12590  
    1515   !!---------------------------------------------------------------------- 
    1616   !!   tra_ldf_iso   : update the tracer trend with the horizontal component of a iso-neutral laplacian operator 
    17    !!                   and with the vertical part of the isopycnal or geopotential s-coord. operator  
     17   !!                   and with the vertical part of the isopycnal or geopotential s-coord. operator 
    1818   !!---------------------------------------------------------------------- 
    1919   USE oce            ! ocean dynamics and active tracers 
     
    4141   !! * Substitutions 
    4242#  include "do_loop_substitute.h90" 
     43#  include "domzgr_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    5455      !!                  ***  ROUTINE tra_ldf_iso  *** 
    5556      !! 
    56       !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive  
    57       !!      trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and  
     57      !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive 
     58      !!      trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and 
    5859      !!      add it to the general trend of tracer equation. 
    5960      !! 
    60       !! ** Method  :   The horizontal component of the lateral diffusive trends  
     61      !! ** Method  :   The horizontal component of the lateral diffusive trends 
    6162      !!      is provided by a 2nd order operator rotated along neural or geopo- 
    6263      !!      tential surfaces to which an eddy induced advection can be added 
     
    6970      !! 
    7071      !!      2nd part :  horizontal fluxes of the lateral mixing operator 
    71       !!      ========     
     72      !!      ======== 
    7273      !!         zftu =  pahu e2u*e3u/e1u di[ tb ] 
    7374      !!               - pahu e2u*uslp    dk[ mi(mk(tb)) ] 
     
    111112      REAL(wp) ::  zcoef0, ze3w_2, zsign, z2dt, z1_2dt   !   -      - 
    112113      REAL(wp), DIMENSION(jpi,jpj)     ::   zdkt, zdk1t, z2d 
    113       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zdit, zdjt, zftu, zftv, ztfw  
     114      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zdit, zdjt, zftu, zftv, ztfw 
    114115      !!---------------------------------------------------------------------- 
    115116      ! 
     
    119120         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 
    120121         ! 
    121          akz     (:,:,:) = 0._wp       
     122         akz     (:,:,:) = 0._wp 
    122123         ah_wslp2(:,:,:) = 0._wp 
    123124      ENDIF 
    124       !    
     125      ! 
    125126      l_hst = .FALSE. 
    126127      l_ptr = .FALSE. 
    127       IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE.  
     128      IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE. 
    128129      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    129130         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE. 
     
    138139      ELSE                    ;   zsign = -1._wp 
    139140      ENDIF 
    140           
     141 
    141142      !!---------------------------------------------------------------------- 
    142143      !!   0 - calculate  ah_wslp2 and akz 
     
    173174               DO_3D_10_10( 2, jpkm1 ) 
    174175                  akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   & 
    175                      &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
     176                     &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk)    & 
     177                     &                        / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
    176178               END_3D 
    177179            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator 
     
    184186           ! 
    185187         ELSE                                    ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit 
    186             akz(:,:,:) = ah_wslp2(:,:,:)       
     188            akz(:,:,:) = ah_wslp2(:,:,:) 
    187189         ENDIF 
    188190      ENDIF 
     
    191193      DO jn = 1, kjpt                                            ! tracer loop 
    192194         !                                                       ! =========== 
    193          !                                                
    194          !!---------------------------------------------------------------------- 
    195          !!   I - masked horizontal derivative  
     195         ! 
     196         !!---------------------------------------------------------------------- 
     197         !!   I - masked horizontal derivative 
    196198         !!---------------------------------------------------------------------- 
    197199!!gm : bug.... why (x,:,:)?   (1,jpj,:) and (jpi,1,:) should be sufficient.... 
     
    200202         !!end 
    201203 
    202          ! Horizontal tracer gradient  
     204         ! Horizontal tracer gradient 
    203205         DO_3D_10_10( 1, jpkm1 ) 
    204206            zdit(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
     
    207209         IF( ln_zps ) THEN      ! botton and surface ocean correction of the horizontal gradient 
    208210            DO_2D_10_10 
    209                zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)           
     211               zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 
    210212               zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 
    211213            END_2D 
    212214            IF( ln_isfcav ) THEN      ! first wet level beneath a cavity 
    213215               DO_2D_10_10 
    214                   IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)           
    215                   IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)      
     216                  IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn) 
     217                  IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn) 
    216218               END_2D 
    217219            ENDIF 
     
    248250               zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   & 
    249251                  &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      & 
    250                   &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk)                   
     252                  &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk) 
    251253            END_2D 
    252254            ! 
     
    254256               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      & 
    255257                  &                                                 + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
    256                   &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     258                  &                                             * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    257259            END_2D 
    258          END DO                                        !   End of slab   
     260         END DO                                        !   End of slab 
    259261 
    260262         !!---------------------------------------------------------------------- 
     
    266268         !                          ! Surface and bottom vertical fluxes set to zero 
    267269         ztfw(:,:, 1 ) = 0._wp      ;      ztfw(:,:,jpk) = 0._wp 
    268           
     270 
    269271         DO_3D_00_00( 2, jpkm1 ) 
    270272            ! 
     
    295297            END_3D 
    296298            ! 
    297          ELSE                                   ! bilaplacian  
     299         ELSE                                   ! bilaplacian 
    298300            SELECT CASE( kpass ) 
    299301            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2 
     
    311313            END SELECT 
    312314         ENDIF 
    313          !          
     315         ! 
    314316         DO_3D_00_00( 1, jpkm1 ) 
    315317            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   & 
Note: See TracChangeset for help on using the changeset viewer.