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/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA – NEMO

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90 @ 10946

Last change on this file since 10946 was 10946, checked in by acc, 5 years ago

2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Convert STO, TRD and USR modules and all knock on effects of these conversions. Note change to USR module may have implications for the TEST CASES (not tested yet). Standard SETTE tested only
  • Property svn:keywords set to Id
File size: 10.8 KB
Line 
1MODULE traadv_cen
2   !!======================================================================
3   !!                     ***  MODULE  traadv_cen  ***
4   !! Ocean  tracers:   advective trend (2nd/4th order centered)
5   !!======================================================================
6   !! History :  3.7  ! 2014-05  (G. Madec)  original code
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
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
12   !!----------------------------------------------------------------------
13   USE dom_oce        ! ocean space and time domain
14   USE eosbn2         ! equation of state
15   USE traadv_fct     ! acces to routine interp_4th_cpt
16   USE trd_oce        ! trends: ocean variables
17   USE trdtra         ! trends manager: tracers
18   USE diaptr         ! poleward transport diagnostics
19   USE diaar5         ! AR5 diagnostics
20   !
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
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   tra_adv_cen   ! called by traadv.F90
30   
31   REAL(wp) ::   r1_6 = 1._wp / 6._wp   ! =1/6
32
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
36
37   !! * Substitutions
38#  include "vectopt_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
41   !! $Id$
42   !! Software governed by the CeCILL license (see ./LICENSE)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE tra_adv_cen( kt, kit000, cdtype, pU, pV, pW,     &
47      &                    Kmm, pt, kjpt, Krhs, kn_cen_h, kn_cen_v ) 
48      !!----------------------------------------------------------------------
49      !!                  ***  ROUTINE tra_adv_cen  ***
50      !!                 
51      !! ** Purpose :   Compute the now trend due to the advection of tracers
52      !!      and add it to the general trend of passive tracer equations.
53      !!
54      !! ** Method  :   The advection is evaluated by a 2nd or 4th order scheme
55      !!               using now fields (leap-frog scheme).
56      !!       kn_cen_h = 2  ==>> 2nd order centered scheme on the horizontal
57      !!                = 4  ==>> 4th order    -        -       -      -
58      !!       kn_cen_v = 2  ==>> 2nd order centered scheme on the vertical
59      !!                = 4  ==>> 4th order COMPACT  scheme     -      -
60      !!
61      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends
62      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T)
63      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T)
64      !!----------------------------------------------------------------------
65      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index
66      INTEGER                                  , INTENT(in   ) ::   Kmm, Krhs       ! ocean time level indices
67      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index
68      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator)
69      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers
70      INTEGER                                  , INTENT(in   ) ::   kn_cen_h        ! =2/4 (2nd or 4th order scheme)
71      INTEGER                                  , INTENT(in   ) ::   kn_cen_v        ! =2/4 (2nd or 4th order scheme)
72      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pU, pV, pW      ! 3 ocean volume flux components
73      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! tracers and RHS of tracer equation
74      !
75      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
76      INTEGER  ::   ierr             ! local integer
77      REAL(wp) ::   zC2t_u, zC4t_u   ! local scalars
78      REAL(wp) ::   zC2t_v, zC4t_v   !   -      -
79      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwx, zwy, zwz, ztu, ztv, ztw
80      !!----------------------------------------------------------------------
81      !
82      IF( kt == kit000 )  THEN
83         IF(lwp) WRITE(numout,*)
84         IF(lwp) WRITE(numout,*) 'tra_adv_cen : centered advection scheme on ', cdtype, ' order h/v =', kn_cen_h,'/', kn_cen_v
85         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~ '
86      ENDIF
87      !                          ! set local switches
88      l_trd = .FALSE.
89      l_hst = .FALSE.
90      l_ptr = .FALSE.
91      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )        l_trd = .TRUE.
92      IF(   cdtype == 'TRA' .AND. ln_diaptr )                                                 l_ptr = .TRUE. 
93      IF(   cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
94         &                          iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE.
95      !
96      !                   
97      zwz(:,:, 1 ) = 0._wp       ! surface & bottom vertical flux set to zero for all tracers
98      zwz(:,:,jpk) = 0._wp
99      !
100      DO jn = 1, kjpt            !==  loop over the tracers  ==!
101         !
102         SELECT CASE( kn_cen_h )       !--  Horizontal fluxes  --!
103         !
104         CASE(  2  )                         !* 2nd order centered
105            DO jk = 1, jpkm1
106               DO jj = 1, jpjm1
107                  DO ji = 1, fs_jpim1   ! vector opt.
108                     zwx(ji,jj,jk) = 0.5_wp * pU(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm) )
109                     zwy(ji,jj,jk) = 0.5_wp * pV(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) )
110                  END DO
111               END DO
112            END DO
113            !
114         CASE(  4  )                         !* 4th order centered
115            ztu(:,:,jpk) = 0._wp                   ! Bottom value : flux set to zero
116            ztv(:,:,jpk) = 0._wp
117            DO jk = 1, jpkm1                       ! masked gradient
118               DO jj = 2, jpjm1
119                  DO ji = fs_2, fs_jpim1   ! vector opt.
120                     ztu(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk)
121                     ztv(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk)
122                  END DO
123               END DO
124            END DO
125            CALL lbc_lnk_multi( 'traadv_cen', ztu, 'U', -1. , ztv, 'V', -1. )   ! Lateral boundary cond.
126            !
127            DO jk = 1, jpkm1                       ! Horizontal advective fluxes
128               DO jj = 2, jpjm1
129                  DO ji = 1, fs_jpim1   ! vector opt.
130                     zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm)   ! C2 interpolation of T at u- & v-points (x2)
131                     zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm)
132                     !                                                  ! C4 interpolation of T at u- & v-points (x2)
133                     zC4t_u =  zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) )
134                     zC4t_v =  zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) )
135                     !                                                  ! C4 fluxes
136                     zwx(ji,jj,jk) =  0.5_wp * pU(ji,jj,jk) * zC4t_u
137                     zwy(ji,jj,jk) =  0.5_wp * pV(ji,jj,jk) * zC4t_v
138                  END DO
139               END DO
140            END DO         
141            !
142         CASE DEFAULT
143            CALL ctl_stop( 'traadv_fct: wrong value for nn_fct' )
144         END SELECT
145         !
146         SELECT CASE( kn_cen_v )       !--  Vertical fluxes  --!   (interior)
147         !
148         CASE(  2  )                         !* 2nd order centered
149            DO jk = 2, jpk
150               DO jj = 2, jpjm1
151                  DO ji = fs_2, fs_jpim1   ! vector opt.
152                     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)
153                  END DO
154               END DO
155            END DO
156            !
157         CASE(  4  )                         !* 4th order compact
158            CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw )      ! ztw = interpolated value of T at w-point
159            DO jk = 2, jpkm1
160               DO jj = 2, jpjm1
161                  DO ji = fs_2, fs_jpim1
162                     zwz(ji,jj,jk) = pW(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)
163                  END DO
164               END DO
165            END DO
166            !
167         END SELECT
168         !
169         IF( ln_linssh ) THEN                !* top value   (linear free surf. only as zwz is multiplied by wmask)
170            IF( ln_isfcav ) THEN                  ! ice-shelf cavities (top of the ocean)
171               DO jj = 1, jpj
172                  DO ji = 1, jpi
173                     zwz(ji,jj, mikt(ji,jj) ) = pW(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm) 
174                  END DO
175               END DO   
176            ELSE                                   ! no ice-shelf cavities (only ocean surface)
177               zwz(:,:,1) = pW(:,:,1) * pt(:,:,1,jn,Kmm)
178            ENDIF
179         ENDIF
180         !               
181         DO jk = 1, jpkm1              !--  Divergence of advective fluxes  --!
182            DO jj = 2, jpjm1
183               DO ji = fs_2, fs_jpim1   ! vector opt.
184                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs)    &
185                     &             - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )    &
186                     &                + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  )    &
187                     &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
188               END DO
189            END DO
190         END DO
191         !                             ! trend diagnostics
192         IF( l_trd ) THEN
193            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_xad, zwx, pU, pt(:,:,:,jn,Kmm) )
194            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_yad, zwy, pV, pt(:,:,:,jn,Kmm) )
195            CALL trd_tra( kt, Kmm, Krhs, cdtype, jn, jptra_zad, zwz, pW, pt(:,:,:,jn,Kmm) )
196         END IF
197         !                                 ! "Poleward" heat and salt transports
198         IF( l_ptr )   CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) )
199         !                                 !  heat and salt transport
200         IF( l_hst )   CALL dia_ar5_hst( jn, 'adv', zwx(:,:,:), zwy(:,:,:) )
201         !
202      END DO
203      !
204   END SUBROUTINE tra_adv_cen
205   
206   !!======================================================================
207END MODULE traadv_cen
Note: See TracBrowser for help on using the repository browser.