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.
step_tam.F90 in branches/2012/dev_r3604_LEGI8_TAM/NEMOGCM/NEMO/OPATAM_SRC – NEMO

source: branches/2012/dev_r3604_LEGI8_TAM/NEMOGCM/NEMO/OPATAM_SRC/step_tam.F90 @ 3641

Last change on this file since 3641 was 3641, checked in by pabouttier, 11 years ago

Add possibility to run in stand-alone Tangent Linear Model - See ticket 1014

  • Property svn:executable set to *
File size: 41.9 KB
Line 
1MODULE step_tam
2#ifdef key_tam
3   !!======================================================================
4   !!                       ***  MODULE step_tam  ***
5   !! Time-stepping    : manager of the adjoint ocean time stepping
6   !!                    Tangent and Adjoint module
7   !!======================================================================
8   !! History of
9   !! the direct:      !  91-03  ()  Original code
10   !!                  !  91-11  (G. Madec)
11   !!                  !  92-06  (M. Imbard)  add a first output record
12   !!                  !  96-04  (G. Madec)  introduction of dynspg
13   !!                  !  96-04  (M.A. Foujols)  introduction of passive tracer
14   !!             8.0  !  97-06  (G. Madec)  new architecture of call
15   !!             8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface
16   !!             8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit
17   !!             8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions
18   !!             9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking
19   !!             " "  !  04-08  (C. Talandier) New trends organization
20   !!             " "  !  05-01  (C. Ethe) Add the KPP closure scheme
21   !!             " "  !  05-11  (V. Garnier) Surface pressure gradient organization
22   !!             " "  !  05-11  (G. Madec)  Reorganisation of tra and dyn calls
23   !!             " "  !  06-01  (L. Debreu, C. Mazauric)  Agrif implementation
24   !!             " "  !  06-07  (S. Masson)  restart using iom
25   !!    "        " "  !  07-04  (K. Mogensen, A. Weaver, M. Martin) Assimilation interface
26   !! History of the TAM
27   !!                  !  08-06  (A. Vidard and A. Weaver) Tangent and Adjoint version of 9.0
28   !!                  !  08-11  (A. Vidard) Nemo v3 update
29   !!                  !  06-09  (F. Vigilant)  Modified to split NEMOVAR / NEMOTAM
30   !!                  !  07-12  (P.-A. Bouttier) Phasing with 3.4 version
31   !!----------------------------------------------------------------------
32   !!   stp_tam        : OPA system time-stepping (tangent linear)
33   !!   stp_adj        : OPA system time-stepping (adjoint)
34   !!----------------------------------------------------------------------
35
36   USE step_oce_tam
37#if defined key_agrif
38#error 'agrif not yet implemented in nemotam'
39#endif
40
41   IMPLICIT NONE
42   PRIVATE
43
44   PUBLIC stp_tan,      &
45      &   stp_adj,      & ! called by simvar.F90
46      &   stp_adj_tst
47
48   !! * Substitutions
49#  include "domzgr_substitute.h90"
50#  include "zdfddm_substitute.h90"
51
52CONTAINS
53   SUBROUTINE stp_tan( kstp )
54      !!----------------------------------------------------------------------
55      !!                     ***  ROUTINE stp_tan  ***
56      !!
57      !! ** Purpose of the direct routine:
58      !!              - Time stepping of OPA (momentum and active tracer eqs.)
59      !!              - Time stepping of LIM (dynamic and thermodynamic eqs.)
60      !!              - Tme stepping  of TRC (passive tracer eqs.)
61      !!
62      !! ** Method of the direct routine:
63      !!              -1- Update forcings and data
64      !!              -2- Update ocean physics
65      !!              -3- Compute the t and s trends
66      !!              -4- Update t and s
67      !!              -5- Compute the momentum trends
68      !!              -6- Update the horizontal velocity
69      !!              -7- Compute the diagnostics variables (rd,N2, div,cur,w)
70      !!              -8- Outputs and diagnostics
71      !!----------------------------------------------------------------------
72      !! * Arguments
73      INTEGER, INTENT( in ) :: kstp   ! ocean time-step index
74
75      !! * local declarations
76      INTEGER ::   indic    ! error indicator if < 0
77      !! ---------------------------------------------------------------------
78
79      indic = 0                    ! reset to no error condition
80      IF ( kstp == nit000 )    CALL tl_trj_wri(nit000-1)
81      IF ( kstp /= nit000 )    CALL day_tam( kstp, 0 )             ! Calendar (day was already called at nit000 in day_init)
82
83                               CALL iom_setkt( kstp )                          ! say to iom that we are at time step kstp
84      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
85      ! Update data, open boundaries, surface boundary condition (including sea-ice)
86      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
87
88      ! update 3D temperature data ... not needed in tangent
89
90      ! update 3D salinity data ... not needed in tangent (to be investigated, see sbc_ssr)
91                             CALL sbc_tan ( kstp ) ! Sea boundary condition (including sea-ice)
92
93      ! Output the initial state and forcings ... not needed in tangent
94
95      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
96      !  Ocean dynamics : ssh, wn, hdiv, rot                                 !
97      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
98
99                             CALL ssh_wzv_tan( kstp )         ! after ssh & vertical velocity
100
101      ! saving direct variables ua,va, ta, sa before entering in tracer
102
103      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
104      ! Ocean physics update                (ua, va, ta, sa used as workspace)
105      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
106                             CALL bn2_tan( tsb, tsb_tl, rn2b_tl )        ! now Brunt-Vaisala frequency
107                             CALL bn2_tan( tsn, tsn_tl, rn2_tl )         ! now Brunt-Vaisala frequency
108
109      !-----------------------------------------------------------------------
110      !  VERTICAL PHYSICS
111      !-----------------------------------------------------------------------
112                             CALL zdf_bfr_tan ( kstp ) ! bottom friction...
113
114      !                                                ! Vertical eddy viscosity and diffusivity coefficients
115      ! Richardson number dependent Kz  ... not available
116      ! TKE closure scheme for Kz ... not available
117      ! KPP closure scheme for Kz ... not available
118      ! Constant Kz (reset avt, avm[uv] to the background value)...
119      ! increase diffusivity at rivers mouths... not needed in tangent
120      ! enhanced vertical eddy diffusivity ... not needed in tangent with lk_zdfcst_tan
121      ! double diffusive mixing ... not needed in tangent with lk_zdfcst_tan
122      ! mixed layer depth... not needed in tangent with lk_zdfcst_tan
123      !-----------------------------------------------------------------------
124      !  LATERAL PHYSICS
125      !-----------------------------------------------------------------------
126      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
127      !-----------------------------------------------------------------------
128      ! before slope of the lateral mixing... not needed in tangent with lk_zdfcst_tan
129#if defined key_traldf_c2d
130      ! eddy induced velocity coefficient... not needed in tangent with lk_zdfcst_tan
131#endif
132      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
133      ! diagnostics and outputs             (ua, va, ta, sa used as workspace)
134      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
135      !      not needed in tangent
136#if defined key_top
137      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
138      ! Passive Tracer Model
139      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
140      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
141      !-----------------------------------------------------------------------
142
143      ! time-stepping... not needed in tangent for the time being
144
145#endif
146
147      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
148      ! Active tracers                              (ua, va used as workspace)
149      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
150                             tsa_tl(:,:,:,:) = 0.0_wp            ! set tracer trends to zero
151
152      ! apply tracer assimilation increment   ... not needed in tangent
153                             CALL tra_sbc_tan( kstp )       ! surface boundary condition
154      IF( ln_traqsr      )   CALL tra_qsr_tan( kstp )       ! penetrative solar radiation qsr
155      IF( ln_trabbc      )   CALL tra_bbc_tan( kstp )       ! bottom heat flux
156      IF( lk_trabbl  )       CALL tra_bbl_tan( kstp )   ! diffusive bottom boundary layer scheme
157      !! advective (and/or diffusive) bottom boundary layer scheme ... currently not available
158      IF( ln_tradmp      )   CALL tra_dmp_tan( kstp )       ! internal damping trends
159                             CALL tra_adv_tan( kstp )       ! horizontal & vertical advection
160                             CALL tra_ldf_tan( kstp )       ! lateral mixing
161                             CALL tra_zdf_tan( kstp )       ! vertical mixing
162
163      IF( ln_dynhpg_imp  ) THEN                             ! semi-implicit hpg (time stepping then eos)
164      ! update the new (t,s) fields by non
165      ! penetrative convective adjustment ... not available
166                             CALL tra_nxt_tan( kstp )                              ! tracer fields at next time step
167                             CALL eos_tan( tsa, tsa_tl,          &                 ! Time-filtered in situ density used in dynhpg module
168                                         & rhd_tl, rhop_tl         )
169         IF( ln_zps    )     CALL zps_hde_tan( kstp, jpts, tsa, tsa_tl, rhd_tl,     &
170            &                                  gtsu_tl, gru_tl, gtsv_tl,    &   ! Partial steps: time filtered hor. gradient
171            &                                  grv_tl )                 ! of t, s, rd at the bottom ocean level
172      ELSE                                              ! centered hpg (default case)
173                             CALL eos_tan( tsn, tsn_tl, &                 ! now (swap=before) in situ density for dynhpg module
174            &                              rhd_tl, rhop_tl )
175
176         IF( ln_zps    )     CALL zps_hde_tan( kstp, jpts, tsn, tsn_tl, rhd_tl,     &
177            &                                  gtsu_tl, gru_tl, gtsv_tl,    &   ! Partial steps: time filtered hor. gradient
178            &                                  grv_tl )
179                             CALL tra_nxt_tan( kstp )                              ! tracer fields at next time step
180      ENDIF
181
182      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
183      ! Dynamics                                    (ta, sa used as workspace)
184      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
185                             ua_tl(:,:,:) = 0.0_wp               ! set dynamics trends to zero
186                             va_tl(:,:,:) = 0.0_wp
187
188      ! apply dynamics assimilation increment ... not needed in tangent
189                             CALL dyn_adv_tan( kstp )       ! advection (vector or flux form)
190                             CALL dyn_vor_tan( kstp )       ! vorticity term including Coriolis
191                             CALL dyn_ldf_tan( kstp )       ! lateral mixing
192                             CALL dyn_hpg_tan( kstp )           ! horizontal gradient of Hydrostatic pressure
193                             CALL dyn_bfr_tan( kstp )           ! bottom friction
194                             CALL dyn_zdf_tan( kstp )           ! vertical diffusion
195                             CALL dyn_spg_tan( kstp, indic )    ! surface pressure gradient
196                             CALL dyn_nxt_tan( kstp )           ! lateral velocity at next time step
197                             CALL ssh_nxt_tan( kstp )           ! sea surface height at next time step
198
199      ! vertical mesh at next time step ... not available
200      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
201      ! Computation of diagnostic variables
202      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
203      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
204      !-----------------------------------------------------------------------
205      CALL tl_trj_wri( kstp )
206      CALL trj_rea( kstp, 1) ! ... Read basic state trajectory at end of current step
207
208      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
209      ! Control, and restarts
210      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
211      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
212      !-----------------------------------------------------------------------
213      !                                                     ! Time loop: control and print
214!*B This fails with cgmod. To be revised      CALL stp_ctl_tan( kstp, indic, 0 )
215      IF( indic < 0          )   CALL ctl_stop( 'step_tan: indic < 0' )
216      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
217      ! diagnostics and outputs
218      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
219      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
220      !-----------------------------------------------------------------------
221      IF ( nstop == 0 ) THEN                                ! Diagnostics
222         ! not needed in tangent
223      ENDIF
224   END SUBROUTINE stp_tan
225
226   SUBROUTINE stp_adj( kstp )
227      !!----------------------------------------------------------------------
228      !!                     ***  ROUTINE stp_adj  ***
229      !!
230      !! ** Purpose of the direct routine:
231      !!              - Time stepping of OPA (momentum and active tracer eqs.)
232      !!              - Time stepping of LIM (dynamic and thermodynamic eqs.)
233      !!              - Tme stepping  of TRC (passive tracer eqs.)
234      !!
235      !! ** Method of the direct routine:
236      !!              -1- Update forcings and data
237      !!              -2- Update ocean physics
238      !!              -3- Compute the t and s trends
239      !!              -4- Update t and s
240      !!              -5- Compute the momentum trends
241      !!              -6- Update the horizontal velocity
242      !!              -7- Compute the diagnostics variables (rd,N2, div,cur,w)
243      !!              -8- Outputs and diagnostics
244      !!----------------------------------------------------------------------
245      !! * Arguments
246      INTEGER, INTENT( in ) :: kstp   ! ocean time-step index
247      !! * local declarations
248      INTEGER ::   indic    ! error indicator if < 0
249      !! ---------------------------------------------------------------------
250
251      indic = 1                    ! reset to no error condition
252                             CALL day_tam( kstp, 1 )            ! Calendar
253
254      ! ... Read basic state trajectory at end of previous step
255                             CALL trj_rea( ( kstp - 1 ), -1 )
256
257      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
258      ! Dynamics                                    (ta, sa used as workspace)
259      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
260      ! apply dynamics assimilation increment ... not needed in adjoint
261                             indic=0
262      ! vertical mesh at next time step ... not available
263                             CALL ssh_nxt_adj( kstp )         ! sea surface height at next time step
264                             CALL dyn_nxt_adj( kstp )           ! lateral velocity at next time step
265                             CALL dyn_spg_adj( kstp, indic )    ! surface pressure gradient
266                             CALL dyn_zdf_adj( kstp )           ! vertical diffusion
267                             CALL dyn_bfr_adj( kstp )           ! bottom friction
268                             CALL dyn_hpg_adj( kstp )           ! horizontal gradient of Hydrostatic pressure
269                             CALL dyn_ldf_adj( kstp )           ! lateral mixing
270                             CALL dyn_vor_adj( kstp )           ! vorticity term including Coriolis
271                             CALL dyn_adv_adj( kstp )           ! advection (vector or flux form)
272
273                             ua_ad(:,:,:)  = 0.0_wp                ! set tracer trends to zero
274                             va_ad(:,:,:)  = 0.0_wp
275
276      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
277      ! Active tracers                              (ua, va used as workspace)
278      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
279
280      IF( ln_dynhpg_imp  ) THEN                             ! semi-implicit hpg
281         IF( ln_zps    )     CALL zps_hde_adj( kstp, jpts, tsa, gtsu, gtsv, tsa_ad,      &   ! Partial steps: time filtered hor. gradient
282            &                                  rhd_ad, gtsu_ad, gru_ad, gtsv_ad,         &   ! of t, s, rd at the bottom ocean level
283            &                                  grv_ad )
284                             CALL eos_adj( tsa, tsa_ad, rhd_ad, rhop_ad )                    ! Time-filtered in situ density used in dynhpg module
285                             CALL tra_nxt_adj( kstp )                                        ! tracer fields at next time step
286      ELSE                                                                                   ! centered hpg (default case)
287                             CALL tra_nxt_adj( kstp )                                        ! tracer fields at next time step
288         IF( ln_zps    )     CALL zps_hde_adj( kstp, jpts, tsn, gtsu, gtsv, tsn_ad,      &   ! Partial steps: time filtered hor. gradient
289            &                                  rhd_ad, gtsu_ad, gru_ad, gtsv_ad,         &   ! of t, s, rd at the bottom ocean level
290            &                                  grv_ad )
291                             CALL eos_adj( tsn, tsn_ad, rhd_ad, rhop_ad )  ! now (swap=before) in situ density for dynhpg module
292      ENDIF
293
294      ! update the new (t,s) fields by non
295      ! penetrative convective adjustment ... not available
296
297                             CALL tra_zdf_adj( kstp )       ! vertical mixing
298                             CALL tra_ldf_adj( kstp )       ! lateral mixing
299                             CALL tra_adv_adj( kstp )      ! horizontal & vertical advection
300      IF( ln_tradmp      )   CALL tra_dmp_adj( kstp )       ! internal damping trends
301      !! advective (and/or diffusive) bottom boundary layer scheme ... currently not available
302      IF( lk_trabbl  )       CALL tra_bbl_adj( kstp )   ! diffusive bottom boundary layer scheme
303      IF( ln_trabbc      )   CALL tra_bbc_adj( kstp )       ! bottom heat flux
304
305      IF( ln_traqsr      )   CALL tra_qsr_adj( kstp )       ! penetrative solar radiation qsr
306
307                             CALL tra_sbc_adj( kstp )       ! surface boundary condition
308
309                             tsa_ad(:,:,:,:) = 0.0_wp            ! set tracer trends to zero
310
311      ! apply tracer assimilation increment ... not needed in adjoint
312#if defined key_top
313      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
314      ! Passive Tracer Model
315      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
316      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
317      !-----------------------------------------------------------------------
318
319      ! time-stepping... not needed in adjoint for the time being
320
321#endif
322      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
323      ! Ocean physics update                (ua, va, ta, sa used as workspace)
324      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
325      !-----------------------------------------------------------------------
326      !  LATERAL PHYSICS
327      !-----------------------------------------------------------------------
328#if defined key_traldf_c2d
329      ! eddy induced velocity coefficient... not needed in tangent with lk_zdfcst_adj
330#endif
331      ! before slope of the lateral mixing... not needed in adjoint with lk_zdfcst_adj
332      !-----------------------------------------------------------------------
333      !  VERTICAL PHYSICS
334      !-----------------------------------------------------------------------
335                             CALL zdf_bfr_adj ( kstp )
336      !! N.B. ua, va, ta, sa arrays are used as workspace in this section
337      !!-----------------------------------------------------------------------
338      !! mixed layer depth... not needed in adjoint with lk_zdfcst_adj
339      !! bottom friction... not needed in adjoint with lk_zdfcst_adj
340      !! double diffusive mixing ... not needed in adjoint with lk_zdfcst_adj
341      !! enhanced vertical eddy diffusivity ... not needed in adjoint with lk_zdfcst_adj
342      !! ! increase diffusivity at rivers mouths... not needed in tangent
343      !! lk_zdfcst_adj:  Constant Kz read from the reference trajectory
344      !! Constant Kz (reset avt, avm[uv] to the background value)... not available
345      !! KPP closure scheme for Kz ... not available
346      !! TKE closure scheme for Kz ... not available
347      !! Richardson number dependent Kz  ... not available
348      !!                                                                 ! Vertical eddy viscosity and diffusivity coefficients
349                             CALL bn2_adj( tsn, tsn_ad, rn2_ad )         ! now Brunt-Vaisala frequency
350                             CALL bn2_adj( tsb, tsb_ad, rn2b_ad )        ! now Brunt-Vaisala frequency
351      !!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
352      !!  Ocean dynamics : ssh, wn, hdiv, rot                                 !
353      !!<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
354                             CALL ssh_wzv_adj( kstp )         ! after ssh & vertical velocity
355      !!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
356      !! Update data, open boundaries and Forcings
357      !!<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
358      !! Output the initial state and forcings ... not needed in adjoint
359                             CALL sbc_adj ( kstp ) ! Sea boundary condition (including sea-ice)
360      ! update 3D salinity data ... not needed in tangent
361      ! update 3D temperature data ... not needed in adjoint
362
363      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
364      ! Control, and restarts
365      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
366      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
367      !-----------------------------------------------------------------------
368      !                                                     ! Time loop: control and print
369                             CALL stp_ctl_adj( kstp, indic, 0 )
370      IF( indic < 0      )   CALL ctl_stop( 'step_adj: indic < 0' )
371
372      ! close input  ocean restart file ... not needed in adjoint
373      ! write output ocean restart file... not needed in adjoint
374      ! write open boundary restart file... not needed in adjoint
375      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
376      ! diagnostics and outputs
377      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
378      ! N.B. ua, va, ta, sa arrays are used as workspace in this section
379      !-----------------------------------------------------------------------
380
381      IF ( nstop == 0 ) THEN                                ! Diagnostics
382         ! not needed in adjoint
383      ENDIF
384
385      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
386      ! Coupled mode
387      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
388
389      ! coupled mode : field exchanges ... not available for the next 20 years
390      !
391      !
392      IF( nn_timing == 1 .AND.  kstp == nit000  )   CALL timing_reset
393      !
394   END SUBROUTINE stp_adj
395
396   SUBROUTINE stp_adj_tst( kumadt )
397      !!-----------------------------------------------------------------------
398      !!
399      !!                  ***  ROUTINE example_adj_tst ***
400      !!
401      !! ** Purpose : Test the adjoint routine.
402      !!
403      !! ** Method  : Verify the scalar product
404      !!
405      !!                 ( L dx )^T W dy  =  dx^T L^T W dy
406      !!
407      !!              where  L   = tangent routine
408      !!                     L^T = adjoint routine
409      !!                     W   = diagonal matrix of scale factors
410      !!                     dx  = input perturbation (random field)
411      !!                     dy  = L dx
412      !!
413      !!
414      !! History :
415      !!        ! 09-03 (F. Vigilant)
416      !!-----------------------------------------------------------------------
417      !! * Modules used
418      USE sbcssr_tam, ONLY: &
419         & qrp_tl, qrp_ad,  &
420         & erp_tl, erp_ad
421      USE sbcfwb_tam, ONLY: &
422         & a_fwb_tl,        & ! for before year.
423         & a_fwb_ad  ! for before year.
424
425      !! * Arguments
426      INTEGER, INTENT(IN) :: &
427         & kumadt               ! Output unit
428
429      !! * Local declarations
430      INTEGER ::             &
431         & ji,               &  ! dummy loop indices
432         & jj,               &
433         & jk,               &
434         & istp
435      INTEGER, DIMENSION(jpi,jpj) :: &
436         & iseed_2d             ! 2D seed for the random number generator
437      REAL(KIND=wp) ::       &
438         & zsp1    ,         &  ! scalar product involving the tangent routine
439         & zsp1_U  ,         &
440         & zsp1_V  ,         &
441         & zsp1_T  ,         &
442         & zsp1_S  ,         &
443         & zsp1_SSH,         &
444         & zsp2    ,         &  ! scalar product involving the adjoint routine
445         & zsp2_U  ,         &
446         & zsp2_V  ,         &
447         & zsp2_T  ,         &
448         & zsp2_S  ,         &
449         & zsp2_SSH
450      REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
451         & zun_tlin ,        & ! Tangent input
452         & zvn_tlin ,        & ! Tangent input
453         & ztn_tlin ,        & ! Tangent input
454         & zsn_tlin ,        & ! Tangent input
455         & zun_tlout,        & ! Tangent output
456         & zvn_tlout,        & ! Tangent output
457         & ztn_tlout,        & ! Tangent output
458         & zsn_tlout,        & ! Tangent outpu
459         & zun_adin ,        & ! Adjoint input
460         & zvn_adin ,        & ! Adjoint input
461         & ztn_adin ,        & ! Adjoint input
462         & zsn_adin ,        & ! Adjoint input
463#if defined key_obc
464         & zub_adin ,        & ! Adjoint input
465         & zvb_adin ,        & ! Adjoint input
466         & ztb_adin ,        & ! Adjoint input
467         & zsb_adin ,        & ! Adjoint input
468         & zub_tlout ,        & ! Adjoint input
469         & zvb_tlout ,        & ! Adjoint input
470         & ztb_tlout ,        & ! Adjoint input
471         & zsb_tlout ,        & ! Adjoint input
472#endif
473         & zun_adout,        & ! Adjoint output
474         & zvn_adout,        & ! Adjoint output
475         & ztn_adout,        & ! Adjoint output
476         & zsn_adout,        & ! Adjoint output
477         & z3r                 ! 3D random field
478      REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
479         & zsshn_tlin  ,     & ! Tangent input
480         & zsshn_tlout ,     & ! Tangent input
481         & zsshn_adin  ,     & ! Adjoint output
482         & zsshn_adout ,     & ! Adjoint output
483#if defined key_obc
484         & zsshb_tlout ,     & ! Tangent input
485         & zsshb_adin  ,     & ! Tangent input
486#endif
487         & z2r                 ! 2D random field
488
489      CHARACTER(LEN=14) ::   &
490         & cl_name
491
492      INTEGER ::             &
493         & jpert
494      INTEGER, PARAMETER :: jpertmax = 6
495
496      ! Allocate memory
497
498      ALLOCATE( &
499         & zun_tlin(   jpi,jpj,jpk), zvn_tlin(   jpi,jpj,jpk), ztn_tlin(  jpi,jpj,jpk), &
500         & zsn_tlin(   jpi,jpj,jpk), zsshn_tlin( jpi,jpj    ), zun_tlout( jpi,jpj,jpk), &
501         & zvn_tlout(  jpi,jpj,jpk), ztn_tlout(  jpi,jpj,jpk), zsn_tlout( jpi,jpj,jpk), &
502         & zsshn_tlout(jpi,jpj    ), zun_adin(   jpi,jpj,jpk), zvn_adin(  jpi,jpj,jpk), &
503         & ztn_adin(   jpi,jpj,jpk), zsn_adin(   jpi,jpj,jpk), zsshn_adin(jpi,jpj    ), &
504         & zun_adout(  jpi,jpj,jpk), zvn_adout(  jpi,jpj,jpk), ztn_adout( jpi,jpj,jpk), &
505         & zsn_adout(  jpi,jpj,jpk), zsshn_adout(jpi,jpj    ), z2r(       jpi,jpj    ), &
506         & z3r(        jpi,jpj,jpk)                                                     &
507         & )
508
509#if defined key_obc
510      ALLOCATE( zub_adin   (jpi,jpj,jpk), zvb_adin  (jpi,jpj,jpk) ,        &
511           &    ztb_adin   (jpi,jpj,jpk), zsb_adin  (jpi,jpj,jpk) ,        &
512           &    zub_tlout  (jpi,jpj,jpk), zvb_tlout (jpi,jpj,jpk) ,        &
513           &    ztb_tlout  (jpi,jpj,jpk), zsb_tlout (jpi,jpj,jpk) ,        &
514           &    zsshb_tlout(jpi,jpj)    , zsshb_adin(jpi,jpj)       )
515#endif
516      !==================================================================
517      ! 1) dx = ( un_tl, vn_tl, tn_tl, sn_tl, sshn_tl ) and
518      !    dy = ( ..... )
519      !==================================================================
520
521
522      !DO jpert = 1, jpertmax
523      jpert = 6
524         !--------------------------------------------------------------------
525         ! Reset the tangent and adjoint variables
526         !--------------------------------------------------------------------
527         zun_tlin   (:,:,:) = 0.0_wp
528         zvn_tlin   (:,:,:) = 0.0_wp
529         ztn_tlin   (:,:,:) = 0.0_wp
530         zsn_tlin   (:,:,:) = 0.0_wp
531         zsshn_tlin (  :,:) = 0.0_wp
532         zun_tlout  (:,:,:) = 0.0_wp
533         zvn_tlout  (:,:,:) = 0.0_wp
534         ztn_tlout  (:,:,:) = 0.0_wp
535         zsn_tlout  (:,:,:) = 0.0_wp
536         zsshn_tlout(  :,:) = 0.0_wp
537         zun_adin   (:,:,:) = 0.0_wp
538         zvn_adin   (:,:,:) = 0.0_wp
539         ztn_adin   (:,:,:) = 0.0_wp
540         zsn_adin   (:,:,:) = 0.0_wp
541         zsshn_adin (  :,:) = 0.0_wp
542         zun_adout  (:,:,:) = 0.0_wp
543         zvn_adout  (:,:,:) = 0.0_wp
544         ztn_adout  (:,:,:) = 0.0_wp
545         zsn_adout  (:,:,:) = 0.0_wp
546         zsshn_adout(  :,:) = 0.0_wp
547         z3r        (:,:,:) = 0.0_wp
548         z2r        (  :,:) = 0.0_wp
549
550         !--------------------------------------------------------------------
551         ! Initialize the tangent input with random noise: dx
552         !--------------------------------------------------------------------
553
554         IF ( (jpert == 1) .OR. (jpert == jpertmax) ) THEN
555            CALL grid_random(  z3r, 'U', 0.0_wp, stdu )
556            DO jk = 1, jpk
557               DO jj = nldj, nlej
558                  DO ji = nldi, nlei
559                     zun_tlin(ji,jj,jk) = z3r(ji,jj,jk)
560                  END DO
561               END DO
562            END DO
563         ENDIF
564         IF ( (jpert == 2) .OR. (jpert == jpertmax) ) THEN
565            CALL grid_random(  z3r, 'V', 0.0_wp, stdv )
566            DO jk = 1, jpk
567               DO jj = nldj, nlej
568                  DO ji = nldi, nlei
569                     zvn_tlin(ji,jj,jk) = z3r(ji,jj,jk)
570                  END DO
571               END DO
572            END DO
573         ENDIF
574         IF ( (jpert == 3) .OR. (jpert == jpertmax) ) THEN
575            CALL grid_random(  z3r, 'T', 0.0_wp, stdt )
576            DO jk = 1, jpk
577               DO jj = nldj, nlej
578                  DO ji = nldi, nlei
579                     ztn_tlin(ji,jj,jk) = z3r(ji,jj,jk)
580                  END DO
581               END DO
582            END DO
583         ENDIF
584         IF ( (jpert == 4) .OR. (jpert == jpertmax) ) THEN
585            CALL grid_random(  z3r, 'T', 0.0_wp, stds )
586            DO jk = 1, jpk
587               DO jj = nldj, nlej
588                  DO ji = nldi, nlei
589                     zsn_tlin(ji,jj,jk) = z3r(ji,jj,jk)
590                  END DO
591               END DO
592            END DO
593         ENDIF
594         IF ( (jpert == 5) .OR. (jpert == jpertmax) ) THEN
595            CALL grid_random(  z2r, 'T', 0.0_wp, stdssh )
596            DO jj = nldj, nlej
597               DO ji = nldi, nlei
598                  zsshn_tlin(ji,jj) = z2r(ji,jj)
599               END DO
600            END DO
601         ENDIF
602
603         IF (jpert>0) THEN
604
605            !CALL     oce_tam_init( 1 )    ! allocate/initialize tl variables
606            !CALL sbc_oce_tam_init( 1 )
607            !CALL sol_oce_tam_init( 1 )
608            !CALL trc_oce_tam_init( 1 )
609
610            qrp_tl = 0.0_wp
611            erp_tl = 0.0_wp
612
613            emp_tl(:,:) = 0.0_wp
614#if defined key_tradmp
615            strdmp_tl = 0.0_wp
616            ttrdmp_tl = 0.0_wp
617#endif
618            a_fwb_tl = 0.0_wp
619
620            istp = nit000 - 1
621            CALL  trj_rea( istp, 1 )
622            !--------------------------------------------------------------------
623            ! Initialize the tangent variables: dy^* = W dy
624            !--------------------------------------------------------------------
625
626            un_tl  (:,:,:) = zun_tlin   (:,:,:)
627            vn_tl  (:,:,:) = zvn_tlin   (:,:,:)
628            tsn_tl  (:,:,:,jp_tem) = ztn_tlin   (:,:,:)
629            tsn_tl  (:,:,:,jp_sal) = zsn_tlin   (:,:,:)
630            sshn_tl(  :,:) = zsshn_tlin (  :,:)
631
632#if defined key_pomme_r025
633            IF ( (jpert == 5) .OR. (jpert == jpertmax) ) THEN
634               !DO ji = 1, jpi
635               !   DO jj = 1, jpj
636               !      sshn_tl(ji,jj) = cos( (2.*rpi)*(FLOAT(ji)-0.5)/FLOAT(jpi) + rpi/2. ) &
637               !                     * sin( (2.*rpi)*(FLOAT(jj)-0.5)/FLOAT(jpj) ) / 100.
638               !   END DO
639               !END DO
640
641               DO ji = 1, jpi
642                  DO jj = 1, jpj
643                     zsshn_tlin(ji,jj) = exp( -((float(ji)-float(jpi)/2.)/(float(jpi)/5.))**2 &
644                          -((float(jj)-float(jpj)/2.)/(float(jpj)/5.))**2 ) / 100. &
645                          * tmask(ji,jj,1)
646                  END DO
647               END DO
648               sshn_tl(:,:) = zsshn_tlin(:,:)
649            ENDIF
650#endif
651
652            !CALL     oce_tam_deallocate( 2 )    ! deallocate adj variables
653            !CALL sbc_oce_tam_deallocate( 2 )
654            !CALL sol_oce_tam_deallocate( 2 )
655
656            !-----------------------------------------------------------------------
657            !  Initialization of the dynamics and tracer fields for the tangent
658            !-----------------------------------------------------------------------
659
660            CALL istate_init_tan
661
662            DO istp = nit000, nitend, 1
663               CALL stp_tan( istp )
664            END DO
665
666            zun_tlout  ( :,:,:) = un_tl   (:,:,:)
667            zvn_tlout  ( :,:,:) = vn_tl   (:,:,:)
668            ztn_tlout  ( :,:,:) = tsn_tl   (:,:,:,jp_tem)
669            zsn_tlout  ( :,:,:) = tsn_tl   (:,:,:,jp_sal)
670            zsshn_tlout(   :,:) = sshn_tl (  :,:)
671
672#if defined key_obc
673            zub_tlout  (:,:,:) = ub_tl  (:,:,:)
674            zvb_tlout  (:,:,:) = vb_tl  (:,:,:)
675            ztb_tlout  (:,:,:) = tsb_tl  (:,:,:,jp_tem)
676            zsb_tlout  (:,:,:) = tsb_tl  (:,:,:,jp_sal)
677            zsshb_tlout(:,:)   = sshb_tl(:,:)
678#endif
679
680            !--------------------------------------------------------------------
681            ! Initialize the adjoint variables: dy^* = W dy
682            !--------------------------------------------------------------------
683
684            DO jk = 1, jpk
685               DO jj = nldj, nlej
686                  DO ji = nldi, nlei
687                     zun_adin(ji,jj,jk) = zun_tlout(ji,jj,jk) &
688                          &               * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) &
689                          &               * umask(ji,jj,jk) * wesp_u
690#if defined key_obc
691                     zub_adin(ji,jj,jk) = zub_tlout(ji,jj,jk) &
692                          &               * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) &
693                          &               * umask(ji,jj,jk) * wesp_u
694#endif
695                  END DO
696               END DO
697            END DO
698
699            DO jk = 1, jpk
700               DO jj = nldj, nlej
701                  DO ji = nldi, nlei
702                     zvn_adin(ji,jj,jk) = zvn_tlout(ji,jj,jk) &
703                          &               * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) &
704                          &               * vmask(ji,jj,jk) * wesp_u
705#if defined key_obc
706                     zvb_adin(ji,jj,jk) = zvb_tlout(ji,jj,jk) &
707                          &               * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) &
708                          &               * vmask(ji,jj,jk) * wesp_u
709#endif
710                  END DO
711               END DO
712            END DO
713
714            DO jk = 1, jpk
715               DO jj = nldj, nlej
716                  DO ji = nldi, nlei
717                     ztn_adin(ji,jj,jk) = ztn_tlout(ji,jj,jk) &
718                          &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
719                          &               * tmask(ji,jj,jk) * wesp_t(jk)
720#if defined key_obc
721                     ztb_adin(ji,jj,jk) = ztb_tlout(ji,jj,jk) &
722                          &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
723                          &               * tmask(ji,jj,jk) * wesp_t(jk)
724#endif
725                  END DO
726               END DO
727            END DO
728
729            DO jk = 1, jpk
730               DO jj = nldj, nlej
731                  DO ji = nldi, nlei
732                     zsn_adin(ji,jj,jk) = zsn_tlout(ji,jj,jk) &
733                          &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
734                          &               * tmask(ji,jj,jk) * wesp_s(jk)
735#if defined key_obc
736                     zsb_adin(ji,jj,jk) = zsb_tlout(ji,jj,jk) &
737                          &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
738                          &               * tmask(ji,jj,jk) * wesp_s(jk)
739#endif
740                  END DO
741               END DO
742            END DO
743
744            DO jj = nldj, nlej
745               DO ji = nldi, nlei
746                  zsshn_adin(ji,jj) = zsshn_tlout(ji,jj) &
747                       &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) &
748                       &               * tmask(ji,jj,1) * wesp_ssh
749#if defined key_obc
750                  zsshb_adin(ji,jj) = zsshb_tlout(ji,jj) &
751                       &               * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) &
752                       &               * tmask(ji,jj,1) * wesp_ssh
753#endif
754               END DO
755            END DO
756
757            !--------------------------------------------------------------------
758            ! Compute the scalar product: ( L dx )^T W dy
759            !--------------------------------------------------------------------
760
761            zsp1_U    = DOT_PRODUCT( zun_tlout  , zun_adin    )
762            zsp1_V    = DOT_PRODUCT( zvn_tlout  , zvn_adin    )
763            zsp1_T    = DOT_PRODUCT( ztn_tlout  , ztn_adin    )
764            zsp1_S    = DOT_PRODUCT( zsn_tlout  , zsn_adin    )
765            zsp1_SSH  = DOT_PRODUCT( zsshn_tlout, zsshn_adin  )
766
767            zsp1      = zsp1_U + zsp1_V + zsp1_T + zsp1_S + zsp1_SSH
768
769#if defined key_obc
770            zsp1_U    = DOT_PRODUCT( zub_tlout  , zub_adin    )
771            zsp1_V    = DOT_PRODUCT( zvb_tlout  , zvb_adin    )
772            zsp1_T    = DOT_PRODUCT( ztb_tlout  , ztb_adin    )
773            zsp1_S    = DOT_PRODUCT( zsb_tlout  , zsb_adin    )
774            zsp1_SSH  = DOT_PRODUCT( zsshb_tlout, zsshb_adin  )
775
776            zsp1      = zsp1 + ( zsp1_U + zsp1_V + zsp1_T + zsp1_S + zsp1_SSH )
777#endif
778            !--------------------------------------------------------------------
779            ! Call the adjoint routine: dx^* = L^T dy^*
780            !--------------------------------------------------------------------
781            !CALL     oce_tam_init( 2 )    ! allocate/initialize adj variables
782            !CALL sbc_oce_tam_init( 2 )
783            !CALL sol_oce_tam_init( 2 )
784            !CALL trc_oce_tam_init( 2 )
785
786            qrp_ad = 0.0_wp
787            erp_ad = 0.0_wp
788            emp_ad(:,:) = 0.0_wp
789#if defined key_tradmp
790            strdmp_ad = 0.0_wp
791            ttrdmp_ad = 0.0_wp
792#endif
793            a_fwb_ad = 0.0_wp
794
795            un_ad  (:,:,:) = zun_adin   (:,:,:)
796            vn_ad  (:,:,:) = zvn_adin   (:,:,:)
797            tsn_ad  (:,:,:,jp_tem) = ztn_adin   (:,:,:)
798            tsn_ad  (:,:,:,jp_sal) = zsn_adin   (:,:,:)
799            sshn_ad(  :,:) = zsshn_adin (  :,:)
800
801#if defined key_obc
802            ub_ad  (:,:,:) = zub_adin   (:,:,:)
803            vb_ad  (:,:,:) = zvb_adin   (:,:,:)
804            tsb_ad  (:,:,:,jp_tem) = ztb_adin   (:,:,:)
805            tsb_ad  (:,:,:,jp_sal) = zsb_adin   (:,:,:)
806            sshb_ad(:,:)   = zsshb_adin (:,:)
807#endif
808
809            !CALL     oce_tam_deallocate( 1 )    !deallocate tl variables
810            !CALL sbc_oce_tam_deallocate( 1 )
811            !CALL sol_oce_tam_deallocate( 1 )
812            !CALL trc_oce_tam_deallocate( 1 )
813
814            DO istp = nitend, nit000, -1
815               CALL stp_adj(istp)
816            END DO
817
818            CALL istate_init_adj
819
820            zun_adout  ( :,:,:) = un_ad   (:,:,:)
821            zvn_adout  ( :,:,:) = vn_ad   (:,:,:)
822            ztn_adout  ( :,:,:) = tsn_ad   (:,:,:,jp_tem)
823            zsn_adout  ( :,:,:) = tsn_ad   (:,:,:,jp_sal)
824            zsshn_adout(   :,:) = sshn_ad (  :,:)
825
826            zsp2_U    = DOT_PRODUCT( zun_tlin  , zun_adout )
827            zsp2_V    = DOT_PRODUCT( zvn_tlin  , zvn_adout )
828            zsp2_T    = DOT_PRODUCT( ztn_tlin  , ztn_adout )
829            zsp2_S    = DOT_PRODUCT( zsn_tlin  , zsn_adout )
830            zsp2_SSH  = DOT_PRODUCT( zsshn_tlin, zsshn_adout )
831
832            zsp2      = zsp2_U + zsp2_V + zsp2_T + zsp2_S + zsp2_SSH
833
834            !CALL     oce_tam_deallocate( 2 )    !deallocate adj variables
835            !CALL sbc_oce_tam_deallocate( 2 )
836            !CALL sol_oce_tam_deallocate( 2 )
837            !CALL trc_oce_tam_deallocate( 2 )
838
839            ! 14 char:'12345678901234'
840            SELECT CASE (jpert)
841            CASE(1)
842               cl_name = 'step_adj     U'
843            CASE(2)
844               cl_name = 'step_adj     V'
845            CASE(3)
846               cl_name = 'step_adj     T'
847            CASE(4)
848               cl_name = 'step_adj     S'
849            CASE(5)
850               cl_name = 'step_adj   SSH'
851            CASE(jpertmax)
852               cl_name = 'step_adj      '
853            END SELECT
854            CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
855         endif
856      !END DO
857
858      DEALLOCATE(                                 &
859         & zun_tlin   , zvn_tlin   , ztn_tlin   , &
860         & zsn_tlin   , zsshn_tlin , zun_tlout  , &
861         & zvn_tlout  , ztn_tlout  , zsn_tlout  , &
862         & zsshn_tlout, zun_adin   , zvn_adin   , &
863         & ztn_adin   , zsn_adin   , zsshn_adin , &
864         & zun_adout  , zvn_adout  , ztn_adout  , &
865         & zsn_adout  , zsshn_adout, z2r        , &
866         & z3r                                    &
867         & )
868
869   END SUBROUTINE stp_adj_tst
870
871      !!----------------------------------------------------------------------
872
873   !!======================================================================
874#endif
875END MODULE step_tam
Note: See TracBrowser for help on using the repository browser.