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.F90 in NEMO/branches/2020/KERNEL-03_Storkey_Coward_RK3_stage2/src/OCE/TRA – NEMO

source: NEMO/branches/2020/KERNEL-03_Storkey_Coward_RK3_stage2/src/OCE/TRA/traadv.F90 @ 12397

Last change on this file since 12397 was 12397, checked in by davestorkey, 4 years ago

2020/KERNEL-03_Storkey_Coward_RK3_stage2 : Consolidation of code to
handle initial Euler timestep in the context of leapfrog
timestepping. This version passes all SETTE tests but fails to bit
compare with the control for several tests (ORCA2_ICE_PISCES, AMM12,
ISOMIP, AGRIF_DEMO, SPITZ12).

  • Property svn:keywords set to Id
File size: 15.4 KB
Line 
1MODULE traadv
2   !!==============================================================================
3   !!                       ***  MODULE  traadv  ***
4   !! Ocean active tracers:  advection trend
5   !!==============================================================================
6   !! History :  2.0  !  2005-11  (G. Madec)  Original code
7   !!            3.3  !  2010-09  (C. Ethe, G. Madec)  merge TRC-TRA + switch from velocity to transport
8   !!            3.6  !  2011-06  (G. Madec)  Addition of Mixed Layer Eddy parameterisation
9   !!            3.7  !  2014-05  (G. Madec)  Add 2nd/4th order cases for CEN and FCT schemes
10   !!             -   !  2014-12  (G. Madec) suppression of cross land advection option
11   !!            3.6  !  2015-06  (E. Clementi) Addition of Stokes drift in case of wave coupling
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   tra_adv       : compute ocean tracer advection trend
16   !!   tra_adv_init  : control the different options of advection scheme
17   !!----------------------------------------------------------------------
18   USE oce            ! ocean dynamics and active tracers
19   USE dom_oce        ! ocean space and time domain
20   USE domvvl         ! variable vertical scale factors
21   USE sbcwave        ! wave module
22   USE sbc_oce        ! surface boundary condition: ocean
23   USE traadv_cen     ! centered scheme            (tra_adv_cen  routine)
24   USE traadv_fct     ! FCT      scheme            (tra_adv_fct  routine)
25   USE traadv_mus     ! MUSCL    scheme            (tra_adv_mus  routine)
26   USE traadv_ubs     ! UBS      scheme            (tra_adv_ubs  routine)
27   USE traadv_qck     ! QUICKEST scheme            (tra_adv_qck  routine)
28   USE tramle         ! Mixed Layer Eddy transport (tra_mle_trp  routine)
29   USE ldftra         ! Eddy Induced transport     (ldf_eiv_trp  routine)
30   USE ldfslp         ! Lateral diffusion: slopes of neutral surfaces
31   USE trd_oce        ! trends: ocean variables
32   USE trdtra         ! trends manager: tracers
33   USE diaptr         ! Poleward heat transport
34   !
35   USE in_out_manager ! I/O manager
36   USE iom            ! I/O module
37   USE prtctl         ! Print control
38   USE lib_mpp        ! MPP library
39   USE timing         ! Timing
40
41   IMPLICIT NONE
42   PRIVATE
43
44   PUBLIC   tra_adv        ! called by step.F90
45   PUBLIC   tra_adv_init   ! called by nemogcm.F90
46
47   !                            !!* Namelist namtra_adv *
48   LOGICAL ::   ln_traadv_OFF    ! no advection on T and S
49   LOGICAL ::   ln_traadv_cen    ! centered scheme flag
50   INTEGER ::      nn_cen_h, nn_cen_v   ! =2/4 : horizontal and vertical choices of the order of CEN scheme
51   LOGICAL ::   ln_traadv_fct    ! FCT scheme flag
52   INTEGER ::      nn_fct_h, nn_fct_v   ! =2/4 : horizontal and vertical choices of the order of FCT scheme
53   LOGICAL ::   ln_traadv_mus    ! MUSCL scheme flag
54   LOGICAL ::      ln_mus_ups           ! use upstream scheme in vivcinity of river mouths
55   LOGICAL ::   ln_traadv_ubs    ! UBS scheme flag
56   INTEGER ::      nn_ubs_v             ! =2/4 : vertical choice of the order of UBS scheme
57   LOGICAL ::   ln_traadv_qck    ! QUICKEST scheme flag
58
59   INTEGER ::   nadv             ! choice of the type of advection scheme
60   !                             ! associated indices:
61   INTEGER, PARAMETER ::   np_NO_adv  = 0   ! no T-S advection
62   INTEGER, PARAMETER ::   np_CEN     = 1   ! 2nd/4th order centered scheme
63   INTEGER, PARAMETER ::   np_FCT     = 2   ! 2nd/4th order Flux Corrected Transport scheme
64   INTEGER, PARAMETER ::   np_MUS     = 3   ! MUSCL scheme
65   INTEGER, PARAMETER ::   np_UBS     = 4   ! 3rd order Upstream Biased Scheme
66   INTEGER, PARAMETER ::   np_QCK     = 5   ! QUICK scheme
67   
68   !!----------------------------------------------------------------------
69   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
70   !! $Id$
71   !! Software governed by the CeCILL license (see ./LICENSE)
72   !!----------------------------------------------------------------------
73CONTAINS
74
75   SUBROUTINE tra_adv( kt, Kbb, Kmm, pts, Krhs )
76      !!----------------------------------------------------------------------
77      !!                  ***  ROUTINE tra_adv  ***
78      !!
79      !! ** Purpose :   compute the ocean tracer advection trend.
80      !!
81      !! ** Method  : - Update (uu(:,:,:,Krhs),vv(:,:,:,Krhs)) with the advection term following nadv
82      !!----------------------------------------------------------------------
83      INTEGER                                  , INTENT(in)    :: kt             ! ocean time-step index
84      INTEGER                                  , INTENT(in)    :: Kbb, Kmm, Krhs ! time level indices
85      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts            ! active tracers and RHS of tracer equation
86      !
87      INTEGER ::   jk   ! dummy loop index
88      REAL(wp), DIMENSION(jpi,jpj,jpk)        :: zuu, zvv, zww   ! 3D workspace
89      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdt, ztrds
90      !!----------------------------------------------------------------------
91      !
92      IF( ln_timing )   CALL timing_start('tra_adv')
93      !
94      !                                         !==  effective transport  ==!
95      zuu(:,:,jpk) = 0._wp
96      zvv(:,:,jpk) = 0._wp
97      zww(:,:,jpk) = 0._wp
98      IF( ln_wave .AND. ln_sdw )  THEN
99         DO jk = 1, jpkm1                                                       ! eulerian transport + Stokes Drift
100            zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) )
101            zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) )
102            zww(:,:,jk) = e1e2t(:,:)                 * ( ww(:,:,jk) + wsd(:,:,jk) )
103         END DO
104      ELSE
105         DO jk = 1, jpkm1
106            zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * uu(:,:,jk,Kmm)               ! eulerian transport only
107            zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * vv(:,:,jk,Kmm)
108            zww(:,:,jk) = e1e2t(:,:)                 * ww(:,:,jk)
109         END DO
110      ENDIF
111      !
112      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN                                ! add z-tilde and/or vvl corrections
113         zuu(:,:,:) = zuu(:,:,:) + un_td(:,:,:)
114         zvv(:,:,:) = zvv(:,:,:) + vn_td(:,:,:)
115      ENDIF
116      !
117      zuu(:,:,jpk) = 0._wp                                                      ! no transport trough the bottom
118      zvv(:,:,jpk) = 0._wp
119      zww(:,:,jpk) = 0._wp
120      !
121      IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   &
122         &              CALL ldf_eiv_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm, Krhs )   ! add the eiv transport (if necessary)
123      !
124      IF( ln_mle    )   CALL tra_mle_trp( kt, nit000, zuu, zvv, zww, 'TRA', Kmm       )   ! add the mle transport (if necessary)
125      !
126      CALL iom_put( "uocetr_eff", zuu )                                        ! output effective transport     
127      CALL iom_put( "vocetr_eff", zvv )
128      CALL iom_put( "wocetr_eff", zww )
129      !
130!!gm ???
131      CALL dia_ptr( kt, Kmm, zvv )                                    ! diagnose the effective MSF
132!!gm ???
133      !
134
135      IF( l_trdtra )   THEN                    !* Save ta and sa trends
136         ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) )
137         ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs)
138         ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs)
139      ENDIF
140      !
141      SELECT CASE ( nadv )                      !==  compute advection trend and add it to general trend  ==!
142      !
143      CASE ( np_CEN )                                 ! Centered scheme : 2nd / 4th order
144         CALL tra_adv_cen    ( kt, nit000, 'TRA',         zuu, zvv, zww, Kmm, pts, jpts, Krhs, nn_cen_h, nn_cen_v )
145      CASE ( np_FCT )                                 ! FCT scheme      : 2nd / 4th order
146         CALL tra_adv_fct    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_fct_h, nn_fct_v )
147      CASE ( np_MUS )                                 ! MUSCL
148         CALL tra_adv_mus    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, ln_mus_ups ) 
149      CASE ( np_UBS )                                 ! UBS
150         CALL tra_adv_ubs    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_ubs_v   )
151      CASE ( np_QCK )                                 ! QUICKEST
152         CALL tra_adv_qck    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs )
153      !
154      END SELECT
155      !
156      IF( l_trdtra )   THEN                      ! save the advective trends for further diagnostics
157         DO jk = 1, jpkm1
158            ztrdt(:,:,jk) = pts(:,:,jk,jp_tem,Krhs) - ztrdt(:,:,jk)
159            ztrds(:,:,jk) = pts(:,:,jk,jp_sal,Krhs) - ztrds(:,:,jk)
160         END DO
161         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_totad, ztrdt )
162         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_totad, ztrds )
163         DEALLOCATE( ztrdt, ztrds )
164      ENDIF
165      !                                              ! print mean trends (used for debugging)
166      IF(sn_cfctl%l_prtctl)   CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' adv  - Ta: ', mask1=tmask,               &
167         &                                  tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
168      !
169      IF( ln_timing )   CALL timing_stop( 'tra_adv' )
170      !
171   END SUBROUTINE tra_adv
172
173
174   SUBROUTINE tra_adv_init
175      !!---------------------------------------------------------------------
176      !!                  ***  ROUTINE tra_adv_init  ***
177      !!               
178      !! ** Purpose :   Control the consistency between namelist options for
179      !!              tracer advection schemes and set nadv
180      !!----------------------------------------------------------------------
181      INTEGER ::   ioptio, ios   ! Local integers
182      !
183      NAMELIST/namtra_adv/ ln_traadv_OFF,                        &   ! No advection
184         &                 ln_traadv_cen , nn_cen_h, nn_cen_v,   &   ! CEN
185         &                 ln_traadv_fct , nn_fct_h, nn_fct_v,   &   ! FCT
186         &                 ln_traadv_mus , ln_mus_ups,           &   ! MUSCL
187         &                 ln_traadv_ubs ,           nn_ubs_v,   &   ! UBS
188         &                 ln_traadv_qck                             ! QCK
189      !!----------------------------------------------------------------------
190      !
191      !                                !==  Namelist  ==!
192      READ  ( numnam_ref, namtra_adv, IOSTAT = ios, ERR = 901)
193901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtra_adv in reference namelist' )
194      !
195      READ  ( numnam_cfg, namtra_adv, IOSTAT = ios, ERR = 902 )
196902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtra_adv in configuration namelist' )
197      IF(lwm) WRITE( numond, namtra_adv )
198      !
199      IF(lwp) THEN                           ! Namelist print
200         WRITE(numout,*)
201         WRITE(numout,*) 'tra_adv_init : choice/control of the tracer advection scheme'
202         WRITE(numout,*) '~~~~~~~~~~~~'
203         WRITE(numout,*) '   Namelist namtra_adv : chose a advection scheme for tracers'
204         WRITE(numout,*) '      No advection on T & S                     ln_traadv_OFF = ', ln_traadv_OFF
205         WRITE(numout,*) '      centered scheme                           ln_traadv_cen = ', ln_traadv_cen
206         WRITE(numout,*) '            horizontal 2nd/4th order               nn_cen_h   = ', nn_fct_h
207         WRITE(numout,*) '            vertical   2nd/4th order               nn_cen_v   = ', nn_fct_v
208         WRITE(numout,*) '      Flux Corrected Transport scheme           ln_traadv_fct = ', ln_traadv_fct
209         WRITE(numout,*) '            horizontal 2nd/4th order               nn_fct_h   = ', nn_fct_h
210         WRITE(numout,*) '            vertical   2nd/4th order               nn_fct_v   = ', nn_fct_v
211         WRITE(numout,*) '      MUSCL scheme                              ln_traadv_mus = ', ln_traadv_mus
212         WRITE(numout,*) '            + upstream scheme near river mouths    ln_mus_ups = ', ln_mus_ups
213         WRITE(numout,*) '      UBS scheme                                ln_traadv_ubs = ', ln_traadv_ubs
214         WRITE(numout,*) '            vertical   2nd/4th order               nn_ubs_v   = ', nn_ubs_v
215         WRITE(numout,*) '      QUICKEST scheme                           ln_traadv_qck = ', ln_traadv_qck
216      ENDIF
217      !
218      !                                !==  Parameter control & set nadv ==!
219      ioptio = 0                       
220      IF( ln_traadv_OFF ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_NO_adv   ;   ENDIF
221      IF( ln_traadv_cen ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_CEN      ;   ENDIF
222      IF( ln_traadv_fct ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_FCT      ;   ENDIF
223      IF( ln_traadv_mus ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_MUS      ;   ENDIF
224      IF( ln_traadv_ubs ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_UBS      ;   ENDIF
225      IF( ln_traadv_qck ) THEN   ;   ioptio = ioptio + 1   ;   nadv = np_QCK      ;   ENDIF
226      !
227      IF( ioptio /= 1 )   CALL ctl_stop( 'tra_adv_init: Choose ONE advection option in namelist namtra_adv' )
228      !
229      IF( ln_traadv_cen .AND. ( nn_cen_h /= 2 .AND. nn_cen_h /= 4 )   &          ! Centered
230                        .AND. ( nn_cen_v /= 2 .AND. nn_cen_v /= 4 )   ) THEN
231        CALL ctl_stop( 'tra_adv_init: CEN scheme, choose 2nd or 4th order' )
232      ENDIF
233      IF( ln_traadv_fct .AND. ( nn_fct_h /= 2 .AND. nn_fct_h /= 4 )   &          ! FCT
234                        .AND. ( nn_fct_v /= 2 .AND. nn_fct_v /= 4 )   ) THEN
235        CALL ctl_stop( 'tra_adv_init: FCT scheme, choose 2nd or 4th order' )
236      ENDIF
237      IF( ln_traadv_ubs .AND. ( nn_ubs_v /= 2 .AND. nn_ubs_v /= 4 )   ) THEN     ! UBS
238        CALL ctl_stop( 'tra_adv_init: UBS scheme, choose 2nd or 4th order' )
239      ENDIF
240      IF( ln_traadv_ubs .AND. nn_ubs_v == 4 ) THEN
241         CALL ctl_warn( 'tra_adv_init: UBS scheme, only 2nd FCT scheme available on the vertical. It will be used' )
242      ENDIF
243      IF( ln_isfcav ) THEN                                                       ! ice-shelf cavities
244         IF(  ln_traadv_cen .AND. nn_cen_v == 4    .OR.   &                            ! NO 4th order with ISF
245            & ln_traadv_fct .AND. nn_fct_v == 4   )   CALL ctl_stop( 'tra_adv_init: 4th order COMPACT scheme not allowed with ISF' )
246      ENDIF
247      !
248      !                                !==  Print the choice  ==! 
249      IF(lwp) THEN
250         WRITE(numout,*)
251         SELECT CASE ( nadv )
252         CASE( np_NO_adv  )   ;   WRITE(numout,*) '   ==>>>   NO T-S advection'
253         CASE( np_CEN     )   ;   WRITE(numout,*) '   ==>>>   CEN      scheme is used. Horizontal order: ', nn_cen_h,   &
254            &                                                                        ' Vertical   order: ', nn_cen_v
255         CASE( np_FCT     )   ;   WRITE(numout,*) '   ==>>>   FCT      scheme is used. Horizontal order: ', nn_fct_h,   &
256            &                                                                        ' Vertical   order: ', nn_fct_v
257         CASE( np_MUS     )   ;   WRITE(numout,*) '   ==>>>   MUSCL    scheme is used'
258         CASE( np_UBS     )   ;   WRITE(numout,*) '   ==>>>   UBS      scheme is used'
259         CASE( np_QCK     )   ;   WRITE(numout,*) '   ==>>>   QUICKEST scheme is used'
260         END SELECT
261      ENDIF
262      !
263      CALL tra_mle_init            !== initialisation of the Mixed Layer Eddy parametrisation (MLE)  ==!
264      !
265   END SUBROUTINE tra_adv_init
266
267  !!======================================================================
268END MODULE traadv
Note: See TracBrowser for help on using the repository browser.