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.

traldf_lap_tam.F90 in branches/TAM_V3_0/NEMOTAM/OPATAM_SRC/TRA – NEMO

Context Navigation

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

Last change on this file since 1885 was 1885, checked in by rblod, 14 years ago
add TAM sources
File size: 63.4 KB

Line
1	MODULE traldf_lap_tam
2	#if defined key_tam
3	!!==============================================================================
4	!! * MODULE traldf_lap *
5	!! Ocean active tracers: horizontal component of the lateral tracer mixing trend
6	!! Tangent and adjoint module
7	!!==============================================================================
8	!! History : 9.0 ! ?????
9	!! History of the T&A module:
10	!! 9.0 ! 09-03 (F. Vigilant) original version
11	!!----------------------------------------------------------------------
12
13	!!----------------------------------------------------------------------
14	!! tra_ldf_lap : update the tracer trend with the horizontal diffusion
15	!! using a iso-level harmonic (laplacien) operator.
16	!!----------------------------------------------------------------------
17	!! * Modules used
18	USE par_kind , ONLY: & ! Precision variables
19	& wp
20	USE par_oce , ONLY: & ! Ocean space and time domain variables
21	& jpiglo, &
22	& jpi, &
23	& jpj, &
24	& jpk, &
25	& jpim1, &
26	& jpjm1, &
27	& jpkm1
28	USE oce_tam , ONLY: & ! ocean dynamics and active tracers
29	& tb_tl, &
30	& sb_tl, &
31	& ta_tl, &
32	& sa_tl, &
33	& gtu_tl, &
34	& gsu_tl, &
35	& gtv_tl, &
36	& gsv_tl, &
37	& tb_ad, &
38	& sb_ad, &
39	& ta_ad, &
40	& sa_ad, &
41	& gtu_ad, &
42	& gsu_ad, &
43	& gtv_ad, &
44	& gsv_ad
45	USE dom_oce , ONLY: & ! ocean space and time domain
46	& e1u, &
47	& e2u, &
48	& e1v, &
49	& e2v, &
50	& e1t, &
51	& e2t, &
52	#if defined key_zco
53	& e3t_0, &
54	#else
55	& e3u, &
56	& e3v, &
57	& e3t, &
58	#endif
59	& tmask, &
60	& umask, &
61	& vmask, &
62	& mbathy, &
63	& ln_zps, &
64	& mig, &
65	& mjg, &
66	& nldi, &
67	& nldj, &
68	& nlei, &
69	& nlej, &
70	& lk_vvl
71	USE ldftra_oce , ONLY: & ! ocean active tracers: lateral physics
72	& ahtv, &
73	& ahtu, &
74	& ahtw, &
75	& ahtb0, &
76	& aht0
77	! USE trdmod ! ocean active tracers trends
78	! USE trdmod_oce ! ocean variables trends
79	USE in_out_manager, ONLY: & ! I/O manager
80	& lwp, &
81	& numout, &
82	& nitend, &
83	& nit000, &
84	& nstop
85	! USE diaptr ! poleward transport diagnostics
86	! USE prtctl ! Print control
87
88	USE gridrandom , ONLY: & ! Random Gaussian noise on grids
89	& grid_random
90	USE dotprodfld , ONLY: & ! Computes dot product for 3D and 2D fields
91	& dot_product
92	USE tstool_tam , ONLY: &
93	& prntst_adj, & !
94	& prntst_tlm, & !
95	& stdt, & ! stdev for u-velocity
96	& stds ! v-velocity
97	USE paresp , ONLY: & ! Weights for an energy-type scalar product
98	& wesp_t, &
99	& wesp_s
100
101	IMPLICIT NONE
102	PRIVATE
103
104	!! * Routine accessibility
105	PUBLIC tra_ldf_lap_tan ! routine called by tradldf_tam.F90
106	PUBLIC tra_ldf_lap_adj ! routine called by tradldf_tam.F90
107	PUBLIC tra_ldf_lap_adj_tst ! routine called by tradldf_tam.F90
108	PUBLIC tra_ldf_lap_tlm_tst
109
110	!! * Substitutions
111	# include "domzgr_substitute.h90"
112	# include "ldftra_substitute.h90"
113	# include "vectopt_loop_substitute.h90"
114	!!----------------------------------------------------------------------
115	!! OPA 9.0 , LOCEAN-IPSL (2005)
116	!! $Id: traldf_lap.F90 1152 2008-06-26 14:11:13Z rblod $
117	!! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
118	!!----------------------------------------------------------------------
119
120	CONTAINS
121
122	SUBROUTINE tra_ldf_lap_tan( kt )
123	!!----------------------------------------------------------------------
124	!! * ROUTINE tra_ldf_lap_tan *
125	!!
126	!! ** Purpose : Compute the before horizontal tracer (t & s) diffusive
127	!! trend and add it to the general trend of tracer equation.
128	!!
129	!! ** Method : Second order diffusive operator evaluated using before
130	!! fields (forward time scheme). The horizontal diffusive trends of
131	!! temperature (idem for salinity) is given by:
132	!! difft = 1/(e1te2te3t) { di-1[ aht e2u*e3u/e1u di(tb) ]
133	!! + dj-1[ aht e1v*e3v/e2v dj(tb) ] }
134	!! Note: key_zco defined, the e3t=e3u=e3v, the trend becomes:
135	!! difft = 1/(e1t*e2t) { di-1[ aht e2u/e1u di(tb) ]
136	!! + dj-1[ aht e1v/e2v dj(tb) ] }
137	!! Add this trend to the general tracer trend (ta,sa):
138	!! (ta,sa) = (ta,sa) + ( difft , diffs )
139	!!
140	!! ** Action : - Update (ta,sa) arrays with the before iso-level
141	!! harmonic mixing trend.
142	!!
143	!! History :
144	!! 1.0 ! 87-06 (P. Andrich, D. L Hostis) Original code
145	!! ! 91-11 (G. Madec)
146	!! ! 95-11 (G. Madec) suppress volumetric scale factors
147	!! ! 96-01 (G. Madec) statement function for e3
148	!! 8.5 ! 02-06 (G. Madec) F90: Free form and module
149	!! 9.0 ! 04-08 (C. Talandier) New trends organization
150	!! ! 05-11 (G. Madec) add zps case
151	!! History of the tangent routine
152	!! 9.0 ! 03-09 (F. Vigilant) tangent of 9.0
153	!!----------------------------------------------------------------------
154
155	!! * Arguments
156	INTEGER, INTENT( in ) :: kt ! ocean time-step index
157
158	!! * Local save
159	REAL(wp), DIMENSION(jpi,jpj), SAVE :: &
160	& ze1ur, ze2vr, zbtr2 ! scale factor coefficients
161
162	!! * Local declarations
163	INTEGER :: ji, jj, jk ! dummy loop indices
164	INTEGER :: iku, ikv ! temporary integers
165	REAL(wp) :: zabe1, zabe2, zbtr ! temporary scalars
166	REAL(wp) :: ztatl, zsatl ! temporary scalars
167	REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztutl, ztvtl ! 3D workspace
168	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsutl, zsvtl ! 3D workspace
169	!!----------------------------------------------------------------------
170
171	IF( kt == nit000 ) THEN
172	IF(lwp) WRITE(numout,*)
173	IF(lwp) WRITE(numout,*) 'tra_ldf_lap_tan : iso-level laplacian diffusion'
174	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
175	ze1ur(:,:) = e2u(:,:) / e1u(:,:)
176	ze2vr(:,:) = e1v(:,:) / e2v(:,:)
177	zbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) )
178	ENDIF
179
180	! ! =============
181	DO jk = 1, jpkm1 ! Vertical slab
182	! ! =============
183	! 1. First derivative (gradient)
184	! -------------------
185	DO jj = 1, jpjm1
186	DO ji = 1, fs_jpim1 ! vector opt.
187	#if defined key_zco
188	zabe1 = fsahtu(ji,jj,jk) * umask(ji,jj,jk) * ze1ur(ji,jj)
189	zabe2 = fsahtv(ji,jj,jk) * vmask(ji,jj,jk) * ze2vr(ji,jj)
190	#else
191	zabe1 = fsahtu(ji,jj,jk) * umask(ji,jj,jk) * ze1ur(ji,jj) * fse3u(ji,jj,jk)
192	zabe2 = fsahtv(ji,jj,jk) * vmask(ji,jj,jk) * ze2vr(ji,jj) * fse3v(ji,jj,jk)
193	#endif
194	ztutl(ji,jj,jk) = zabe1 * ( tb_tl(ji+1,jj ,jk) - tb_tl(ji,jj,jk) )
195	zsutl(ji,jj,jk) = zabe1 * ( sb_tl(ji+1,jj ,jk) - sb_tl(ji,jj,jk) )
196	ztvtl(ji,jj,jk) = zabe2 * ( tb_tl(ji ,jj+1,jk) - tb_tl(ji,jj,jk) )
197	zsvtl(ji,jj,jk) = zabe2 * ( sb_tl(ji ,jj+1,jk) - sb_tl(ji,jj,jk) )
198	END DO
199	END DO
200	IF( ln_zps ) THEN ! set gradient at partial step level
201	DO jj = 1, jpjm1
202	DO ji = 1, fs_jpim1 ! vector opt.
203	! last level
204	iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj ) ) - 1
205	ikv = MIN ( mbathy(ji,jj), mbathy(ji ,jj+1) ) - 1
206	IF( iku == jk ) THEN
207	zabe1 = fsahtu(ji,jj,iku) * umask(ji,jj,iku) * ze1ur(ji,jj) * fse3u(ji,jj,iku)
208	ztutl(ji,jj,jk) = zabe1 * gtu_tl(ji,jj)
209	zsutl(ji,jj,jk) = zabe1 * gsu_tl(ji,jj)
210	ENDIF
211	IF( ikv == jk ) THEN
212	zabe2 = fsahtv(ji,jj,ikv) * vmask(ji,jj,ikv) * ze2vr(ji,jj) * fse3v(ji,jj,ikv)
213	ztvtl(ji,jj,jk) = zabe2 * gtv_tl(ji,jj)
214	zsvtl(ji,jj,jk) = zabe2 * gsv_tl(ji,jj)
215	ENDIF
216	END DO
217	END DO
218	ENDIF
219
220
221	! 2. Second derivative (divergence)
222	! --------------------
223	DO jj = 2, jpjm1
224	DO ji = fs_2, fs_jpim1 ! vector opt.
225	#if defined key_zco
226	zbtr = zbtr2(ji,jj)
227	#else
228	zbtr = zbtr2(ji,jj) / fse3t(ji,jj,jk)
229	#endif
230	! horizontal diffusive trends
231	ztatl = zbtr * ( ztutl(ji,jj,jk) - ztutl(ji-1,jj ,jk) &
232	& + ztvtl(ji,jj,jk) - ztvtl(ji ,jj-1,jk) )
233	zsatl = zbtr * ( zsutl(ji,jj,jk) - zsutl(ji-1,jj ,jk) &
234	& + zsvtl(ji,jj,jk) - zsvtl(ji ,jj-1,jk) )
235	! add it to the general tracer trends
236	ta_tl(ji,jj,jk) = ta_tl(ji,jj,jk) + ztatl
237	sa_tl(ji,jj,jk) = sa_tl(ji,jj,jk) + zsatl
238	END DO
239	END DO
240	! ! =============
241	END DO ! End of slab
242	! ! =============
243
244	! "zonal" mean lateral diffusive heat and salt transport
245	! Warning: no update of pht_ldf and /pst_ldf/ when /ln_diaptr/ is activated
246
247	END SUBROUTINE tra_ldf_lap_tan
248
249
250	SUBROUTINE tra_ldf_lap_adj( kt )
251	!!----------------------------------------------------------------------
252	!! * ROUTINE tra_ldf_lap_adj *
253	!!
254	!! ** Purpose : Compute the before horizontal tracer (t & s) diffusive
255	!! trend and add it to the general trend of tracer equation.
256	!!
257	!! ** Method : Second order diffusive operator evaluated using before
258	!! fields (forward time scheme). The horizontal diffusive trends of
259	!! temperature (idem for salinity) is given by:
260	!! difft = 1/(e1te2te3t) { di-1[ aht e2u*e3u/e1u di(tb) ]
261	!! + dj-1[ aht e1v*e3v/e2v dj(tb) ] }
262	!! Note: key_zco defined, the e3t=e3u=e3v, the trend becomes:
263	!! difft = 1/(e1t*e2t) { di-1[ aht e2u/e1u di(tb) ]
264	!! + dj-1[ aht e1v/e2v dj(tb) ] }
265	!! Add this trend to the general tracer trend (ta,sa):
266	!! (ta,sa) = (ta,sa) + ( difft , diffs )
267	!!
268	!! ** Action : - Update (ta,sa) arrays with the before iso-level
269	!! harmonic mixing trend.
270	!!
271	!! History :
272	!! 1.0 ! 87-06 (P. Andrich, D. L Hostis) Original code
273	!! ! 91-11 (G. Madec)
274	!! ! 95-11 (G. Madec) suppress volumetric scale factors
275	!! ! 96-01 (G. Madec) statement function for e3
276	!! 8.5 ! 02-06 (G. Madec) F90: Free form and module
277	!! 9.0 ! 04-08 (C. Talandier) New trends organization
278	!! ! 05-11 (G. Madec) add zps case
279	!! History of the tangent routine
280	!! 9.0 ! 03-09 (F. Vigilant) adjoint of 9.0
281	!!----------------------------------------------------------------------
282
283	!! * Arguments
284	INTEGER, INTENT( in ) :: kt ! ocean time-step index
285
286	!! * Local save
287	REAL(wp), DIMENSION(jpi,jpj), SAVE :: &
288	& ze1ur, ze2vr, zbtr2 ! scale factor coefficients
289
290	!! * Local declarations
291	INTEGER :: ji, jj, jk ! dummy loop indices
292	INTEGER :: iku, ikv ! temporary integers
293	REAL(wp) :: zabe1, zabe2, zbtr ! temporary scalars
294	REAL(wp) :: ztaad, zsaad ! temporary scalars
295	REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztuad, ztvad ! 3D workspace
296	REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsuad, zsvad ! 3D workspace
297	!!----------------------------------------------------------------------
298	ztvad(:,:,:) = 0.0_wp
299	zsvad(:,:,:) = 0.0_wp
300	ztuad(:,:,:) = 0.0_wp
301	zsuad(:,:,:) = 0.0_wp
302	ztaad = 0.0_wp
303	zsaad = 0.0_wp
304
305	IF( kt == nit000 ) THEN
306	IF(lwp) WRITE(numout,*)
307	IF(lwp) WRITE(numout,*) 'tra_ldf_lap_adj : iso-level laplacian diffusion'
308	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
309	ze1ur(:,:) = e2u(:,:) / e1u(:,:)
310	ze2vr(:,:) = e1v(:,:) / e2v(:,:)
311	zbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) )
312	ENDIF
313
314
315
316	! ! =============
317	DO jk = 1, jpkm1 ! Vertical slab
318	! ! =============
319
320	! 2. Second derivative (divergence)
321	! --------------------
322	DO jj = jpjm1, 2, -1
323	DO ji = fs_jpim1, fs_2, -1 ! vector opt.
324	#if defined key_zco
325	zbtr = zbtr2(ji,jj)
326	#else
327	zbtr = zbtr2(ji,jj) / fse3t(ji,jj,jk)
328	#endif
329	! horizontal diffusive trends
330	ztaad = zbtr * ta_ad(ji,jj,jk)
331	zsaad = zbtr * sa_ad(ji,jj,jk)
332
333	ztuad(ji ,jj ,jk) = ztuad(ji ,jj ,jk) + ztaad
334	ztuad(ji-1,jj ,jk) = ztuad(ji-1,jj ,jk) - ztaad
335	ztvad(ji ,jj ,jk) = ztvad(ji ,jj ,jk) + ztaad
336	ztvad(ji ,jj-1,jk) = ztvad(ji ,jj-1,jk) - ztaad
337
338	zsuad(ji ,jj ,jk) = zsuad(ji ,jj ,jk) + zsaad
339	zsuad(ji-1,jj ,jk) = zsuad(ji-1,jj ,jk) - zsaad
340	zsvad(ji ,jj ,jk) = zsvad(ji ,jj ,jk) + zsaad
341	zsvad(ji ,jj-1,jk) = zsvad(ji ,jj-1,jk) - zsaad
342
343	!ztaad = 0.0_wp
344	!zsaad = 0.0_wp
345
346	END DO
347	END DO
348
349	! 1. First derivative (gradient)
350	! -------------------
351
352	IF( ln_zps ) THEN ! set gradient at partial step level
353	DO jj = jpjm1, 1, -1
354	DO ji = fs_jpim1, 1, -1 ! vector opt.?
355	! last level
356	iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj ) ) - 1
357	ikv = MIN ( mbathy(ji,jj), mbathy(ji ,jj+1) ) - 1
358	IF( iku == jk ) THEN
359	zabe1 = fsahtu(ji,jj,iku) * umask(ji,jj,iku) * ze1ur(ji,jj) * fse3u(ji,jj,iku)
360
361	gtu_ad(ji,jj) = zabe1 * ztuad(ji,jj,jk) + gtu_ad(ji,jj)
362	gsu_ad(ji,jj) = zabe1 * zsuad(ji,jj,jk) + gsu_ad(ji,jj)
363
364	ztuad(ji,jj,jk) = 0.0_wp
365	zsuad(ji,jj,jk) = 0.0_wp
366	ENDIF
367	IF( ikv == jk ) THEN
368	zabe2 = fsahtv(ji,jj,ikv) * vmask(ji,jj,ikv) * ze2vr(ji,jj) * fse3v(ji,jj,ikv)
369
370	gtv_ad(ji,jj) = zabe2 * ztvad(ji,jj,jk) + gtv_ad(ji,jj)
371	gsv_ad(ji,jj) = zabe2 * zsvad(ji,jj,jk) + gsv_ad(ji,jj)
372
373	ztvad(ji,jj,jk) = 0.0_wp
374	zsvad(ji,jj,jk) = 0.0_wp
375	ENDIF
376	END DO
377	END DO
378	ENDIF
379
380	DO jj = jpjm1, 1, -1
381	DO ji = fs_jpim1, 1, -1 ! vector opt. ?
382	#if defined key_zco
383	zabe1 = fsahtu(ji,jj,jk) * umask(ji,jj,jk) * ze1ur(ji,jj)
384	zabe2 = fsahtv(ji,jj,jk) * vmask(ji,jj,jk) * ze2vr(ji,jj)
385	#else
386	zabe1 = fsahtu(ji,jj,jk) * umask(ji,jj,jk) * ze1ur(ji,jj) * fse3u(ji,jj,jk)
387	zabe2 = fsahtv(ji,jj,jk) * vmask(ji,jj,jk) * ze2vr(ji,jj) * fse3v(ji,jj,jk)
388	#endif
389
390	tb_ad(ji ,jj ,jk) = tb_ad(ji ,jj ,jk) - (zabe1 * ztuad(ji,jj,jk) + zabe2 * ztvad(ji,jj,jk))
391	tb_ad(ji+1,jj ,jk) = tb_ad(ji+1,jj ,jk) + zabe1 * ztuad(ji,jj,jk)
392	tb_ad(ji ,jj+1,jk) = tb_ad(ji ,jj+1,jk) + zabe2 * ztvad(ji,jj,jk)
393
394	sb_ad(ji ,jj ,jk) = sb_ad(ji ,jj ,jk) - (zabe1 * zsuad(ji,jj,jk) + zabe2 * zsvad(ji,jj,jk))
395	sb_ad(ji+1,jj ,jk) = sb_ad(ji+1,jj ,jk) + zabe1 * zsuad(ji,jj,jk)
396	sb_ad(ji ,jj+1,jk) = sb_ad(ji ,jj+1,jk) + zabe2 * zsvad(ji,jj,jk)
397
398	ztuad(ji,jj,jk) = 0.0_wp
399	zsuad(ji,jj,jk) = 0.0_wp
400	ztvad(ji,jj,jk) = 0.0_wp
401	zsvad(ji,jj,jk) = 0.0_wp
402
403	END DO
404	END DO
405
406	! ! =============
407	END DO ! End of slab
408	! ! =============
409
410	END SUBROUTINE tra_ldf_lap_adj
411
412
413	SUBROUTINE tra_ldf_lap_adj_tst ( kumadt )
414	!!-----------------------------------------------------------------------
415	!!
416	!! * ROUTINE example_adj_tst *
417	!!
418	!! ** Purpose : Test the adjoint routine.
419	!!
420	!! ** Method : Verify the scalar product
421	!!
422	!! ( L dx )^T W dy = dx^T L^T W dy
423	!!
424	!! where L = tangent routine
425	!! L^T = adjoint routine
426	!! W = diagonal matrix of scale factors
427	!! dx = input perturbation (random field)
428	!! dy = L dx
429	!!
430	!! History :
431	!! ! 08-08 (A. Vidard)
432	!!-----------------------------------------------------------------------
433	!! * Modules used
434
435	!! * Arguments
436	INTEGER, INTENT(IN) :: &
437	& kumadt ! Output unit
438
439	!! * Local declarations
440	INTEGER :: &
441	& ji, & ! dummy loop indices
442	& jj, &
443	& jk
444	INTEGER, DIMENSION(jpi,jpj) :: &
445	& iseed_2d ! 2D seed for the random number generator
446	REAL(KIND=wp) :: &
447	& zsp1, & ! scalar product involving the tangent routine
448	& zsp1_T, &
449	& zsp1_S, &
450	& zsp2, & ! scalar product involving the adjoint routine
451	& zsp2_1, &
452	& zsp2_2, &
453	& zsp2_3, &
454	& zsp2_4, &
455	& zsp2_5, &
456	& zsp2_6, &
457	& zsp2_7, &
458	& zsp2_8, &
459	& zsp2_T, &
460	& zsp2_S
461	REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
462	& ztb_tlin , & ! Tangent input
463	& zsb_tlin , & ! Tangent input
464	& zta_tlin , & ! Tangent input
465	& zsa_tlin , & ! Tangent input
466	& zta_tlout, & ! Tangent output
467	& zsa_tlout, & ! Tangent output
468	& zta_adin, & ! Adjoint input
469	& zsa_adin, & ! Adjoint input
470	& ztb_adout , & ! Adjoint output
471	& zsb_adout , & ! Adjoint output
472	& zta_adout , & ! Adjoint output
473	& zsa_adout , & ! Adjoint output
474	& z3r ! 3D random field
475	REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
476	& zgtu_tlin , & ! Tangent input
477	& zgsu_tlin , & ! Tangent input
478	& zgtv_tlin , & ! Tangent input
479	& zgsv_tlin , & ! Tangent input
480	& zgtu_adout , & ! Adjoint output
481	& zgsu_adout , & ! Adjoint output
482	& zgtv_adout , & ! Adjoint output
483	& zgsv_adout , & ! Adjoint output
484	& z2r ! 2D random field
485	CHARACTER(LEN=14) :: cl_name
486	! Allocate memory
487
488	ALLOCATE( &
489	& ztb_tlin(jpi,jpj,jpk), &
490	& zsb_tlin(jpi,jpj,jpk), &
491	& zta_tlin(jpi,jpj,jpk), &
492	& zsa_tlin(jpi,jpj,jpk), &
493	& zgtu_tlin(jpi,jpj), &
494	& zgsu_tlin(jpi,jpj), &
495	& zgtv_tlin(jpi,jpj), &
496	& zgsv_tlin(jpi,jpj), &
497	& zta_tlout(jpi,jpj,jpk), &
498	& zsa_tlout(jpi,jpj,jpk), &
499	& zta_adin(jpi,jpj,jpk), &
500	& zsa_adin(jpi,jpj,jpk), &
501	& ztb_adout(jpi,jpj,jpk), &
502	& zsb_adout(jpi,jpj,jpk), &
503	& zta_adout(jpi,jpj,jpk), &
504	& zsa_adout(jpi,jpj,jpk), &
505	& zgtu_adout(jpi,jpj), &
506	& zgsu_adout(jpi,jpj), &
507	& zgtv_adout(jpi,jpj), &
508	& zgsv_adout(jpi,jpj), &
509	& z3r(jpi,jpj,jpk), &
510	& z2r(jpi,jpj) &
511	& )
512
513
514	!=======================================================================
515	! 1) dx = ( tb_tl, ta_tl, sb_tl, sa_tl, gtu_tl, gtv_tl, gsu_tl, gsv_tl )
516	! dy = ( ta_tl, sa_tl )
517	!=======================================================================
518
519	!--------------------------------------------------------------------
520	! Reset the tangent and adjoint variables
521	!--------------------------------------------------------------------
522
523	ztb_tlin(:,:,:) = 0.0_wp
524	zsb_tlin(:,:,:) = 0.0_wp
525	zta_tlin(:,:,:) = 0.0_wp
526	zsa_tlin(:,:,:) = 0.0_wp
527	zgtu_tlin(:,:) = 0.0_wp
528	zgsu_tlin(:,:) = 0.0_wp
529	zgtv_tlin(:,:) = 0.0_wp
530	zgsv_tlin(:,:) = 0.0_wp
531	zta_tlout(:,:,:) = 0.0_wp
532	zsa_tlout(:,:,:) = 0.0_wp
533	zta_adin(:,:,:) = 0.0_wp
534	zsa_adin(:,:,:) = 0.0_wp
535	ztb_adout(:,:,:) = 0.0_wp
536	zsb_adout(:,:,:) = 0.0_wp
537	zta_adout(:,:,:) = 0.0_wp
538	zsa_adout(:,:,:) = 0.0_wp
539	zgtu_adout(:,:) = 0.0_wp
540	zgsu_adout(:,:) = 0.0_wp
541	zgtv_adout(:,:) = 0.0_wp
542	zgsv_adout(:,:) = 0.0_wp
543
544	tb_tl(:,:,:) = 0.0_wp
545	sb_tl(:,:,:) = 0.0_wp
546	ta_tl(:,:,:) = 0.0_wp
547	sa_tl(:,:,:) = 0.0_wp
548	gtu_tl(:,:) = 0.0_wp
549	gsu_tl(:,:) = 0.0_wp
550	gtv_tl(:,:) = 0.0_wp
551	gsv_tl(:,:) = 0.0_wp
552	tb_ad(:,:,:) = 0.0_wp
553	sb_ad(:,:,:) = 0.0_wp
554	ta_ad(:,:,:) = 0.0_wp
555	sa_ad(:,:,:) = 0.0_wp
556	gtu_ad(:,:) = 0.0_wp
557	gsu_ad(:,:) = 0.0_wp
558	gtv_ad(:,:) = 0.0_wp
559	gsv_ad(:,:) = 0.0_wp
560
561	!--------------------------------------------------------------------
562	! Initialize the tangent input with random noise: dx
563	!--------------------------------------------------------------------
564
565	DO jj = 1, jpj
566	DO ji = 1, jpi
567	iseed_2d(ji,jj) = - ( 596035 + &
568	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
569	END DO
570	END DO
571	CALL grid_random( iseed_2d, z3r, 'T', 0.0_wp, stdt )
572	ztb_tlin(:,:,:) = z3r(:,:,:)
573
574	DO jj = 1, jpj
575	DO ji = 1, jpi
576	iseed_2d(ji,jj) = - ( 727391 + &
577	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
578	END DO
579	END DO
580	CALL grid_random( iseed_2d, z3r, 'T', 0.0_wp, stds )
581	zsb_tlin(:,:,:) = z3r(:,:,:)
582
583	DO jj = 1, jpj
584	DO ji = 1, jpi
585	iseed_2d(ji,jj) = - ( 249741 + &
586	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
587	END DO
588	END DO
589	CALL grid_random( iseed_2d, z3r, 'T', 0.0_wp, stdt )
590	zta_tlin(:,:,:) = z3r(:,:,:)
591
592	DO jj = 1, jpj
593	DO ji = 1, jpi
594	iseed_2d(ji,jj) = - ( 182029 + &
595	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
596	END DO
597	END DO
598	CALL grid_random( iseed_2d, z3r, 'T', 0.0_wp, stds )
599	zsa_tlin(:,:,:) = z3r(:,:,:)
600
601	DO jj = 1, jpj
602	DO ji = 1, jpi
603	iseed_2d(ji,jj) = - ( 596035 + &
604	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
605	END DO
606	END DO
607	CALL grid_random( iseed_2d, z2r, 'U', 0.0_wp, stds )
608	zgtu_tlin(:,:) = z2r(:,:)
609
610	DO jj = 1, jpj
611	DO ji = 1, jpi
612	iseed_2d(ji,jj) = - ( 727391 + &
613	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
614	END DO
615	END DO
616	CALL grid_random( iseed_2d, z2r, 'U', 0.0_wp, stds )
617	zgsu_tlin(:,:) = z2r(:,:)
618
619	DO jj = 1, jpj
620	DO ji = 1, jpi
621	iseed_2d(ji,jj) = - ( 249741 + &
622	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
623	END DO
624	END DO
625	CALL grid_random( iseed_2d, z2r, 'V', 0.0_wp, stds )
626	zgtv_tlin(:,:) = z2r(:,:)
627
628	DO jj = 1, jpj
629	DO ji = 1, jpi
630	iseed_2d(ji,jj) = - ( 182029 + &
631	& mig(ji) + ( mjg(jj) - 1 ) * jpiglo )
632	END DO
633	END DO
634	CALL grid_random( iseed_2d, z2r, 'V', 0.0_wp, stds )
635	zgsv_tlin(:,:) = z2r(:,:)
636
637	tb_tl(:,:,:) = ztb_tlin(:,:,:)
638	sb_tl(:,:,:) = zsb_tlin(:,:,:)
639	ta_tl(:,:,:) = zta_tlin(:,:,:)
640	sa_tl(:,:,:) = zsa_tlin(:,:,:)
641	gtu_tl(:,:) = zgtu_tlin(:,:)
642	gsu_tl(:,:) = zgsu_tlin(:,:)
643	gtv_tl(:,:) = zgtv_tlin(:,:)
644	gsv_tl(:,:) = zgsv_tlin(:,:)
645
646	CALL tra_ldf_lap_tan( nit000 )
647
648	zta_tlout(:,:,:) = ta_tl(:,:,:)
649	zsa_tlout(:,:,:) = sa_tl(:,:,:)
650
651	!--------------------------------------------------------------------
652	! Initialize the adjoint variables: dy^* = W dy
653	!--------------------------------------------------------------------
654
655	DO jk = 1, jpk
656	DO jj = nldj, nlej
657	DO ji = nldi, nlei
658	zsa_adin(ji,jj,jk) = zsa_tlout(ji,jj,jk) &
659	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
660	& * tmask(ji,jj,jk) * wesp_s(jk)
661	zta_adin(ji,jj,jk) = zta_tlout(ji,jj,jk) &
662	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
663	& * tmask(ji,jj,jk) * wesp_t(jk)
664	END DO
665	END DO
666	END DO
667
668	!--------------------------------------------------------------------
669	! Compute the scalar product: ( L dx )^T W dy
670	!-------------------------------------------------------------------
671
672	zsp1_T = DOT_PRODUCT( zta_tlout, zta_adin )
673	zsp1_S = DOT_PRODUCT( zsa_tlout, zsa_adin )
674	zsp1 = zsp1_T + zsp1_S
675
676	!--------------------------------------------------------------------
677	! Call the adjoint routine: dx^* = L^T dy^*
678	!--------------------------------------------------------------------
679
680	ta_ad(:,:,:) = zta_adin(:,:,:)
681	sa_ad(:,:,:) = zsa_adin(:,:,:)
682
683	CALL tra_ldf_lap_adj( nit000 )
684
685	ztb_adout(:,:,:) = tb_ad(:,:,:)
686	zsb_adout(:,:,:) = sb_ad(:,:,:)
687	zta_adout(:,:,:) = ta_ad(:,:,:)
688	zsa_adout(:,:,:) = sa_ad(:,:,:)
689	zgtu_adout(:,:) = gtu_ad(:,:)
690	zgsu_adout(:,:) = gsu_ad(:,:)
691	zgtv_adout(:,:) = gtv_ad(:,:)
692	zgsv_adout(:,:) = gsv_ad(:,:)
693
694	zsp2_1 = DOT_PRODUCT( ztb_tlin , ztb_adout )
695	zsp2_2 = DOT_PRODUCT( zta_tlin , zta_adout )
696	zsp2_3 = DOT_PRODUCT( zgtu_tlin, zgtu_adout )
697	zsp2_4 = DOT_PRODUCT( zgtv_tlin, zgtv_adout )
698	zsp2_5 = DOT_PRODUCT( zsb_tlin , zsb_adout )
699	zsp2_6 = DOT_PRODUCT( zsa_tlin , zsa_adout )
700	zsp2_7 = DOT_PRODUCT( zgsu_tlin, zgsu_adout )
701	zsp2_8 = DOT_PRODUCT( zgsv_tlin, zgsv_adout )
702
703	zsp2_T = zsp2_1 + zsp2_2 + zsp2_3 + zsp2_4
704	zsp2_S = zsp2_5 + zsp2_6 + zsp2_7 + zsp2_8
705	zsp2 = zsp2_T + zsp2_S
706
707	cl_name = 'tra_ldf_lap_ad'
708	CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
709
710	DEALLOCATE( &
711	& ztb_tlin, & ! Tangent input
712	& zsb_tlin, & ! Tangent input
713	& zta_tlin, & ! Tangent input
714	& zsa_tlin, & ! Tangent input
715	& zgtu_tlin, & ! Tangent input
716	& zgsu_tlin, & ! Tangent input
717	& zgtv_tlin, & ! Tangent input
718	& zgsv_tlin, & ! Tangent input
719	& zta_tlout, & ! Tangent output
720	& zsa_tlout, & ! Tangent output
721	& zta_adin, & ! Adjoint input
722	& zsa_adin, & ! Adjoint input
723	& ztb_adout, & ! Adjoint output
724	& zsb_adout, & ! Adjoint output
725	& zta_adout, & ! Adjoint output
726	& zsa_adout, & ! Adjoint output
727	& zgtu_adout, & ! Adjoint output
728	& zgsu_adout, & ! Adjoint output
729	& zgtv_adout, & ! Adjoint output
730	& zgsv_adout, & ! Adjoint output
731	& z3r, & ! 3D random field
732	& z2r &
733	& )
734
735
736	END SUBROUTINE tra_ldf_lap_adj_tst
737
738	SUBROUTINE tra_ldf_lap_tlm_tst ( kumadt )
739	!!-----------------------------------------------------------------------
740	!!
741	!! * ROUTINE tra_ldf_lap_tlm_tst *
742	!!
743	!! ** Purpose : Test the tangent routine.
744	!!
745	!! ** Method : Verify the tangent with Taylor expansion
746	!!
747	!! M(x+hdx) = M(x) + L(hdx) + O(h^2)
748	!!
749	!! where L = Linear routine
750	!! M = direct routine
751	!! dx = input perturbation (random field)
752	!! h = ration on perturbation
753	!!
754	!! In the tangent test we verify that:
755	!! M(x+h*dx) - M(x)
756	!! g(h) = ------------------ ---> 1 as h ---> 0
757	!! L(h*dx)
758	!! and
759	!! g(h) - 1
760	!! f(h) = ---------- ---> k (costant) as h ---> 0
761	!! p
762	!!
763	!!
764	!! History :
765	!! ! 09-07 (A. Vigilant)
766	!!-----------------------------------------------------------------------
767	!! * Modules used
768	USE step_c1d ! Time stepping loop for the 1D configuration
769	!! USE asminc ! assimilation increments (asm_inc_init routine)
770	USE tamtrj ! writing out state trajectory
771	USE lib_mpp, ONLY: & ! distributed memory computing
772	& lk_mpp, &
773	& mpp_max
774	USE c1d, ONLY: & ! 1D configuration
775	& lk_c1d
776	USE par_tlm, ONLY: &
777	& cur_loop, &
778	& h_ratio
779	USE istate_mod
780	USE wzvmod ! vertical velocity
781	USE gridrandom, ONLY: &
782	& grid_rd_sd
783	USE trj_tam
784	USE oce , ONLY: & ! ocean dynamics and tracers variables
785	& tb, sb, tn, sn, ta, &
786	& sa, gtu, gsu, gtv, &
787	& gsv
788	USE traldf_lap ! lateral mixing (tra_ldf routine)
789	USE opatam_tst_ini, ONLY: &
790	& tlm_namrd
791	USE tamctl, ONLY: & ! Control parameters
792	& numtan, numtan_sc
793	!! * Arguments
794	INTEGER, INTENT(IN) :: &
795	& kumadt ! Output unit
796
797	!! * Local declarations
798	INTEGER :: &
799	& ji, & ! dummy loop indices
800	& jj, &
801	& jk, &
802	& istp
803	INTEGER, DIMENSION(jpi,jpj) :: &
804	& iseed_2d ! 2D seed for the random number generator
805	REAL(KIND=wp) :: &
806	& zsp1, & ! scalar product involving the tangent routine
807	& zsp1_Tb, &
808	& zsp1_Sb, &
809	& zsp1_Ta, &
810	& zsp1_Sa, &
811	& zsp1_Gtu, &
812	& zsp1_Gsu, &
813	& zsp1_Gtv, &
814	& zsp1_Gsv, &
815	& zsp2, & ! scalar product involving the tangent routine
816	& zsp2_Tb, &
817	& zsp2_Sb, &
818	& zsp2_Ta, &
819	& zsp2_Sa, &
820	& zsp2_Gtu, &
821	& zsp2_Gsu, &
822	& zsp2_Gtv, &
823	& zsp2_Gsv, &
824	& zsp3, & ! scalar product involving the tangent routine
825	& zsp3_Tb, &
826	& zsp3_Sb, &
827	& zsp3_Ta, &
828	& zsp3_Sa, &
829	& zsp3_Gtu, &
830	& zsp3_Gsu, &
831	& zsp3_Gtv, &
832	& zsp3_Gsv, &
833	& zzsp, & ! scalar product involving the tangent routine
834	& zzsp_Tb, &
835	& zzsp_Sb, &
836	& zzsp_Ta, &
837	& zzsp_Sa, &
838	& zzsp_Gtu, &
839	& zzsp_Gsu, &
840	& zzsp_Gtv, &
841	& zzsp_Gsv, &
842	& gamma, &
843	& zgsp1, &
844	& zgsp2, &
845	& zgsp3, &
846	& zgsp4, &
847	& zgsp5, &
848	& zgsp6, &
849	& zgsp7
850	REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
851	& ztb_tlin , & ! Tangent input
852	& zsb_tlin , & ! Tangent input
853	& zta_tlin , & ! Tangent input
854	& zsa_tlin , & ! Tangent input
855	& ztb_out , & ! Direct output
856	& zsb_out , & ! Direct output
857	& zta_out , & ! Direct output
858	& zsa_out , & ! Direct output
859	& ztb_wop , & !
860	& zsb_wop , & !
861	& zta_wop , & !
862	& zsa_wop , & !
863	& z3r ! 3D random field
864	REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
865	& zgtu_tlin , & ! Tangent input
866	& zgsu_tlin , & ! Tangent input
867	& zgtv_tlin , & ! Tangent input
868	& zgsv_tlin , & ! Tangent input
869	& zgtu_out , & ! Direct output
870	& zgsu_out , & ! Direct output
871	& zgtv_out , & ! Direct output
872	& zgsv_out , & ! Direct output
873	& zgtu_wop , & !
874	& zgsu_wop , & !
875	& zgtv_wop , & !
876	& zgsv_wop , & !
877	& z2r ! 2D random field
878	CHARACTER(LEN=14) :: cl_name
879	CHARACTER (LEN=128) :: file_out_sc, file_wop, file_out
880	CHARACTER (LEN=90) :: &
881	& FMT
882	REAL(KIND=wp), DIMENSION(100):: &
883	& zsctb, zscta, &
884	& zscsb, zscsa, &
885	& zscgtu, zscgsu, &
886	& zscgtv, zscgsv, &
887	& zscerrtb, &
888	& zscerrta, &
889	& zscerrsb, &
890	& zscerrsa, &
891	& zscerrgtu, &
892	& zscerrgsu, &
893	& zscerrgtv, &
894	& zscerrgsv
895	INTEGER, DIMENSION(100):: &
896	& iipostb, iiposta, iipossb, iipossa, &
897	& iiposgtu, iiposgsu, iiposgtv, iiposgsv, &
898	& ijpostb, ijposta, ijpossb, ijpossa, &
899	& ijposgtu, ijposgsu, ijposgtv, ijposgsv, &
900	& ikpostb, ikposta, ikpossb, ikpossa
901	INTEGER:: &
902	& ii, &
903	& isamp=40, &
904	& jsamp=40, &
905	& ksamp=10, &
906	& numsctlm, &
907	& numtlm
908	REAL(KIND=wp), DIMENSION(jpi,jpj,jpk) :: &
909	& zerrtb, zerrta, &
910	& zerrsb, zerrsa
911
912	REAL(KIND=wp), DIMENSION(jpi,jpj) :: &
913	& zerrgtu, zerrgsu, &
914	& zerrgtv, zerrgsv
915
916	! Allocate memory
917	ALLOCATE( &
918	& ztb_tlin(jpi,jpj,jpk), &
919	& zsb_tlin(jpi,jpj,jpk), &
920	& zta_tlin(jpi,jpj,jpk), &
921	& zsa_tlin(jpi,jpj,jpk), &
922	& ztb_out(jpi,jpj,jpk), &
923	& zsb_out(jpi,jpj,jpk), &
924	& zta_out(jpi,jpj,jpk), &
925	& zsa_out(jpi,jpj,jpk), &
926	& ztb_wop(jpi,jpj,jpk), &
927	& zsb_wop(jpi,jpj,jpk), &
928	& zta_wop(jpi,jpj,jpk), &
929	& zsa_wop(jpi,jpj,jpk), &
930	& zgtu_tlin(jpi,jpj), &
931	& zgsu_tlin(jpi,jpj), &
932	& zgtv_tlin(jpi,jpj), &
933	& zgsv_tlin(jpi,jpj), &
934	& zgtu_out(jpi,jpj), &
935	& zgsu_out(jpi,jpj), &
936	& zgtv_out(jpi,jpj), &
937	& zgsv_out(jpi,jpj), &
938	& zgtu_wop(jpi,jpj), &
939	& zgsu_wop(jpi,jpj), &
940	& zgtv_wop(jpi,jpj), &
941	& zgsv_wop(jpi,jpj), &
942	& z3r(jpi,jpj,jpk), &
943	& z2r(jpi,jpj) &
944	& )
945
946	!--------------------------------------------------------------------
947	! Reset variables
948	!--------------------------------------------------------------------
949
950	ztb_tlin(:,:,:) = 0.0_wp
951	zsb_tlin(:,:,:) = 0.0_wp
952	zta_tlin(:,:,:) = 0.0_wp
953	zsa_tlin(:,:,:) = 0.0_wp
954	ztb_out(:,:,:) = 0.0_wp
955	zsb_out(:,:,:) = 0.0_wp
956	zta_out(:,:,:) = 0.0_wp
957	zsa_out(:,:,:) = 0.0_wp
958	ztb_wop(:,:,:) = 0.0_wp
959	zsb_wop(:,:,:) = 0.0_wp
960	zta_wop(:,:,:) = 0.0_wp
961	zsa_wop(:,:,:) = 0.0_wp
962	zgtu_tlin(:,:) = 0.0_wp
963	zgsu_tlin(:,:) = 0.0_wp
964	zgtv_tlin(:,:) = 0.0_wp
965	zgsv_tlin(:,:) = 0.0_wp
966	zgtu_out(:,:) = 0.0_wp
967	zgsu_out(:,:) = 0.0_wp
968	zgtv_out(:,:) = 0.0_wp
969	zgsv_out(:,:) = 0.0_wp
970	zgtu_wop(:,:) = 0.0_wp
971	zgsu_wop(:,:) = 0.0_wp
972	zgtv_wop(:,:) = 0.0_wp
973	zgsv_wop(:,:) = 0.0_wp
974
975	tb_tl(:,:,:) = 0.0_wp
976	sb_tl(:,:,:) = 0.0_wp
977	ta_tl(:,:,:) = 0.0_wp
978	sa_tl(:,:,:) = 0.0_wp
979	gtu_tl(:,:) = 0.0_wp
980	gsu_tl(:,:) = 0.0_wp
981	gtv_tl(:,:) = 0.0_wp
982	gsv_tl(:,:) = 0.0_wp
983
984	zsctb(:) = 0.0_wp
985	zscta(:) = 0.0_wp
986	zscsb(:) = 0.0_wp
987	zscsa(:) = 0.0_wp
988	zscgtu(:) = 0.0_wp
989	zscgsu(:) = 0.0_wp
990	zscgtv(:) = 0.0_wp
991	zscgsv(:) = 0.0_wp
992	zscerrtb(:) = 0.0_wp
993	zscerrta(:) = 0.0_wp
994	zscerrsb(:) = 0.0_wp
995	zscerrsa(:) = 0.0_wp
996	zscerrgtu(:) = 0.0_wp
997	zscerrgsu(:) = 0.0_wp
998	zscerrgtv(:) = 0.0_wp
999	zscerrgsv(:) = 0.0_wp
1000
1001	!--------------------------------------------------------------------
1002	! Output filename Xn=F(X0)
1003	!--------------------------------------------------------------------
1004	file_wop='trj_wop_tldf_lap'
1005	CALL tlm_namrd
1006	gamma = h_ratio
1007	!--------------------------------------------------------------------
1008	! Initialize the tangent input with random noise: dx
1009	!--------------------------------------------------------------------
1010	IF (( cur_loop .NE. 0) .OR. ( gamma .NE. 0.0_wp) )THEN
1011	CALL grid_rd_sd( 596035, z3r, 'T', 0.0_wp, stdt)
1012	DO jk = 1, jpk
1013	DO jj = nldj, nlej
1014	DO ji = nldi, nlei
1015	ztb_tlin(ji,jj,jk) = z3r(ji,jj,jk)
1016	END DO
1017	END DO
1018	END DO
1019	CALL grid_rd_sd( 727391, z3r, 'S', 0.0_wp, stds)
1020	DO jk = 1, jpk
1021	DO jj = nldj, nlej
1022	DO ji = nldi, nlei
1023	zsb_tlin(ji,jj,jk) = z3r(ji,jj,jk)
1024	END DO
1025	END DO
1026	END DO
1027	CALL grid_rd_sd( 249741, z3r, 'T', 0.0_wp, stdt)
1028	DO jk = 1, jpk
1029	DO jj = nldj, nlej
1030	DO ji = nldi, nlei
1031	zta_tlin(ji,jj,jk) = z3r(ji,jj,jk)
1032	END DO
1033	END DO
1034	END DO
1035	CALL grid_rd_sd( 182029, z3r, 'S', 0.0_wp, stds)
1036	DO jk = 1, jpk
1037	DO jj = nldj, nlej
1038	DO ji = nldi, nlei
1039	zsa_tlin(ji,jj,jk) = z3r(ji,jj,jk)
1040	END DO
1041	END DO
1042	END DO
1043	CALL grid_rd_sd( 596035, z2r, 'U', 0.0_wp, stds)
1044	DO jj = nldj, nlej
1045	DO ji = nldi, nlei
1046	zgtu_tlin(ji,jj) = z2r(ji,jj)
1047	END DO
1048	END DO
1049	CALL grid_rd_sd( 727391, z2r, 'U', 0.0_wp, stds)
1050	DO jj = nldj, nlej
1051	DO ji = nldi, nlei
1052	zgsu_tlin(ji,jj) = z2r(ji,jj)
1053	END DO
1054	END DO
1055	CALL grid_rd_sd( 249741, z2r, 'V', 0.0_wp, stds)
1056	DO jj = nldj, nlej
1057	DO ji = nldi, nlei
1058	zgtv_tlin(ji,jj) = z2r(ji,jj)
1059	END DO
1060	END DO
1061	CALL grid_rd_sd( 182029, z2r, 'V', 0.0_wp, stds)
1062	DO jj = nldj, nlej
1063	DO ji = nldi, nlei
1064	zgsv_tlin(ji,jj) = z2r(ji,jj)
1065	END DO
1066	END DO
1067	ENDIF
1068	!--------------------------------------------------------------------
1069	! Complete Init for Direct
1070	!-------------------------------------------------------------------
1071	CALL istate_p
1072
1073	! *** initialize the reference trajectory
1074	! ------------
1075	CALL trj_rea( nit000-1, 1 )
1076	CALL trj_rea( nit000, 1 )
1077
1078	IF (( cur_loop .NE. 0) .OR. ( gamma .NE. 0.0_wp) )THEN
1079	ztb_tlin(:,:,:) = gamma * ztb_tlin(:,:,:)
1080	tb(:,:,:) = tb(:,:,:) + ztb_tlin(:,:,:)
1081
1082	zsb_tlin(:,:,:) = gamma * zsb_tlin(:,:,:)
1083	sb(:,:,:) = sb(:,:,:) + zsb_tlin(:,:,:)
1084
1085	zta_tlin(:,:,:) = gamma * zta_tlin(:,:,:)
1086	ta(:,:,:) = ta(:,:,:) + zta_tlin(:,:,:)
1087
1088	zsa_tlin(:,:,:) = gamma * zsa_tlin(:,:,:)
1089	sa(:,:,:) = sa(:,:,:) + zsa_tlin(:,:,:)
1090
1091	zgtu_tlin(:,:) = gamma * zgtu_tlin(:,:)
1092	gtu(:,:) = gtu(:,:) + zgtu_tlin(:,:)
1093
1094	zgsu_tlin(:,:) = gamma * zgsu_tlin(:,:)
1095	gsu(:,:) = gsu(:,:) + zgsu_tlin(:,:)
1096
1097	zgtv_tlin(:,:) = gamma * zgtv_tlin(:,:)
1098	gtv(:,:) = gtv(:,:) + zgtv_tlin(:,:)
1099
1100	zgsv_tlin(:,:) = gamma * zgsv_tlin(:,:)
1101	gsv(:,:) = gsv(:,:) + zgsv_tlin(:,:)
1102	ENDIF
1103	IF( .NOT. lk_vvl ) CALL wzv(nit000)
1104	!--------------------------------------------------------------------
1105	! Compute the direct model F(X0,t=n) = Xn
1106	!--------------------------------------------------------------------
1107	CALL tra_ldf_lap( nit000 )
1108	IF ( cur_loop .EQ. 0) CALL trj_wri_spl(file_wop)
1109
1110	!--------------------------------------------------------------------
1111	! Compute the Tangent
1112	!--------------------------------------------------------------------
1113	IF ( cur_loop .NE. 0) THEN
1114	!--------------------------------------------------------------------
1115	! Storing data
1116	!--------------------------------------------------------------------
1117	ztb_out (:,:,:) = tb (:,:,:)
1118	zsb_out (:,:,:) = sb (:,:,:)
1119	zta_out (:,:,:) = ta (:,:,:)
1120	zsa_out (:,:,:) = sa (:,:,:)
1121	zgtu_out (:,: ) = gtu (:,: )
1122	zgsu_out (:,: ) = gsu (:,: )
1123	zgtv_out (:,: ) = gtv (:,: )
1124	zgsv_out (:,: ) = gsv (:,: )
1125
1126	!--------------------------------------------------------------------
1127	! Initialize the tangent variables: dy^* = W dy
1128	!--------------------------------------------------------------------
1129	CALL trj_rea( nit000-1, 1 )
1130	CALL trj_rea( nit000, 1 )
1131	tb_tl (:,:,:) = ztb_tlin (:,:,:)
1132	sb_tl (:,:,:) = zsb_tlin (:,:,:)
1133	ta_tl (:,:,:) = zta_tlin (:,:,:)
1134	sa_tl (:,:,:) = zsa_tlin (:,:,:)
1135	gtu_tl ( :,:) = zgtu_tlin ( :,:)
1136	gsu_tl ( :,:) = zgsu_tlin ( :,:)
1137	gtv_tl ( :,:) = zgtv_tlin ( :,:)
1138	gsv_tl ( :,:) = zgsv_tlin ( :,:)
1139
1140	!-----------------------------------------------------------------------
1141	! Initialization of the dynamics and tracer fields for the tangent
1142	!-----------------------------------------------------------------------
1143	CALL tra_ldf_lap_tan(nit000)
1144
1145	!--------------------------------------------------------------------
1146	! Compute the scalar product: ( L(t0,tn) gamma dx0 ) )
1147	!--------------------------------------------------------------------
1148
1149	zsp2_Tb = DOT_PRODUCT( tb_tl , tb_tl )
1150	zsp2_Sb = DOT_PRODUCT( sb_tl , sb_tl )
1151	zsp2_Ta = DOT_PRODUCT( ta_tl , ta_tl )
1152	zsp2_Sa = DOT_PRODUCT( sa_tl , sa_tl )
1153	zsp2_Gtu = DOT_PRODUCT( gtu_tl, gtu_tl )
1154	zsp2_Gsu = DOT_PRODUCT( gsu_tl, gsu_tl )
1155	zsp2_Gtv = DOT_PRODUCT( gtv_tl, gtv_tl )
1156	zsp2_Gsv = DOT_PRODUCT( gsv_tl, gsv_tl )
1157
1158	zsp2 = zsp2_Tb + zsp2_Sb + zsp2_Ta + zsp2_Sa + &
1159	& zsp2_Gtu + Zsp2_Gsu + zsp2_Gtv + zsp2_Gsv
1160
1161	!--------------------------------------------------------------------
1162	! Storing data
1163	!--------------------------------------------------------------------
1164	CALL trj_rd_spl(file_wop)
1165	ztb_wop (:,:,:) = tb (:,:,:)
1166	zsb_wop (:,:,:) = sb (:,:,:)
1167	zta_wop (:,:,:) = ta (:,:,:)
1168	zsa_wop (:,:,:) = sa (:,:,:)
1169	zgtu_wop (:,: ) = gtu (:,: )
1170	zgsu_wop (:,: ) = gsu (:,: )
1171	zgtv_wop (:,: ) = gtv (:,: )
1172	zgsv_wop (:,: ) = gsv (:,: )
1173
1174	!--------------------------------------------------------------------
1175	! Compute the Linearization Error
1176	! Nn = M( X0+gamma.dX0, t0,tn) - M(X0, t0,tn)
1177	! and
1178	! Compute the Linearization Error
1179	! En = Nn -TL(gamma.dX0, t0,tn)
1180	!--------------------------------------------------------------------
1181	! Warning: Here we re-use local variables z()_out and z()_wop
1182	ii=0
1183	DO jk = 1, jpk
1184	DO jj = 1, jpj
1185	DO ji = 1, jpi
1186	ztb_out (ji,jj,jk) = ztb_out (ji,jj,jk) - ztb_wop (ji,jj,jk)
1187	ztb_wop (ji,jj,jk) = ztb_out (ji,jj,jk) - tb_tl (ji,jj,jk)
1188	IF ( tb_tl(ji,jj,jk) .NE. 0.0_wp ) &
1189	& zerrtb(ji,jj,jk) = ztb_out(ji,jj,jk)/tb_tl(ji,jj,jk)
1190	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1191	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1192	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1193	ii = ii+1
1194	iipostb(ii) = ji
1195	ijpostb(ii) = jj
1196	ikpostb(ii) = jk
1197	IF ( INT(tmask(ji,jj,jk)) .NE. 0) THEN
1198	! zsctb (ii) = ztb_out(ji,jj,jk)
1199	zsctb (ii) = ztb_wop(ji,jj,jk)
1200	zscerrtb (ii) = ( zerrtb(ji,jj,jk) - 1.0_wp ) / gamma
1201	ENDIF
1202	ENDIF
1203	END DO
1204	END DO
1205	END DO
1206	ii=0
1207	DO jk = 1, jpk
1208	DO jj = 1, jpj
1209	DO ji = 1, jpi
1210	zsb_out (ji,jj,jk) = zsb_out (ji,jj,jk) - zsb_wop (ji,jj,jk)
1211	zsb_wop (ji,jj,jk) = zsb_out (ji,jj,jk) - sb_tl (ji,jj,jk)
1212	IF ( sb_tl(ji,jj,jk) .NE. 0.0_wp ) &
1213	& zerrsb(ji,jj,jk) = zsb_out(ji,jj,jk)/sb_tl(ji,jj,jk)
1214	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1215	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1216	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1217	ii = ii+1
1218	iipossb(ii) = ji
1219	ijpossb(ii) = jj
1220	ikpossb(ii) = jk
1221	IF ( INT(tmask(ji,jj,jk)) .NE. 0) THEN
1222	! zscsb (ii) = zsb_out(ji,jj,jk)
1223	zscsb (ii) = zsb_wop(ji,jj,jk)
1224	zscerrsb (ii) = ( zerrsb(ji,jj,jk) - 1.0_wp ) / gamma
1225	ENDIF
1226	ENDIF
1227	END DO
1228	END DO
1229	END DO
1230	ii=0
1231	DO jk = 1, jpk
1232	DO jj = 1, jpj
1233	DO ji = 1, jpi
1234	zta_out (ji,jj,jk) = zta_out (ji,jj,jk) - zta_wop (ji,jj,jk)
1235	zta_wop (ji,jj,jk) = zta_out (ji,jj,jk) - ta_tl (ji,jj,jk)
1236	IF ( ta_tl(ji,jj,jk) .NE. 0.0_wp ) &
1237	& zerrta(ji,jj,jk) = zta_out(ji,jj,jk)/ta_tl(ji,jj,jk)
1238	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1239	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1240	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1241	ii = ii+1
1242	iiposta(ii) = ji
1243	ijposta(ii) = jj
1244	ikposta(ii) = jk
1245	IF ( INT(tmask(ji,jj,jk)) .NE. 0) THEN
1246	! zscta (ii) = zta_out(ji,jj,jk)
1247	zscta (ii) = zta_wop(ji,jj,jk)
1248	zscerrta (ii) = ( zerrta(ji,jj,jk) - 1.0_wp ) / gamma
1249	ENDIF
1250	ENDIF
1251	END DO
1252	END DO
1253	END DO
1254	ii=0
1255	DO jk = 1, jpk
1256	DO jj = 1, jpj
1257	DO ji = 1, jpi
1258	zsa_out (ji,jj,jk) = zsa_out (ji,jj,jk) - zsa_wop (ji,jj,jk)
1259	zsa_wop (ji,jj,jk) = zsa_out (ji,jj,jk) - sa_tl (ji,jj,jk)
1260	IF ( sa_tl(ji,jj,jk) .NE. 0.0_wp ) &
1261	& zerrsa(ji,jj,jk) = zsa_out(ji,jj,jk)/sa_tl(ji,jj,jk)
1262	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1263	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1264	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1265	ii = ii+1
1266	iipossa(ii) = ji
1267	ijpossa(ii) = jj
1268	ikpossa(ii) = jk
1269	IF ( INT(tmask(ji,jj,jk)) .NE. 0) THEN
1270	! zscsa (ii) = zsa_out(ji,jj,jk)
1271	zscsa (ii) = zsa_wop(ji,jj,jk)
1272	zscerrsa (ii) = ( zerrsa(ji,jj,jk) - 1.0_wp ) / gamma
1273	ENDIF
1274	ENDIF
1275	END DO
1276	END DO
1277	END DO
1278	ii=0
1279	jk=1
1280	DO jj = 1, jpj
1281	DO ji = 1, jpi
1282	zgtu_out (ji,jj)= zgtu_out (ji,jj) - zgtu_wop (ji,jj)
1283	zgtu_wop (ji,jj) = zgtu_out (ji,jj) - gtu_tl (ji,jj)
1284	IF ( gtu_tl(ji,jj) .NE. 0.0_wp ) &
1285	& zerrgtu(ji,jj) = zgtu_out(ji,jj)/gtu_tl(ji,jj)
1286	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1287	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1288	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1289	ii = ii+1
1290	iiposgtu(ii) = ji
1291	ijposgtu(ii) = jj
1292	IF ( INT(umask(ji,jj,jk)) .NE. 0) THEN
1293	! zscgtu (ii) = zgtu_out(ji,jj)
1294	zscgtu (ii) = zgtu_wop(ji,jj)
1295	zscerrgtu (ii) = ( zerrgtu(ji,jj) - 1.0_wp ) / gamma
1296	ENDIF
1297	ENDIF
1298	END DO
1299	END DO
1300	ii=0
1301	jk=1
1302	DO jj = 1, jpj
1303	DO ji = 1, jpi
1304	zgtv_out (ji,jj) = zgtv_out (ji,jj) - zgtv_wop (ji,jj)
1305	zgtv_wop (ji,jj) = zgtv_out (ji,jj) - gtv_tl (ji,jj)
1306	IF ( gtv_tl(ji,jj) .NE. 0.0_wp ) &
1307	& zerrgtv(ji,jj) = zgtv_out(ji,jj)/gtv_tl(ji,jj)
1308	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1309	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1310	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1311	ii = ii+1
1312	iiposgtv(ii) = ji
1313	ijposgtv(ii) = jj
1314	IF ( INT(vmask(ji,jj,jk)) .NE. 0) THEN
1315	! zscgtv (ii) = zgtv_out(ji,jj)
1316	zscgtv (ii) = zgtv_wop(ji,jj)
1317	zscerrgtv (ii) = ( zerrgtv(ji,jj) - 1.0_wp ) /gamma
1318	ENDIF
1319	ENDIF
1320	END DO
1321	END DO
1322	ii=0
1323	jk=1
1324	DO jj = 1, jpj
1325	DO ji = 1, jpi
1326	zgsu_out (ji,jj) = zgsu_out (ji,jj) - zgsu_wop (ji,jj)
1327	zgsu_wop (ji,jj) = zgsu_out (ji,jj) - gsu_tl (ji,jj)
1328	IF ( gsu_tl(ji,jj) .NE. 0.0_wp ) &
1329	& zerrgsu(ji,jj) = zgsu_out(ji,jj)/gsu_tl(ji,jj)
1330	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1331	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1332	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1333	ii = ii+1
1334	iiposgsu(ii) = ji
1335	ijposgsu(ii) = jj
1336	IF ( INT(umask(ji,jj,jk)) .NE. 0) THEN
1337	! zscgsu (ii) = zgsu_out(ji,jj)
1338	zscgsu (ii) = zgsu_wop(ji,jj)
1339	zscerrgsu (ii) = ( zerrgsu(ji,jj) - 1.0_wp ) / gamma
1340	ENDIF
1341	ENDIF
1342	END DO
1343	END DO
1344	ii=0
1345	jk=1
1346	DO jj = 1, jpj
1347	DO ji = 1, jpi
1348	zgsv_out (ji,jj) = zgsv_out (ji,jj) - zgsv_wop (ji,jj)
1349	zgsv_wop (ji,jj) = zgsv_out (ji,jj) - gsv_tl (ji,jj)
1350	IF ( gsv_tl(ji,jj) .NE. 0.0_wp ) &
1351	& zerrgsv(ji,jj) = zgsv_out(ji,jj)/gsv_tl(ji,jj)
1352	IF( (MOD(ji, isamp) .EQ. 0) .AND. &
1353	& (MOD(jj, jsamp) .EQ. 0) .AND. &
1354	& (MOD(jk, ksamp) .EQ. 0) ) THEN
1355	ii = ii+1
1356	iiposgsv(ii) = ji
1357	ijposgsv(ii) = jj
1358	IF ( INT(vmask(ji,jj,jk)) .NE. 0) THEN
1359	! zscgsv (ii) = zgsv_out(ji,jj)
1360	zscgsv (ii) = zgsv_wop(ji,jj)
1361	zscerrgsv (ii) = ( zerrgsv(ji,jj) - 1.0_wp ) /gamma
1362	ENDIF
1363	ENDIF
1364	END DO
1365	END DO
1366
1367	zsp1_Tb = DOT_PRODUCT( ztb_out , ztb_out )
1368	zsp1_Sb = DOT_PRODUCT( zsb_out , zsb_out )
1369	zsp1_Ta = DOT_PRODUCT( zta_out , zta_out )
1370	zsp1_Sa = DOT_PRODUCT( zsa_out , zsa_out )
1371	zsp1_Gtu = DOT_PRODUCT( zgtu_out , zgtu_out )
1372	zsp1_Gsu = DOT_PRODUCT( zgsu_out , zgsu_out )
1373	zsp1_Gtv = DOT_PRODUCT( zgtv_out , zgtv_out )
1374	zsp1_Gsv = DOT_PRODUCT( zgsv_out , zgsv_out )
1375
1376	zsp1 = zsp1_Tb + zsp1_Sb + zsp1_Ta + zsp1_Sa + &
1377	& zsp1_Gtu + Zsp1_Gsu + zsp1_Gtv + zsp1_Gsv
1378
1379	zsp3_Tb = DOT_PRODUCT( ztb_wop , ztb_wop )
1380	zsp3_Sb = DOT_PRODUCT( zsb_wop , zsb_wop )
1381	zsp3_Ta = DOT_PRODUCT( zta_wop , zta_wop )
1382	zsp3_Sa = DOT_PRODUCT( zsa_wop , zsa_wop )
1383	zsp3_Gtu = DOT_PRODUCT( zgtu_wop , zgtu_wop )
1384	zsp3_Gsu = DOT_PRODUCT( zgsu_wop , zgsu_wop )
1385	zsp3_Gtv = DOT_PRODUCT( zgtv_wop , zgtv_wop )
1386	zsp3_Gsv = DOT_PRODUCT( zgsv_wop , zgsv_wop )
1387
1388	zsp3 = zsp3_Tb + zsp3_Sb + zsp3_Ta + zsp3_Sa + &
1389	& zsp3_Gtu + Zsp3_Gsu + zsp3_Gtv + zsp3_Gsv
1390
1391	!--------------------------------------------------------------------
1392	! Print the linearization error En - norme 2
1393	!--------------------------------------------------------------------
1394	! 14 char:'12345678901234'
1395	cl_name = 'traldf_tam:En '
1396	zzsp = SQRT(zsp3)
1397	zzsp_Tb = SQRT(zsp3_Tb)
1398	zzsp_Sb = SQRT(zsp3_Sb)
1399	zzsp_Ta = SQRT(zsp3_Ta)
1400	zzsp_Sa = SQRT(zsp3_Sa)
1401	zzsp_Gtu = SQRT(zsp3_Gtu)
1402	zzsp_Gsu = SQRT(zsp3_Gsu)
1403	zzsp_Gtv = SQRT(zsp3_Gtv)
1404	zzsp_Gsv = SQRT(zsp3_Gsv)
1405
1406	CALL prntst_tlm( cl_name, kumadt, zzsp, h_ratio )
1407
1408	!--------------------------------------------------------------------
1409	! Compute the relative error Er of the linear model
1410	! Er = norm(En) / norm(TL(gamma.dX0, t0,tn))
1411	!--------------------------------------------------------------------
1412	cl_name = 'traldf:Er=En/L'
1413	zzsp = SQRT( zsp3/zsp2 )
1414	zzsp_Tb = SQRT( zsp3_Tb/zsp2_Tb)
1415	zzsp_Sb = SQRT( zsp3_Sb/zsp2_Sb)
1416	zzsp_Ta = SQRT(zsp3_Ta/zsp2_Ta)
1417	zzsp_Sa = SQRT(zsp3_Sa/zsp2_Sa)
1418	zzsp_Gtu = SQRT(zsp3_Gtu/zsp2_Gtu)
1419	zzsp_Gsu = SQRT(zsp3_Gsu/zsp2_Gsu)
1420	zzsp_Gtv = SQRT(zsp3_Gtv/zsp2_Gtv)
1421	zzsp_Gsv = SQRT(zsp3_Gsv/zsp2_Gsv)
1422
1423	zgsp3 = zzsp
1424	zgsp7 = zzsp/gamma
1425	zgsp4 = SQRT( zsp2 )
1426	zgsp5 = SQRT( zsp3 )
1427	CALL prntst_tlm( cl_name, kumadt, zzsp, h_ratio )
1428
1429	!--------------------------------------------------------------------
1430	! Compute TLM norm2
1431	!--------------------------------------------------------------------
1432	zzsp = SQRT(zsp2)
1433	zzsp_Tb = SQRT(zsp2_Tb)
1434	zzsp_Sb = SQRT(zsp2_Sb)
1435	zzsp_Ta = SQRT(zsp2_Ta)
1436	zzsp_Sa = SQRT(zsp2_Sa)
1437	zzsp_Gtu = SQRT(zsp2_Gtu)
1438	zzsp_Gsu = SQRT(zsp2_Gsu)
1439	zzsp_Gtv = SQRT(zsp2_Gtv)
1440	zzsp_Gsv = SQRT(zsp2_Gsv)
1441
1442	cl_name = 'traldf:L_n2 '
1443	CALL prntst_tlm( cl_name, kumadt, zzsp, h_ratio )
1444
1445	!--------------------------------------------------------------------
1446	! Print the linearization error Nn - norme 2
1447	!--------------------------------------------------------------------
1448	! 14 char:'12345678901234'
1449	cl_name = 'traldf:Mhdx-Mx'
1450	zzsp = SQRT(zsp1)
1451	zzsp_Tb = SQRT(zsp1_Tb)
1452	zzsp_Sb = SQRT(zsp1_Sa)
1453	zzsp_Ta = SQRT(zsp1_Ta)
1454	zzsp_Sa = SQRT(zsp1_Sa)
1455	zzsp_Gtu = SQRT(zsp1_Gtu)
1456	zzsp_Gsu = SQRT(zsp1_Gsu)
1457	zzsp_Gtv = SQRT(zsp1_Gtv)
1458	zzsp_Gsv = SQRT(zsp1_Gsv)
1459
1460	zgsp1 = zzsp
1461	zgsp2 = zgsp1 / SQRT(zsp2)
1462	zgsp6 = (zgsp2 - 1.0_wp)/gamma
1463	CALL prntst_tlm( cl_name, kumadt, zzsp, h_ratio )
1464
1465	FMT = "(A8,2X,I4.4,2X,E6.1,2X,E20.13,2X,E20.13,2X,E20.13,2X,E20.13,2X,E20.13,2X,E20.13,2X,E20.13)"
1466	WRITE(numtan,FMT) 'tldf_lap', cur_loop, h_ratio, zgsp1, zgsp2, zgsp3, zgsp4, zgsp5, zgsp6, zgsp7
1467
1468	!--------------------------------------------------------------------
1469	! Unitary calculus
1470	!--------------------------------------------------------------------
1471	FMT = "(A8,2X,A8,2X,I4.4,2X,E6.1,2X,I4.4,2X,I4.4,2X,I4.4,2X,E20.13,1X)"
1472	cl_name = 'tldf_lap'
1473	IF(lwp) THEN
1474	DO ii=1, 100, 1
1475	IF ( zsctb(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zsctb ', &
1476	& cur_loop, h_ratio, ii, iipostb(ii), ijpostb(ii), zsctb(ii)
1477	ENDDO
1478	DO ii=1, 100, 1
1479	IF ( zscsb(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscsb ', &
1480	& cur_loop, h_ratio, ii, iipossb(ii), ijpossb(ii), zscsb(ii)
1481	ENDDO
1482	DO ii=1, 100, 1
1483	IF ( zscta(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscta ', &
1484	& cur_loop, h_ratio, ii, iiposta(ii), ijposta(ii), zscta(ii)
1485	ENDDO
1486	DO ii=1, 100, 1
1487	IF ( zscsa(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscsa ', &
1488	& cur_loop, h_ratio, ii, iipossa(ii), ijpossa(ii), zscsa(ii)
1489	ENDDO
1490	DO ii=1, 100, 1
1491	IF ( zscerrtb(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscerrtb ', &
1492	& cur_loop, h_ratio, ii, iipostb(ii), ijpostb(ii), zscerrtb(ii)
1493	ENDDO
1494	DO ii=1, 100, 1
1495	IF ( zscerrsb(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscerrsb ', &
1496	& cur_loop, h_ratio, ii, iipossb(ii), ijpossb(ii), zscerrsb(ii)
1497	ENDDO
1498	DO ii=1, 100, 1
1499	IF ( zscerrta(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscerrta ', &
1500	& cur_loop, h_ratio, ii, iiposta(ii), ijposta(ii), zscerrta(ii)
1501	ENDDO
1502	DO ii=1, 100, 1
1503	IF ( zscerrsa(ii) .NE. 0.0_wp ) WRITE(numtan_sc,FMT) cl_name, 'zscerrsa ', &
1504	& cur_loop, h_ratio, ii, iipossa(ii), ijpossa(ii), zscerrsa(ii)
1505	ENDDO
1506	! write separator
1507	WRITE(numtan_sc,"(A4)") '===='
1508	ENDIF
1509
1510	ENDIF
1511
1512	DEALLOCATE( &
1513	& ztb_tlin, & ! Tangent input
1514	& zsb_tlin, & ! Tangent input
1515	& zta_tlin, & ! Tangent input
1516	& zsa_tlin, & ! Tangent input
1517	& ztb_out, &
1518	& zsb_out, &
1519	& zta_out, &
1520	& zsa_out, &
1521	& ztb_wop, &
1522	& zsb_wop, &
1523	& zta_wop, &
1524	& zsa_wop, &
1525	& zgtu_tlin, & ! Tangent input
1526	& zgsu_tlin, & ! Tangent input
1527	& zgtv_tlin, & ! Tangent input
1528	& zgsv_tlin, & ! Tangent input
1529	& zgtu_out, &
1530	& zgsu_out, &
1531	& zgtv_out, &
1532	& zgsv_out, &
1533	& zgtu_wop, &
1534	& zgsu_wop, &
1535	& zgtv_wop, &
1536	& zgsv_wop, &
1537	& z3r, & ! 3D random field
1538	& z2r &
1539	& )
1540
1541	END SUBROUTINE tra_ldf_lap_tlm_tst
1542
1543
1544
1545	SUBROUTINE asm_trj_wop_wri
1546	!!-----------------------------------------------------------------------
1547	!!
1548	!! * ROUTINE asm_trj__wop_wri *
1549	!!
1550	!! ** Purpose : Write to file the model state trajectory for use with
1551	!! 4D-Var.
1552	!!
1553	!! ** Method :
1554	!!
1555	!! ** Action :
1556	!!
1557	!! History :
1558	!! ! 09-07 (F. Vigilant)
1559	!!-----------------------------------------------------------------------
1560	!! *Module udes
1561	USE iom
1562	USE oce , ONLY: & ! ocean dynamics and tracers variables
1563	& tb, sb, ta, sa, gtu, gsu, gtv, gsv
1564	!! * Arguments
1565	!! * Local declarations
1566	INTEGER :: &
1567	& inum ! File unit number
1568	CHARACTER (LEN=50) :: &
1569	& filename
1570	! Define the output file
1571	filename='trj_wop'
1572	WRITE(filename, FMT='(A,A)' ) TRIM( filename ), '.nc'
1573	filename = TRIM( filename )
1574	CALL iom_open( filename, inum, ldwrt = .TRUE., kiolib = jprstlib)
1575
1576	! Output trajectory fields
1577	CALL iom_rstput( 1, 1, inum, 'tb' , tb )
1578	CALL iom_rstput( 1, 1, inum, 'sb' , sb )
1579	CALL iom_rstput( 1, 1, inum, 'ta' , ta )
1580	CALL iom_rstput( 1, 1, inum, 'sa' , sa )
1581	CALL iom_rstput( 1, 1, inum, 'gtu' , gtu )
1582	CALL iom_rstput( 1, 1, inum, 'gsu' , gsu )
1583	CALL iom_rstput( 1, 1, inum, 'gtv' , gtv )
1584	CALL iom_rstput( 1, 1, inum, 'gsv' , gsv )
1585	CALL iom_close( inum )
1586	END SUBROUTINE asm_trj_wop_wri
1587
1588	SUBROUTINE asm_trj_wop_rd
1589	!!-----------------------------------------------------------------------
1590	!!
1591	!! * ROUTINE asm_trj__wop_rd *
1592	!!
1593	!! ** Purpose : Read file the model state trajectory for use with
1594	!! 4D-Var.
1595	!!
1596	!! ** Method :
1597	!!
1598	!! ** Action :
1599	!!
1600	!! History :
1601	!! ! 09-07 (F. Vigilant)
1602	!!-----------------------------------------------------------------------
1603	!! *Module udes
1604	USE iom ! I/O module
1605	USE oce , ONLY: & ! ocean dynamics and tracers variables
1606	& tb, sb, ta, sa, gtu, gsu, gtv, gsv
1607	!! * Arguments
1608	!! * Local declarations
1609	INTEGER :: &
1610	& inum ! File unit number
1611	CHARACTER (LEN=50) :: &
1612	& filename
1613	! Define the output file
1614	filename='trj_wop'
1615	WRITE(filename, FMT='(A,A)' ) TRIM( filename ), '.nc'
1616	filename = TRIM( filename )
1617	CALL iom_open( filename, inum)
1618
1619	! Output trajectory fields
1620	CALL iom_get( inum, jpdom_data, 'tb' , tb, 1 )
1621	CALL iom_get( inum, jpdom_data, 'sb' , sb, 1 )
1622	CALL iom_get( inum, jpdom_data, 'ta' , ta, 1 )
1623	CALL iom_get( inum, jpdom_data, 'sa' , sa, 1 )
1624	CALL iom_get( inum, jpdom_data, 'gtu' , gtu, 1 )
1625	CALL iom_get( inum, jpdom_data, 'gsu' , gsu, 1 )
1626	CALL iom_get( inum, jpdom_data, 'gtv' , gtv, 1 )
1627	CALL iom_get( inum, jpdom_data, 'gsv' , gsv, 1 )
1628	CALL iom_close( inum )
1629	END SUBROUTINE asm_trj_wop_rd
1630
1631	#endif
1632	!!==============================================================================
1633	END MODULE traldf_lap_tam

Note: See TracBrowser for help on using the repository browser.

Download in other formats: