Context Navigation

traldf_lap_tam.F90 @ 2587

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

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

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

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

Download in other formats: