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.
trazdf_tam.F90 in branches/TAM_V3_0/NEMOTAM/OPATAM_SRC/TRA – NEMO

source: branches/TAM_V3_0/NEMOTAM/OPATAM_SRC/TRA/trazdf_tam.F90 @ 1885

Last change on this file since 1885 was 1885, checked in by rblod, 14 years ago

add TAM sources
File size: 11.8 KB
Line 
1MODULE trazdf_tam
2#ifdef key_tam
3   !!==============================================================================
4   !!                 ***  MODULE  trazdf_zdf  ***
5   !! Ocean active tracers:  vertical component of the tracer mixing trend
6   !!                 Tangent and Adjoint Module
7   !!==============================================================================
8   !! History of the direct module: 
9   !!            9.0  !  05-11  (G. Madec)  Original code
10   !! History of the TAM module:
11   !!            9.0  !  08-06  (A. Vidard) Skeleton
12   !!            9.0  !  09-01  (A. Vidard) TAM of the 05-11 version
13   !!----------------------------------------------------------------------
14
15   !!----------------------------------------------------------------------
16   !!   tra_zdf_tan  : Update the tracer trend with the vertical diffusion (tan)
17   !!   tra_zdf_adj  : Update the tracer trend with the vertical diffusion (adj)
18   !!       zdf_ctl  : ???
19   !!----------------------------------------------------------------------
20   USE par_kind      , ONLY: & ! Precision variables
21      & wp
22   USE par_oce       , ONLY: & ! Ocean space and time domain variables
23      & jpk,                 &
24      & lk_esopa
25   USE oce_tam       , ONLY: & ! ocean dynamics and active tracers
26      & ta_tl,               &
27      & sa_tl,               &
28      & ta_ad,               &
29      & sa_ad
30   USE dom_oce       , ONLY: & ! ocean space and time domain
31      & tmask,               &
32      & lk_vvl,              &
33      & neuler,              &
34      & rdttra,              &
35      & ln_sco
36   USE ldftra_oce    , ONLY: &
37      & ln_traldf_iso,       &
38      & ln_traldf_hor
39   USE zdf_oce      , ONLY : &
40      & ln_zdfexp
41   USE trazdf_exp_tam, ONLY: & ! vertical diffusion: explicit (tra_zdf_exp     routine)
42      & tra_zdf_exp_tan,     &
43      & tra_zdf_exp_adj,     &
44      & tra_zdf_exp_adj_tst
45   USE trazdf_imp_tam, ONLY: & ! vertical diffusion: implicit (tra_zdf_imp     routine)
46      & tra_zdf_imp_tan,     &
47      & tra_zdf_imp_adj,     &
48      & tra_zdf_imp_adj_tst, &
49      & tra_zdf_imp_tlm_tst
50   USE in_out_manager, ONLY: & ! I/O manager
51      & lwp,                 &
52      & numout,              & 
53      & nitend,              & 
54      & nit000,              &
55      & ctl_stop
56   USE prtctl        , ONLY: &
57      & prt_ctl
58
59   IMPLICIT NONE
60   PRIVATE
61
62   PUBLIC  &
63      & tra_zdf_tan, &
64      & tra_zdf_adj         ! routines called by step_tam.F90
65   PUBLIC  &
66      & tra_zdf_adj_tst, &  ! routine called by tst.F90
67      & tra_zdf_tlm_tst     ! routine called by tst.F90
68   INTEGER ::   nzdf = 0               ! type vertical diffusion algorithm used
69      !                                ! defined from ln_zdf...  namlist logicals)
70
71   REAL(wp), DIMENSION(jpk) ::   r2dt  ! vertical profile time-step, = 2 rdttra
72      !                                ! except at nit000 (=rdttra) if neuler=0
73
74   !! * Substitutions
75#  include "domzgr_substitute.h90"
76#  include "zdfddm_substitute.h90"
77#  include "vectopt_loop_substitute.h90"
78
79CONTAINS
80   
81   SUBROUTINE tra_zdf_tan( kt )
82      !!----------------------------------------------------------------------
83      !!                  ***  ROUTINE tra_zdf_tan  ***
84      !!
85      !! ** Purpose of the direct routine:   
86      !!            compute the vertical ocean tracer physics.
87      !!---------------------------------------------------------------------
88      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
89
90      !!---------------------------------------------------------------------
91
92      IF( kt == nit000 )   CALL zdf_ctl_tam          ! initialisation & control of options
93
94      !                                          ! set time step
95      IF( neuler == 0 .AND. kt == nit000 ) THEN     ! at nit000
96         r2dt(:) =  rdttra(:)                          ! = rdtra (restarting with Euler time stepping)
97      ELSEIF( kt <= nit000 + 1) THEN                ! at nit000 or nit000+1
98         r2dt(:) = 2. * rdttra(:)                      ! = 2 rdttra (leapfrog)
99      ENDIF
100
101
102      IF( lk_vvl ) THEN !CALL dom_vvl_ssh( kt )      ! compute ssha field
103         IF (lwp) WRITE(numout,*) 'key_vvl not available in tangent yet'
104         CALL abort
105      END IF
106         
107      SELECT CASE ( nzdf )                       ! compute lateral mixing trend and add it to the general trend
108      CASE ( -1 )                                       ! esopa: test all possibility with control print
109         CALL tra_zdf_exp_tan    ( kt, r2dt )
110         CALL prt_ctl( tab3d_1=ta_tl, clinfo1=' zdf0 - Ta_tl: ', mask1=tmask,               &
111            &          tab3d_2=sa_tl, clinfo2=       ' Sa_tl: ', mask2=tmask, clinfo3='tra' )
112         CALL tra_zdf_imp_tan    ( kt, r2dt )
113         CALL prt_ctl( tab3d_1=ta_tl, clinfo1=' zdf1 - Ta_tl: ', mask1=tmask,               &
114            &          tab3d_2=sa_tl, clinfo2=       ' Sa_tl: ', mask2=tmask, clinfo3='tra' )
115
116      CASE ( 0 )                                       ! explicit scheme
117         CALL tra_zdf_exp_tan    ( kt, r2dt )
118
119      CASE ( 1 )                                       ! implicit scheme
120         CALL tra_zdf_imp_tan    ( kt, r2dt )
121
122      END SELECT
123
124   END SUBROUTINE tra_zdf_tan
125   SUBROUTINE tra_zdf_adj( kt )
126      !!----------------------------------------------------------------------
127      !!                  ***  ROUTINE tra_zdf_adj  ***
128      !!
129      !! ** Purpose of the direct routine:   
130      !!            compute the vertical ocean tracer physics.
131      !!---------------------------------------------------------------------
132      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
133
134      !!---------------------------------------------------------------------
135
136      IF( kt == nitend )   CALL zdf_ctl_tam          ! initialisation & control of options
137      !                                          ! set time step
138      IF( neuler == 0 .AND. kt == nit000 ) THEN     ! at nit000
139         r2dt(:) =  rdttra(:)                          ! = rdtra (restarting with Euler time stepping)
140      ELSEIF( kt <= nit000 + 1) THEN                ! at nit000 or nit000+1
141         r2dt(:) = 2. * rdttra(:)                      ! = 2 rdttra (leapfrog)
142      ELSEIF( kt == nitend) THEN
143         r2dt(:) = 2. * rdttra(:)
144      ENDIF
145         
146      SELECT CASE ( nzdf )                       ! compute lateral mixing trend and add it to the general trend
147      CASE ( -1 )                                       ! esopa: test all possibility with control print
148         CALL tra_zdf_exp_adj    ( kt, r2dt )
149         CALL prt_ctl( tab3d_1=ta_ad, clinfo1=' zdf0 - Ta_ad: ', mask1=tmask,               &
150            &          tab3d_2=sa_ad, clinfo2=       ' Sa_ad: ', mask2=tmask, clinfo3='tra' )
151         CALL tra_zdf_imp_adj    ( kt, r2dt )
152         CALL prt_ctl( tab3d_1=ta_ad, clinfo1=' zdf1 - Ta_ad: ', mask1=tmask,               &
153            &          tab3d_2=sa_ad, clinfo2=       ' Sa_ad: ', mask2=tmask, clinfo3='tra' )
154
155      CASE ( 0 )                                       ! explicit scheme
156         CALL tra_zdf_exp_adj    ( kt, r2dt )
157
158      CASE ( 1 )                                       ! implicit scheme
159         CALL tra_zdf_imp_adj    ( kt, r2dt )
160
161      END SELECT
162
163
164      IF( lk_vvl ) THEN !CALL dom_vvl_ssh( kt )      ! compute ssha field
165         IF (lwp) WRITE(numout,*) 'key_vvl not available in tangent yet'
166         CALL abort
167      END IF
168
169   END SUBROUTINE tra_zdf_adj
170   SUBROUTINE tra_zdf_adj_tst( kumadt )
171      !! ** Purpose : Test the adjoint routines.
172      !!
173      !! ** Method  : Verify the scalar product
174      !!           
175      !!                 ( L dx )^T W dy  =  dx^T L^T W dy
176      !!
177      !!              where  L   = tangent routine
178      !!                     L^T = adjoint routine
179      !!                     W   = diagonal matrix of scale factors
180      !!                     dx  = input perturbation (random field)
181      !!                     dy  = L dx
182      !!
183      !!       
184      !!-----------------------------------------------------------------------
185      !! * Modules used
186
187      !! * Arguments
188      INTEGER, INTENT(IN) :: &
189         & kumadt             ! Output unit
190 
191      !! * Local declarations
192      ! init
193      CALL zdf_ctl_tam 
194      ! Test the explicit formulation           
195      CALL tra_zdf_exp_adj_tst    ( kumadt )
196      ! Test the implicit formulation           
197      CALL tra_zdf_imp_adj_tst    ( kumadt )
198   END SUBROUTINE tra_zdf_adj_tst
199   !!==============================================================================
200   SUBROUTINE zdf_ctl_tam
201      !!----------------------------------------------------------------------
202      !!                 ***  ROUTINE zdf_ctl_tam  ***
203      !!
204      !! ** Purpose :   Choose the vertical mixing scheme
205      !!
206      !! ** Method  :   Set nzdf from ln_zdfexp
207      !!      nzdf = 0   explicit (time-splitting) scheme (ln_zdfexp=T)
208      !!           = 1   implicit (euler backward) scheme (ln_zdfexp=F)
209      !!      NB: rotation of lateral mixing operator or TKE or KPP scheme,
210      !!      the implicit scheme is required.
211      !!----------------------------------------------------------------------
212      USE zdftke
213      USE zdfkpp
214      !!----------------------------------------------------------------------
215
216      !  Define the vertical tracer physics scheme
217      ! ==========================================
218
219      ! Choice from ln_zdfexp already read in namelist in zdfini module
220      IF( ln_zdfexp ) THEN               ! use explicit scheme
221         nzdf = 0
222      ELSE                               ! use implicit scheme
223         nzdf = 1
224      ENDIF
225
226      ! Force implicit schemes
227      IF( lk_zdftke .OR. lk_zdfkpp   )   nzdf = 1      ! TKE or KPP physics
228      IF( ln_traldf_iso              )   nzdf = 1      ! iso-neutral lateral physics
229      IF( ln_traldf_hor .AND. ln_sco )   nzdf = 1      ! horizontal lateral physics in s-coordinate
230
231      IF( ln_zdfexp .AND. nzdf == 1 )   THEN
232         CALL ctl_stop( 'tra_zdf_tam : If using the rotation of lateral mixing operator or TKE ', &
233            &           '            or KPP scheme, the implicit scheme is required, set ln_zdfexp = .false.' )
234      ENDIF
235
236      ! Test: esopa
237      IF( lk_esopa )    nzdf = -1                      ! All schemes used
238
239      IF(lwp) THEN
240         WRITE(numout,*)
241         WRITE(numout,*) 'tra:zdf_ctl_tam : vertical tracer physics scheme'
242         WRITE(numout,*) '~~~~~~~~~~~~~~~'
243         IF( nzdf == -1 )   WRITE(numout,*) '              ESOPA test All scheme used'
244         IF( nzdf ==  0 )   WRITE(numout,*) '              Explicit time-splitting scheme'
245         IF( nzdf ==  1 )   WRITE(numout,*) '              Implicit (euler backward) scheme'
246      ENDIF
247
248   END SUBROUTINE zdf_ctl_tam
249
250   SUBROUTINE tra_zdf_tlm_tst( kumadt )
251      !!-----------------------------------------------------------------------
252      !!
253      !!                  ***  ROUTINE tra_zdf_tlm_tst ***
254      !!
255      !! ** Purpose : Test the tangent routine.
256      !!
257      !! ** Method  : Verify the tangent with Taylor expansion
258      !!           
259      !!                 M(x+hdx) = M(x) + L(hdx) + O(h^2)
260      !!
261      !!              where  L   = tangent routine
262      !!                     M   = direct routine
263      !!                     dx  = input perturbation (random field)
264      !!                     h   = ration on perturbation
265      !!                   
266      !! History :
267      !!        ! 09-08 (A. Vigilant)
268      !!-----------------------------------------------------------------------
269      !! * Modules used
270
271      !! * Arguments
272      INTEGER, INTENT(IN) :: &
273         & kumadt             ! Output unit
274 
275      !! * Local declarations
276      ! init
277      CALL zdf_ctl_tam 
278      ! Test the explicit formulation           
279!      CALL tra_zdf_exp_adj_tst    ( kumadt )
280      ! Test the implicit formulation           
281      CALL tra_zdf_imp_tlm_tst    ( kumadt )
282   END SUBROUTINE tra_zdf_tlm_tst
283#endif
284END MODULE trazdf_tam
Note: See TracBrowser for help on using the repository browser.