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 branches/UKMO/dev_r5518_GO6_package_asm_surf_bgc/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: branches/UKMO/dev_r5518_GO6_package_asm_surf_bgc/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90 @ 8485

Last change on this file since 8485 was 7771, checked in by frrh, 7 years ago

Apply optimisations to various areas of code replacing the use of
allocated pointers with straightforward direct ALLOCATE and DEALLOCATE
operations.

These optimisations largely have an impact in models featuring MEDUSA,
i.e. those with significant numbers of tracers, although they are
expected to have a small impact in all configurations.

Code developed and tested in NEMO branch branches/UKMO/dev_r5518_optim_GO6_alloc
Tested in stand-alone GO6-GSI8, GO6-GSI8-MEDUSA and UKESM coupled models.
NEMO ticket #1821 documents this change further.

File size: 15.0 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   !!            4.0  !  2011-06  (G. Madec)  Addition of Mixed Layer Eddy parameterisation
9   !!----------------------------------------------------------------------
10
11   !!----------------------------------------------------------------------
12   !!   tra_adv      : compute ocean tracer advection trend
13   !!   tra_adv_ctl  : control the different options of advection scheme
14   !!----------------------------------------------------------------------
15   USE oce             ! ocean dynamics and active tracers
16   USE dom_oce         ! ocean space and time domain
17   USE domvvl          ! variable vertical scale factors
18   USE traadv_cen2     ! 2nd order centered scheme (tra_adv_cen2   routine)
19   USE traadv_tvd      ! TVD      scheme           (tra_adv_tvd    routine)
20   USE traadv_muscl    ! MUSCL    scheme           (tra_adv_muscl  routine)
21   USE traadv_muscl2   ! MUSCL2   scheme           (tra_adv_muscl2 routine)
22   USE traadv_ubs      ! UBS      scheme           (tra_adv_ubs    routine)
23   USE traadv_qck      ! QUICKEST scheme           (tra_adv_qck    routine)
24   USE traadv_eiv      ! eddy induced velocity     (tra_adv_eiv    routine)
25   USE traadv_mle      ! ML eddy induced velocity  (tra_adv_mle    routine)
26   USE cla             ! cross land advection      (cla_traadv     routine)
27   USE ldftra_oce      ! lateral diffusion coefficient on tracers
28   USE trd_oce         ! trends: ocean variables
29   USE trdtra          ! trends manager: tracers
30   !
31   USE in_out_manager  ! I/O manager
32   USE iom             ! I/O module
33   USE prtctl          ! Print control
34   USE lib_mpp         ! MPP library
35   USE wrk_nemo        ! Memory Allocation
36   USE timing          ! Timing
37   USE sbc_oce
38   USE diaptr          ! Poleward heat transport
39
40
41   IMPLICIT NONE
42   PRIVATE
43
44   PUBLIC   tra_adv        ! routine called by step module
45   PUBLIC   tra_adv_init   ! routine called by opa module
46
47   !                              !!* Namelist namtra_adv *
48   LOGICAL ::   ln_traadv_cen2     ! 2nd order centered scheme flag
49   LOGICAL ::   ln_traadv_tvd      ! TVD scheme flag
50   LOGICAL ::   ln_traadv_tvd_zts  ! TVD scheme flag with vertical sub time-stepping
51   LOGICAL ::   ln_traadv_muscl    ! MUSCL scheme flag
52   LOGICAL ::   ln_traadv_muscl2   ! MUSCL2 scheme flag
53   LOGICAL ::   ln_traadv_ubs      ! UBS scheme flag
54   LOGICAL ::   ln_traadv_qck      ! QUICKEST scheme flag
55   LOGICAL ::   ln_traadv_msc_ups  ! use upstream scheme within muscl
56
57
58   INTEGER ::   nadv   ! choice of the type of advection scheme
59
60   !! * Substitutions
61#  include "domzgr_substitute.h90"
62#  include "vectopt_loop_substitute.h90"
63   !!----------------------------------------------------------------------
64   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
65   !! $Id$
66   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
67   !!----------------------------------------------------------------------
68CONTAINS
69
70   SUBROUTINE tra_adv( kt )
71      !!----------------------------------------------------------------------
72      !!                  ***  ROUTINE tra_adv  ***
73      !!
74      !! ** Purpose :   compute the ocean tracer advection trend.
75      !!
76      !! ** Method  : - Update (ua,va) with the advection term following nadv
77      !!----------------------------------------------------------------------
78      !
79      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
80      !
81      INTEGER ::   jk   ! dummy loop index
82      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zun, zvn, zwn
83      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds   ! 3D workspace
84      !!----------------------------------------------------------------------
85      !
86      IF( nn_timing == 1 )  CALL timing_start('tra_adv')
87      !
88      ALLOCATE(zun(1:jpi, 1:jpj, 1:jpk))
89      ALLOCATE(zvn(1:jpi, 1:jpj, 1:jpk))
90      ALLOCATE(zwn(1:jpi, 1:jpj, 1:jpk))
91      !                                          ! set time step
92      IF( neuler == 0 .AND. kt == nit000 ) THEN     ! at nit000
93         r2dtra(:) =  rdttra(:)                          ! = rdtra (restarting with Euler time stepping)
94      ELSEIF( kt <= nit000 + 1) THEN                ! at nit000 or nit000+1
95         r2dtra(:) = 2._wp * rdttra(:)                   ! = 2 rdttra (leapfrog)
96      ENDIF
97      !
98      IF( nn_cla == 1 .AND. cp_cfg == 'orca' .AND. jp_cfg == 2 )   CALL cla_traadv( kt )       !==  Cross Land Advection  ==! (hor. advection)
99      !
100      !                                               !==  effective transport  ==!
101      DO jk = 1, jpkm1
102         zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)                  ! eulerian transport only
103         zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)
104         zwn(:,:,jk) = e1t(:,:) * e2t(:,:)      * wn(:,:,jk)
105      END DO
106      !
107      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN
108         zun(:,:,:) = zun(:,:,:) + un_td(:,:,:)
109         zvn(:,:,:) = zvn(:,:,:) + vn_td(:,:,:)
110      ENDIF
111      !
112      zun(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom
113      zvn(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom
114      zwn(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom
115      !
116      IF( lk_traldf_eiv .AND. .NOT. ln_traldf_grif )   &
117         &              CALL tra_adv_eiv( kt, nit000, zun, zvn, zwn, 'TRA' )    ! add the eiv transport (if necessary)
118      !
119      IF( ln_mle    )   CALL tra_adv_mle( kt, nit000, zun, zvn, zwn, 'TRA' )    ! add the mle transport (if necessary)
120      !
121      CALL iom_put( "uocetr_eff", zun )                                         ! output effective transport     
122      CALL iom_put( "vocetr_eff", zvn )
123      CALL iom_put( "wocetr_eff", zwn )
124      !
125      IF( ln_diaptr )   CALL dia_ptr( zvn )                                     ! diagnose the effective MSF
126      !
127      IF( l_trdtra )   THEN                    !* Save ta and sa trends
128         ALLOCATE(ztrdt( 1:jpi, 1:jpj, 1:jpk) )
129         ALLOCATE(ztrds( 1:jpi, 1:jpj, 1:jpk) )
130         ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
131         ztrds(:,:,:) = tsa(:,:,:,jp_sal)
132      ENDIF
133      !
134      SELECT CASE ( nadv )                            !==  compute advection trend and add it to general trend  ==!
135      CASE ( 1 )   ;    CALL tra_adv_cen2   ( kt, nit000, 'TRA',         zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  2nd order centered
136      CASE ( 2 )   ;    CALL tra_adv_tvd    ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  TVD
137      CASE ( 3 )   ;    CALL tra_adv_muscl  ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb,      tsa, jpts, ln_traadv_msc_ups )   !  MUSCL
138      CASE ( 4 )   ;    CALL tra_adv_muscl2 ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  MUSCL2
139      CASE ( 5 )   ;    CALL tra_adv_ubs    ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  UBS
140      CASE ( 6 )   ;    CALL tra_adv_qck    ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  QUICKEST
141      CASE ( 7 )   ;    CALL tra_adv_tvd_zts( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )   !  TVD ZTS
142      !
143      CASE (-1 )                                      !==  esopa: test all possibility with control print  ==!
144         CALL tra_adv_cen2  ( kt, nit000, 'TRA',         zun, zvn, zwn, tsb, tsn, tsa, jpts )         
145         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv0 - Ta: ', mask1=tmask,               &
146            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
147         CALL tra_adv_tvd   ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )         
148         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv1 - Ta: ', mask1=tmask,               &
149            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
150         CALL tra_adv_muscl ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb,      tsa, jpts, ln_traadv_msc_ups )         
151         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv3 - Ta: ', mask1=tmask,               &
152            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
153         CALL tra_adv_muscl2( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )         
154         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv4 - Ta: ', mask1=tmask,               &
155            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
156         CALL tra_adv_ubs   ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )         
157         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv5 - Ta: ', mask1=tmask,               &
158            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
159         CALL tra_adv_qck   ( kt, nit000, 'TRA', r2dtra, zun, zvn, zwn, tsb, tsn, tsa, jpts )         
160         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv6 - Ta: ', mask1=tmask,               &
161            &          tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
162      END SELECT
163      !
164      IF( l_trdtra )   THEN                      ! save the advective trends for further diagnostics
165         DO jk = 1, jpkm1
166            ztrdt(:,:,jk) = tsa(:,:,jk,jp_tem) - ztrdt(:,:,jk)
167            ztrds(:,:,jk) = tsa(:,:,jk,jp_sal) - ztrds(:,:,jk)
168         END DO
169         CALL trd_tra( kt, 'TRA', jp_tem, jptra_totad, ztrdt )
170         CALL trd_tra( kt, 'TRA', jp_sal, jptra_totad, ztrds )
171         DEALLOCATE (ztrdt)
172         DEALLOCATE (ztrds)
173      ENDIF
174      !                                              ! print mean trends (used for debugging)
175      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' adv  - Ta: ', mask1=tmask,               &
176         &                       tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
177      !
178      IF( nn_timing == 1 )  CALL timing_stop( 'tra_adv' )
179      !
180      DEALLOCATE(zun)
181      DEALLOCATE(zvn)
182      DEALLOCATE(zwn)
183      !                                         
184   END SUBROUTINE tra_adv
185
186
187   SUBROUTINE tra_adv_init
188      !!---------------------------------------------------------------------
189      !!                  ***  ROUTINE tra_adv_init  ***
190      !!               
191      !! ** Purpose :   Control the consistency between namelist options for
192      !!              tracer advection schemes and set nadv
193      !!----------------------------------------------------------------------
194      INTEGER ::   ioptio
195      INTEGER ::   ios                 ! Local integer output status for namelist read
196      !!
197      NAMELIST/namtra_adv/ ln_traadv_cen2 , ln_traadv_tvd,     &
198         &                 ln_traadv_muscl, ln_traadv_muscl2,  &
199         &                 ln_traadv_ubs  , ln_traadv_qck,     &
200         &                 ln_traadv_msc_ups, ln_traadv_tvd_zts
201      !!----------------------------------------------------------------------
202
203      REWIND( numnam_ref )              ! Namelist namtra_adv in reference namelist : Tracer advection scheme
204      READ  ( numnam_ref, namtra_adv, IOSTAT = ios, ERR = 901)
205901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtra_adv in reference namelist', lwp )
206
207      REWIND( numnam_cfg )              ! Namelist namtra_adv in configuration namelist : Tracer advection scheme
208      READ  ( numnam_cfg, namtra_adv, IOSTAT = ios, ERR = 902 )
209902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtra_adv in configuration namelist', lwp )
210      IF(lwm) WRITE ( numond, namtra_adv )
211
212      IF(lwp) THEN                    ! Namelist print
213         WRITE(numout,*)
214         WRITE(numout,*) 'tra_adv_init : choice/control of the tracer advection scheme'
215         WRITE(numout,*) '~~~~~~~~~~~'
216         WRITE(numout,*) '   Namelist namtra_adv : chose a advection scheme for tracers'
217         WRITE(numout,*) '      2nd order advection scheme     ln_traadv_cen2    = ', ln_traadv_cen2
218         WRITE(numout,*) '      TVD advection scheme           ln_traadv_tvd     = ', ln_traadv_tvd
219         WRITE(numout,*) '      MUSCL  advection scheme        ln_traadv_muscl   = ', ln_traadv_muscl
220         WRITE(numout,*) '      MUSCL2 advection scheme        ln_traadv_muscl2  = ', ln_traadv_muscl2
221         WRITE(numout,*) '      UBS    advection scheme        ln_traadv_ubs     = ', ln_traadv_ubs
222         WRITE(numout,*) '      QUICKEST advection scheme      ln_traadv_qck     = ', ln_traadv_qck
223         WRITE(numout,*) '      upstream scheme within muscl   ln_traadv_msc_ups = ', ln_traadv_msc_ups
224         WRITE(numout,*) '      TVD advection scheme with zts  ln_traadv_tvd_zts = ', ln_traadv_tvd_zts
225      ENDIF
226
227      ioptio = 0                      ! Parameter control
228      IF( ln_traadv_cen2   )   ioptio = ioptio + 1
229      IF( ln_traadv_tvd    )   ioptio = ioptio + 1
230      IF( ln_traadv_muscl  )   ioptio = ioptio + 1
231      IF( ln_traadv_muscl2 )   ioptio = ioptio + 1
232      IF( ln_traadv_ubs    )   ioptio = ioptio + 1
233      IF( ln_traadv_qck    )   ioptio = ioptio + 1
234      IF( ln_traadv_tvd_zts)   ioptio = ioptio + 1
235      IF( lk_esopa         )   ioptio =          1
236
237      IF( ( ln_traadv_muscl .OR. ln_traadv_muscl2 .OR. ln_traadv_ubs .OR. ln_traadv_qck .OR. ln_traadv_tvd_zts )   &
238         .AND. ln_isfcav )   CALL ctl_stop( 'Only traadv_cen2 and traadv_tvd is compatible with ice shelf cavity')
239
240      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE advection scheme in namelist namtra_adv' )
241
242      !                              ! Set nadv
243      IF( ln_traadv_cen2   )   nadv =  1
244      IF( ln_traadv_tvd    )   nadv =  2
245      IF( ln_traadv_muscl  )   nadv =  3
246      IF( ln_traadv_muscl2 )   nadv =  4
247      IF( ln_traadv_ubs    )   nadv =  5
248      IF( ln_traadv_qck    )   nadv =  6
249      IF( ln_traadv_tvd_zts)   nadv =  7
250      IF( lk_esopa         )   nadv = -1
251
252      IF(lwp) THEN                   ! Print the choice
253         WRITE(numout,*)
254         IF( nadv ==  1 )   WRITE(numout,*) '         2nd order scheme is used'
255         IF( nadv ==  2 )   WRITE(numout,*) '         TVD       scheme is used'
256         IF( nadv ==  3 )   WRITE(numout,*) '         MUSCL     scheme is used'
257         IF( nadv ==  4 )   WRITE(numout,*) '         MUSCL2    scheme is used'
258         IF( nadv ==  5 )   WRITE(numout,*) '         UBS       scheme is used'
259         IF( nadv ==  6 )   WRITE(numout,*) '         QUICKEST  scheme is used'
260         IF( nadv ==  7 )   WRITE(numout,*) '         TVD ZTS   scheme is used'
261         IF( nadv == -1 )   WRITE(numout,*) '         esopa test: use all advection scheme'
262      ENDIF
263      !
264      CALL tra_adv_mle_init          ! initialisation of the Mixed Layer Eddy parametrisation (MLE)
265      !
266   END SUBROUTINE tra_adv_init
267
268  !!======================================================================
269END MODULE traadv
Note: See TracBrowser for help on using the repository browser.