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.
traadv_cen.F90 in NEMO/trunk/src/OCE/TRA – NEMO

source: NEMO/trunk/src/OCE/TRA/traadv_cen.F90 @ 14587

Last change on this file since 14587 was 14433, checked in by smasson, 3 years ago

trunk: merge dev_r14312_MPI_Interface into the trunk, #2598
  • Property svn:keywords set to Id
File size: 10.3 KB
RevLine 
[5770]1MODULE traadv_cen
2   !!======================================================================
3   !!                     ***  MODULE  traadv_cen  ***
[6140]4   !! Ocean  tracers:   advective trend (2nd/4th order centered)
[5770]5   !!======================================================================
6   !! History :  3.7  ! 2014-05  (G. Madec)  original code
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
[6140]10   !!   tra_adv_cen   : update the tracer trend with the advection trends using a centered or scheme (2nd or 4th order)
11   !!                   NB: on the vertical it is actually a 4th order COMPACT scheme which is used
[5770]12   !!----------------------------------------------------------------------
[6140]13   USE dom_oce        ! ocean space and time domain
14   USE eosbn2         ! equation of state
[14072]15   USE traadv_fct     ! acces to routine interp_4th_cpt
[6140]16   USE trd_oce        ! trends: ocean variables
[14072]17   USE trdtra         ! trends manager: tracers
[6140]18   USE diaptr         ! poleward transport diagnostics
[7646]19   USE diaar5         ! AR5 diagnostics
[5770]20   !
[6140]21   USE in_out_manager ! I/O manager
22   USE iom            ! IOM library
23   USE trc_oce        ! share passive tracers/Ocean variables
24   USE lib_mpp        ! MPP library
[5770]25
26   IMPLICIT NONE
27   PRIVATE
28
[9019]29   PUBLIC   tra_adv_cen   ! called by traadv.F90
[14072]30
[5770]31   REAL(wp) ::   r1_6 = 1._wp / 6._wp   ! =1/6
32
[9019]33   LOGICAL ::   l_trd   ! flag to compute trends
34   LOGICAL ::   l_ptr   ! flag to compute poleward transport
35   LOGICAL ::   l_hst   ! flag to compute heat/salt transport
[7646]36
[5770]37   !! * Substitutions
[12377]38#  include "do_loop_substitute.h90"
[13237]39#  include "domzgr_substitute.h90"
[5770]40   !!----------------------------------------------------------------------
[9598]41   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
[10068]42   !! $Id$
43   !! Software governed by the CeCILL license (see ./LICENSE)
[5770]44   !!----------------------------------------------------------------------
45CONTAINS
46
[12377]47   SUBROUTINE tra_adv_cen( kt, kit000, cdtype, pU, pV, pW,     &
[14072]48      &                    Kmm, pt, kjpt, Krhs, kn_cen_h, kn_cen_v )
[5770]49      !!----------------------------------------------------------------------
50      !!                  ***  ROUTINE tra_adv_cen  ***
[14072]51      !!
[5770]52      !! ** Purpose :   Compute the now trend due to the advection of tracers
53      !!      and add it to the general trend of passive tracer equations.
54      !!
55      !! ** Method  :   The advection is evaluated by a 2nd or 4th order scheme
[14072]56      !!               using now fields (leap-frog scheme).
[5770]57      !!       kn_cen_h = 2  ==>> 2nd order centered scheme on the horizontal
58      !!                = 4  ==>> 4th order    -        -       -      -
59      !!       kn_cen_v = 2  ==>> 2nd order centered scheme on the vertical
60      !!                = 4  ==>> 4th order COMPACT  scheme     -      -
61      !!
[12377]62      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends
[6140]63      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T)
[12377]64      !!             - poleward advective heat and salt transport (l_diaptr=T)
[5770]65      !!----------------------------------------------------------------------
[12377]66      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index
67      INTEGER                                  , INTENT(in   ) ::   Kmm, Krhs       ! ocean time level indices
68      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index
69      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator)
70      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers
71      INTEGER                                  , INTENT(in   ) ::   kn_cen_h        ! =2/4 (2nd or 4th order scheme)
72      INTEGER                                  , INTENT(in   ) ::   kn_cen_v        ! =2/4 (2nd or 4th order scheme)
[13982]73      ! TEMP: [tiling] This can be A2D(nn_hls) if using XIOS (subdomain support)
[12377]74      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pU, pV, pW      ! 3 ocean volume flux components
75      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! tracers and RHS of tracer equation
[5770]76      !
77      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
78      INTEGER  ::   ierr             ! local integer
79      REAL(wp) ::   zC2t_u, zC4t_u   ! local scalars
80      REAL(wp) ::   zC2t_v, zC4t_v   !   -      -
[13982]81      REAL(wp), DIMENSION(A2D(nn_hls),jpk) ::   zwx, zwy, zwz, ztu, ztv, ztw
[5770]82      !!----------------------------------------------------------------------
83      !
[13982]84      IF( ntile == 0 .OR. ntile == 1 )  THEN                       ! Do only on the first tile
85         IF( kt == kit000 )  THEN
86            IF(lwp) WRITE(numout,*)
87            IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v
88            IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ '
89         ENDIF
90         !                          ! set local switches
91         l_trd = .FALSE.
92         l_hst = .FALSE.
93         l_ptr = .FALSE.
94         IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )       l_trd = .TRUE.
95         IF(   cdtype == 'TRA' .AND. ( iom_use( 'sophtadv' ) .OR. iom_use( 'sophtadv' ) )  )    l_ptr = .TRUE.
96         IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
97            &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE.
[5770]98      ENDIF
[7646]99      !
[14072]100      !
[6140]101      zwz(:,:, 1 ) = 0._wp       ! surface & bottom vertical flux set to zero for all tracers
102      zwz(:,:,jpk) = 0._wp
[5770]103      !
104      DO jn = 1, kjpt            !==  loop over the tracers  ==!
105         !
[6140]106         SELECT CASE( kn_cen_h )       !--  Horizontal fluxes  --!
[5770]107         !
[6140]108         CASE(  2  )                         !* 2nd order centered
[13295]109            DO_3D( 1, 0, 1, 0, 1, jpkm1 )
[12377]110               zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm) )
111               zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) )
112            END_3D
[5770]113            !
[6140]114         CASE(  4  )                         !* 4th order centered
115            ztu(:,:,jpk) = 0._wp                   ! Bottom value : flux set to zero
[5770]116            ztv(:,:,jpk) = 0._wp
[13982]117            DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )          ! masked gradient
[12377]118               ztu(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk)
119               ztv(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk)
120            END_3D
[14433]121            IF (nn_hls.EQ.1) CALL lbc_lnk( 'traadv_cen', ztu, 'U', -1.0_wp , ztv, 'V', -1.0_wp )   ! Lateral boundary cond.
[5770]122            !
[13982]123            DO_3D( nn_hls-1, 0, nn_hls-1, 0, 1, jpkm1 )           ! Horizontal advective fluxes
[12377]124               zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm)   ! C2 interpolation of T at u- & v-points (x2)
125               zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm)
126               !                                                  ! C4 interpolation of T at u- & v-points (x2)
127               zC4t_u =  zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) )
128               zC4t_v =  zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) )
129               !                                                  ! C4 fluxes
130               zwx(ji,jj,jk) =  0.5_wp * pU(ji,jj,jk) * zC4t_u
131               zwy(ji,jj,jk) =  0.5_wp * pV(ji,jj,jk) * zC4t_v
132            END_3D
[14433]133            IF (nn_hls.EQ.1) CALL lbc_lnk( 'traadv_cen', zwx, 'U', -1. , zwy, 'V', -1. )
[5770]134            !
135         CASE DEFAULT
[13457]136            CALL ctl_stop( 'traadv_cen: wrong value for nn_cen' )
[5770]137         END SELECT
138         !
[6140]139         SELECT CASE( kn_cen_v )       !--  Vertical fluxes  --!   (interior)
[5770]140         !
[6140]141         CASE(  2  )                         !* 2nd order centered
[13295]142            DO_3D( 0, 0, 0, 0, 2, jpk )
[12377]143               zwz(ji,jj,jk) = 0.5 * pW(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) * wmask(ji,jj,jk)
144            END_3D
[5770]145            !
[6140]146         CASE(  4  )                         !* 4th order compact
[12377]147            CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw )      ! ztw = interpolated value of T at w-point
[13295]148            DO_3D( 0, 0, 0, 0, 2, jpkm1 )
[12377]149               zwz(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)
150            END_3D
[5770]151            !
152         END SELECT
153         !
[6140]154         IF( ln_linssh ) THEN                !* top value   (linear free surf. only as zwz is multiplied by wmask)
[5770]155            IF( ln_isfcav ) THEN                  ! ice-shelf cavities (top of the ocean)
[13295]156               DO_2D( 1, 1, 1, 1 )
[14072]157                  zwz(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm)
[12377]158               END_2D
[5770]159            ELSE                                   ! no ice-shelf cavities (only ocean surface)
[13982]160               DO_2D( 1, 1, 1, 1 )
161                  zwz(ji,jj,1) = pW(ji,jj,1) * pt(ji,jj,1,jn,Kmm)
162               END_2D
[5770]163            ENDIF
164         ENDIF
[14072]165         !
[13497]166         DO_3D( 0, 0, 0, 0, 1, jpkm1 )   !--  Divergence of advective fluxes  --!
[12377]167            pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs)    &
168               &             - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )    &
169               &                + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  )    &
[13237]170               &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) &
171               &                * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
[12377]172         END_3D
[13497]173         !                               ! trend diagnostics
[7646]174         IF( l_trd ) THEN
[12377]175            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
176            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
177            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
[13982]178         ENDIF
179         !                                 ! "Poleward" heat and salt transports
[9019]180         IF( l_ptr )   CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) )
[7646]181         !                                 !  heat and salt transport
[9019]182         IF( l_hst )   CALL dia_ar5_hst( jn, 'adv', zwx(:,:,:), zwy(:,:,:) )
[5770]183         !
184      END DO
185      !
186   END SUBROUTINE tra_adv_cen
[14072]187
[5770]188   !!======================================================================
189END MODULE traadv_cen
Note: See TracBrowser for help on using the repository browser.