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.

cla_tam.F90 in branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/LBC – NEMO

Context Navigation

source: branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/LBC/cla_tam.F90 @ 4391

Last change on this file since 4391 was 3668, checked in by pabouttier, 12 years ago
Missing allocation for sbc_ice variables in TAM - See Ticket #1026
File size: 115.5 KB

Line
1	MODULE cla_tam
2	#if defined key_tam
3	!!======================================================================
4	!! * MODULE cla *
5	!! Cross Land Advection : specific update of the horizontal divergence,
6	!! tracer trends and after velocity
7	!!
8	!! --- Specific to ORCA_R2 ---
9	!!
10	!!======================================================================
11	!! History : 1.0 ! 2002-11 (A. Bozec) Original code
12	!! 3.2 ! 2009-07 (G. Madec) merge cla, cla_div, tra_cla, cla_dynspg
13	!! ! and correct a mpp bug reported by A.R. Porter
14	!! History of TAM :
15	!! 3.4 ! 2012-07 (P.-A. Bouttier)
16	!!----------------------------------------------------------------------
17	#if defined key_orca_r2
18	!!----------------------------------------------------------------------
19	!! 'key_orca_r2' global ocean model R2
20	!!----------------------------------------------------------------------
21	!! cla_div : update of horizontal divergence at cla straits
22	!! tra_cla : update of tracers at cla straits
23	!! cla_dynspg : update of after horizontal velocities at cla straits
24	!! cla_init : initialisation - control check
25	!! cla_bab_el_mandeb : cross land advection for Bab-el-mandeb strait
26	!! cla_gibraltar : cross land advection for Gibraltar strait
27	!! cla_hormuz : cross land advection for Hormuz strait
28	!!----------------------------------------------------------------------
29	USE oce ! ocean dynamics and tracers
30	USE dom_oce ! ocean space and time domain
31	USE sbc_oce ! surface boundary condition: ocean
32	USE dynspg_oce ! ocean dynamics: surface pressure gradient variables
33	USE oce_tam ! ocean dynamics and tracers
34	USE sbc_oce_tam ! surface boundary condition: ocean
35	!USE dynspg_oce_tam ! ocean dynamics: surface pressure gradient variables
36	USE in_out_manager ! I/O manager
37	USE lib_mpp ! distributed memory computing library
38	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
39	USE dotprodfld
40	USE tstool_tam
41	USE paresp
42	USE gridrandom
43
44	IMPLICIT NONE
45	PRIVATE
46
47	PUBLIC cla_init_tam ! routine called by opatam.F90
48	PUBLIC cla_div_tan ! routine called by divcur_tan.F90
49	PUBLIC cla_traadv_tan ! routine called by traadv_tan.F90
50	PUBLIC cla_dynspg_tan ! routine called by dynspg_flt_tan.F90
51	PUBLIC cla_div_adj ! routine called by divcur_adj.F90
52	PUBLIC cla_traadv_adj ! routine called by traadv_adj.F90
53	PUBLIC cla_dynspg_adj ! routine called by dynspg_flt_adj.F90
54	PUBLIC cla_div_adj_tst ! routine called by tamtst.F90
55	PUBLIC cla_traadv_adj_tst ! routine called by traadv_adj.F90
56	PUBLIC cla_dynspg_adj_tst ! routine called by dynspg_flt_adj.F90
57
58	INTEGER :: nbab, ngib, nhor ! presence or not of required grid-point on local domain
59	! ! for Bab-el-Mandeb, Gibraltar, and Hormuz straits
60
61	! ! fixed part ! time evolving !!! profile of hdiv for some straits
62	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_101, hdiv_139_101_kt ! Gibraltar (i,j)=(172,101)
63	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_102 ! Gibraltar (i,j)=(139,102)
64	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_141_102, hdiv_141_102_kt ! Gibraltar (i,j)=(141,102)
65	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_88 , hdiv_161_88_kt ! Bab-el-Mandeb (i,j)=(161,88)
66	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_87 ! Bab-el-Mandeb (i,j)=(161,87)
67	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_160_89 , hdiv_160_89_kt ! Bab-el-Mandeb (i,j)=(160,89)
68	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_172_94 ! Hormuz (i,j)=(172, 94)
69
70	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: t_171_94_hor, s_171_94_hor ! Temperature, salinity in Hormuz strait
71
72	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_101_tl, hdiv_139_101_kt_tl ! Gibraltar (i,j)=(172,101)
73	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_102_tl ! Gibraltar (i,j)=(139,102)
74	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_141_102_tl, hdiv_141_102_kt_tl ! Gibraltar (i,j)=(141,102)
75	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_88_tl , hdiv_161_88_kt_tl ! Bab-el-Mandeb (i,j)=(161,88)
76	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_87_tl ! Bab-el-Mandeb (i,j)=(161,87)
77	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_160_89_tl , hdiv_160_89_kt_tl ! Bab-el-Mandeb (i,j)=(160,89)
78	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_172_94_tl ! Hormuz (i,j)=(172, 94)
79
80	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: t_171_94_hor_tl, s_171_94_hor_tl ! Temperature, salinity in Hormuz strait
81
82	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_101_ad, hdiv_139_101_kt_ad ! Gibraltar (i,j)=(172,101)
83	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_139_102_ad ! Gibraltar (i,j)=(139,102)
84	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_141_102_ad, hdiv_141_102_kt_ad ! Gibraltar (i,j)=(141,102)
85	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_88_ad , hdiv_161_88_kt_ad ! Bab-el-Mandeb (i,j)=(161,88)
86	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_161_87_ad ! Bab-el-Mandeb (i,j)=(161,87)
87	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_160_89_ad , hdiv_160_89_kt_ad ! Bab-el-Mandeb (i,j)=(160,89)
88	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: hdiv_172_94_ad ! Hormuz (i,j)=(172, 94)
89
90	REAL(wp), ALLOCATABLE, SAVE, DIMENSION (:) :: t_171_94_hor_ad, s_171_94_hor_ad ! Temperature, salinity in Hormuz strait
91	!! * Substitutions
92	# include "domzgr_substitute.h90"
93	!!----------------------------------------------------------------------
94	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
95	!! $Id$
96	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
97	!!----------------------------------------------------------------------
98	CONTAINS
99
100	SUBROUTINE cla_div_tan( kt )
101	!!----------------------------------------------------------------------
102	!! * ROUTINE div_cla_tan *
103	!!
104	!! ** Purpose : update the horizontal divergence of the velocity field
105	!! at some straits ( Gibraltar, Bab el Mandeb and Hormuz ).
106	!!
107	!! ** Method : - first time-step: initialisation of cla
108	!! - all time-step: using imposed transport at each strait,
109	!! the now horizontal divergence is updated
110	!!
111	!! ** Action : phdivn updted now horizontal divergence at cla straits
112	!!----------------------------------------------------------------------
113	INTEGER, INTENT( in ) :: kt ! ocean time-step index
114	!!----------------------------------------------------------------------
115	!
116	IF( kt == nit000 ) THEN
117	!
118	IF(lwp) WRITE(numout,*)
119	IF(lwp) WRITE(numout,*) 'div_cla_tan : cross land advection on hdiv '
120	IF(lwp) WRITE(numout,*) '~~~~~~~~'
121	!
122	IF( nbab == 1 ) CALL cla_bab_el_mandeb_tan('ini') ! Bab el Mandeb ( Red Sea - Indian ocean )
123	IF( ngib == 1 ) CALL cla_gibraltar_tan ('ini') ! Gibraltar strait (Med Sea - Atlantic ocean)
124	IF( nhor == 1 ) CALL cla_hormuz_tan ('ini') ! Hormuz Strait ( Persian Gulf - Indian ocean )
125	!
126	ENDIF
127	!
128	IF( nbab == 1 ) CALL cla_bab_el_mandeb_tan('div') ! Bab el Mandeb ( Red Sea - Indian ocean )
129	IF( ngib == 1 ) CALL cla_gibraltar_tan ('div') ! Gibraltar strait (Med Sea - Atlantic ocean)
130	IF( nhor == 1 ) CALL cla_hormuz_tan ('div') ! Hormuz Strait ( Persian Gulf - Indian ocean )
131	!
132	!!gm lbc useless here, no?
133	!!gm CALL lbc_lnk( hdivn, 'T', 1. ) ! Lateral boundary conditions on hdivn
134	!
135	END SUBROUTINE cla_div_tan
136
137	SUBROUTINE cla_div_adj( kt )
138	!!----------------------------------------------------------------------
139	!! * ROUTINE div_cla_adj *
140	!!
141	!! ** Purpose : update the horizontal divergence of the velocity field
142	!! at some straits ( Gibraltar, Bab el Mandeb and Hormuz ).
143	!!
144	!! ** Method : - first time-step: initialisation of cla
145	!! - all time-step: using imposed transport at each strait,
146	!! the now horizontal divergence is updated
147	!!
148	!! ** Action : phdivn updted now horizontal divergence at cla straits
149	!!----------------------------------------------------------------------
150	INTEGER, INTENT( in ) :: kt ! ocean time-step index
151	!!----------------------------------------------------------------------
152	!
153	IF( kt == nitend ) THEN
154	!
155	IF(lwp) WRITE(numout,*)
156	IF(lwp) WRITE(numout,*) 'div_cla_adj : cross land advection on hdiv '
157	IF(lwp) WRITE(numout,*) '~~~~~~~~'
158	!
159	IF( nbab == 1 ) CALL cla_bab_el_mandeb_adj('ini') ! Bab el Mandeb ( Red Sea - Indian ocean )
160	IF( ngib == 1 ) CALL cla_gibraltar_adj ('ini') ! Gibraltar strait (Med Sea - Atlantic ocean)
161	IF( nhor == 1 ) CALL cla_hormuz_adj ('ini') ! Hormuz Strait ( Persian Gulf - Indian ocean )
162	!
163	ENDIF
164	!
165	IF( nbab == 1 ) CALL cla_bab_el_mandeb_adj('div') ! Bab el Mandeb ( Red Sea - Indian ocean )
166	IF( ngib == 1 ) CALL cla_gibraltar_adj ('div') ! Gibraltar strait (Med Sea - Atlantic ocean)
167	IF( nhor == 1 ) CALL cla_hormuz_adj ('div') ! Hormuz Strait ( Persian Gulf - Indian ocean )
168	!
169	!!gm lbc useless here, no?
170	!!gm CALL lbc_lnk( hdivn, 'T', 1. ) ! Lateral boundary conditions on hdivn
171	!
172	END SUBROUTINE cla_div_adj
173
174	SUBROUTINE cla_traadv_tan( kt )
175	!!----------------------------------------------------------------------
176	!! * ROUTINE tra_cla_tan *
177	!!
178	!! ** Purpose : Update the now trend due to the advection of tracers
179	!! and add it to the general trend of passive tracer equations
180	!! at some straits ( Bab el Mandeb, Gibraltar, Hormuz ).
181	!!
182	!! ** Method : using both imposed transport at each strait and T & S
183	!! budget, the now tracer trends is updated
184	!!
185	!! ** Action : (ta,sa) updated now tracer trends at cla straits
186	!!----------------------------------------------------------------------
187	INTEGER, INTENT( in ) :: kt ! ocean time-step index
188	!!----------------------------------------------------------------------
189	!
190	IF( kt == nit000 ) THEN
191	IF(lwp) WRITE(numout,*)
192	IF(lwp) WRITE(numout,*) 'tra_cla_tan : cross land advection on tracers '
193	IF(lwp) WRITE(numout,*) '~~~~~~~~'
194	ENDIF
195	!
196	IF( nbab == 1 ) CALL cla_bab_el_mandeb_tan('tra') ! Bab el Mandeb strait
197	IF( ngib == 1 ) CALL cla_gibraltar_tan ('tra') ! Gibraltar strait
198	IF( nhor == 1 ) CALL cla_hormuz_tan ('tra') ! Hormuz Strait ( Persian Gulf)
199	!
200	END SUBROUTINE cla_traadv_tan
201
202	SUBROUTINE cla_traadv_adj( kt )
203	!!----------------------------------------------------------------------
204	!! * ROUTINE tra_cla_adj *
205	!!
206	!! ** Purpose : Update the now trend due to the advection of tracers
207	!! and add it to the general trend of passive tracer equations
208	!! at some straits ( Bab el Mandeb, Gibraltar, Hormuz ).
209	!!
210	!! ** Method : using both imposed transport at each strait and T & S
211	!! budget, the now tracer trends is updated
212	!!
213	!! ** Action : (ta,sa) updated now tracer trends at cla straits
214	!!----------------------------------------------------------------------
215	INTEGER, INTENT( in ) :: kt ! ocean time-step index
216	!!----------------------------------------------------------------------
217	!
218	IF( kt == nitend ) THEN
219	IF(lwp) WRITE(numout,*)
220	IF(lwp) WRITE(numout,*) 'tra_cla_adj : cross land advection on tracers '
221	IF(lwp) WRITE(numout,*) '~~~~~~~~'
222	ENDIF
223	!
224	IF( nbab == 1 ) CALL cla_bab_el_mandeb_adj('tra') ! Bab el Mandeb strait
225	IF( ngib == 1 ) CALL cla_gibraltar_adj ('tra') ! Gibraltar strait
226	IF( nhor == 1 ) CALL cla_hormuz_adj ('tra') ! Hormuz Strait ( Persian Gulf)
227	!
228	END SUBROUTINE cla_traadv_adj
229
230	SUBROUTINE cla_dynspg_tan( kt )
231	!!----------------------------------------------------------------------
232	!! * ROUTINE cla_dynspg_tan *
233	!!
234	!! ** Purpose : Update the after velocity at some straits
235	!! (Bab el Mandeb, Gibraltar, Hormuz).
236	!!
237	!! ** Method : required to compute the filtered surface pressure gradient
238	!!
239	!! ** Action : (ua,va) after velocity at the cla straits
240	!!----------------------------------------------------------------------
241	INTEGER, INTENT( in ) :: kt ! ocean time-step index
242	!!----------------------------------------------------------------------
243	!
244	IF( kt == nit000 ) THEN
245	IF(lwp) WRITE(numout,*)
246	IF(lwp) WRITE(numout,*) 'cla_dynspg_tan : cross land advection on (ua,va) '
247	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
248	ENDIF
249	!
250	IF( nbab == 1 ) CALL cla_bab_el_mandeb_tan('spg') ! Bab el Mandeb strait
251	IF( ngib == 1 ) CALL cla_gibraltar_tan ('spg') ! Gibraltar strait
252	IF( nhor == 1 ) CALL cla_hormuz_tan ('spg') ! Hormuz Strait ( Persian Gulf)
253	!
254	!!gm lbc is needed here, not?
255	!!gm CALL lbc_lnk( hdivn, 'U', -1. ) ; CALL lbc_lnk( hdivn, 'V', -1. ) ! Lateral boundary conditions
256	!
257	END SUBROUTINE cla_dynspg_tan
258
259	SUBROUTINE cla_dynspg_adj( kt )
260	!!----------------------------------------------------------------------
261	!! * ROUTINE cla_dynspg_adj *
262	!!
263	!! ** Purpose : Update the after velocity at some straits
264	!! (Bab el Mandeb, Gibraltar, Hormuz).
265	!!
266	!! ** Method : required to compute the filtered surface pressure gradient
267	!!
268	!! ** Action : (ua,va) after velocity at the cla straits
269	!!----------------------------------------------------------------------
270	INTEGER, INTENT( in ) :: kt ! ocean time-step index
271	!!----------------------------------------------------------------------
272	!
273	IF( kt == nitend ) THEN
274	IF(lwp) WRITE(numout,*)
275	IF(lwp) WRITE(numout,*) 'cla_dynspg_adj : cross land advection on (ua,va) '
276	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
277	ENDIF
278	!
279	IF( nbab == 1 ) CALL cla_bab_el_mandeb_adj('spg') ! Bab el Mandeb strait
280	IF( ngib == 1 ) CALL cla_gibraltar_adj ('spg') ! Gibraltar strait
281	IF( nhor == 1 ) CALL cla_hormuz_adj ('spg') ! Hormuz Strait ( Persian Gulf)
282	!
283	!!gm lbc is needed here, not?
284	!!gm CALL lbc_lnk( hdivn, 'U', -1. ) ; CALL lbc_lnk( hdivn, 'V', -1. ) ! Lateral boundary conditions
285	!
286	END SUBROUTINE cla_dynspg_adj
287
288	SUBROUTINE cla_init_tam
289	!! -------------------------------------------------------------------
290	!! * ROUTINE cla_init_tam *
291	!!
292	!! ** Purpose : control check for mpp computation
293	!!
294	!! ** Method : - All the strait grid-points must be inside one of the
295	!! local domain interior for the cla advection to work
296	!! properly in mpp (i.e. inside (2:jpim1,2:jpjm1) ).
297	!! Define the corresponding indicators (nbab, ngib, nhor)
298	!! - The profiles of cross-land fluxes are currently hard
299	!! coded for L31 levels. Stop if jpk/=31
300	!!
301	!! ** Action : nbab, ngib, nhor strait inside the local domain or not
302	!!---------------------------------------------------------------------
303	REAL(wp) :: ztemp ! local scalar
304	INTEGER :: ierr ! local integer
305	!!---------------------------------------------------------------------
306	!
307	IF(lwp) WRITE(numout,*)
308	IF(lwp) WRITE(numout,*) 'cla_init_tam : cross land advection initialisation '
309	IF(lwp) WRITE(numout,*) '~~~~~~~~~'
310	!
311	! ! Allocate arrays for this module
312	ALLOCATE( hdiv_139_101(jpk) , hdiv_139_101_kt(jpk) , & ! Gibraltar
313	& hdiv_139_102(jpk) , &
314	& hdiv_141_102(jpk) , hdiv_141_102_kt(jpk) , &
315	& hdiv_161_88 (jpk) , hdiv_161_88_kt (jpk) , & ! Bab-el-Mandeb
316	& hdiv_161_87 (jpk) , &
317	& hdiv_160_89 (jpk) , hdiv_160_89_kt (jpk) , & ! Hormuz
318	& hdiv_172_94 (jpk) , &
319	& t_171_94_hor(jpk) , s_171_94_hor (jpk) , STAT=ierr )
320	ALLOCATE( hdiv_139_101_tl(jpk) , hdiv_139_101_kt_tl(jpk) , & ! Gibraltar
321	& hdiv_139_102_tl(jpk) , &
322	& hdiv_141_102_tl(jpk) , hdiv_141_102_kt_tl(jpk) , &
323	& hdiv_161_88_tl (jpk) , hdiv_161_88_kt_tl (jpk) , & ! Bab-el-Mandeb
324	& hdiv_161_87_tl (jpk) , &
325	& hdiv_160_89_tl (jpk) , hdiv_160_89_kt_tl (jpk) , & ! Hormuz
326	& hdiv_172_94_tl (jpk) , &
327	& t_171_94_hor_tl(jpk) , s_171_94_hor_tl (jpk) , STAT=ierr )
328	ALLOCATE( hdiv_139_101_ad(jpk) , hdiv_139_101_kt_ad(jpk) , & ! Gibraltar
329	& hdiv_139_102_ad(jpk) , &
330	& hdiv_141_102_ad(jpk) , hdiv_141_102_kt_ad(jpk) , &
331	& hdiv_161_88_ad (jpk) , hdiv_161_88_kt_ad (jpk) , & ! Bab-el-Mandeb
332	& hdiv_161_87_ad (jpk) , &
333	& hdiv_160_89_ad (jpk) , hdiv_160_89_kt_ad (jpk) , & ! Hormuz
334	& hdiv_172_94_ad (jpk) , &
335	& t_171_94_hor_ad(jpk) , s_171_94_hor_ad (jpk) , STAT=ierr )
336	IF( lk_mpp ) CALL mpp_sum( ierr )
337	IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'cla_init_tam: unable to allocate arrays' )
338	!
339	IF( .NOT.lk_dynspg_flt ) CALL ctl_stop( 'cla_init_tam: Cross Land Advection works only with lk_dynspg_flt=T ' )
340	!
341	IF( lk_vvl ) CALL ctl_stop( 'cla_init_tam: Cross Land Advection does not work with lk_vvl=T option' )
342	!
343	IF( jpk /= 31 ) CALL ctl_stop( 'cla_init_tam: Cross Land Advection hard coded for ORCA_R2_L31' )
344	!
345	! _\|_______\|_______\|_
346	! 89 \| \|///////\|
347	! _\|_______\|_______\|_
348	! ------------------------ ! 88 \|///////\| \|
349	! Bab el Mandeb strait ! _\|_______\|_______\|_
350	! ------------------------ ! 87 \|///////\| \|
351	! _\|_______\|_______\|_
352	! \| 160 \| 161 \|
353	!
354	! The 6 Bab el Mandeb grid-points must be inside one of the interior of the
355	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
356	nbab = 0
357	IF( ( 1 <= mj0( 88) .AND. mj1( 89) <= jpj ) .AND. & !* (161,89), (161,88) and (161,88) on the local pocessor
358	& ( 1 <= mi0(160) .AND. mi1(161) <= jpi ) ) nbab = 1
359	!
360	! test if there is no local domain that includes all required grid-points
361	ztemp = REAL( nbab )
362	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
363	IF( ztemp == 0 ) THEN ! Only 2 points in each direction, this should never be a problem
364	CALL ctl_stop( ' cross land advection at Bab-el_Mandeb does not work with your processor cutting: change it' )
365	ENDIF
366	! ___________________________
367	! ------------------------ ! 102 \| \|///////\| \|
368	! Gibraltar strait ! _\|_______\|_______\|_______\|_
369	! ------------------------ ! 101 \| \|///////\| \|
370	! _\|_______\|_______\|_______\|_
371	! \| 139 \| 140 \| 141 \|
372	!
373	! The 6 Gibraltar grid-points must be inside one of the interior of the
374	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
375	ngib = 0
376	IF( ( 2 <= mj0(101) .AND. mj1(102) <= jpjm1 ) .AND. & !* (139:141,101:102) on the local pocessor
377	& ( 2 <= mi0(139) .AND. mi1(141) <= jpim1 ) ) ngib = 1
378	!
379	! test if there is no local domain that includes all required grid-points
380	ztemp = REAL( ngib )
381	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
382	IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting
383	CALL ctl_stop( ' cross land advection at Gibraltar does not work with your processor cutting: change it' )
384	ENDIF
385	! _______________
386	! ------------------------ ! 94 \|/////\| \|
387	! Hormuz strait ! _\|_____\|_____\|_
388	! ------------------------ ! 171 172
389	!
390	! The 2 Hormuz grid-points must be inside one of the interior of the
391	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
392	nhor = 0
393	IF( 2 <= mj0( 94) .AND. mj1( 94) <= jpjm1 .AND. &
394	& 2 <= mi0(171) .AND. mi1(172) <= jpim1 ) nhor = 1
395	!
396	! test if there is no local domain that includes all required grid-points
397	ztemp = REAL( nhor )
398	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
399	IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting
400	CALL ctl_stop( ' cross land advection at Hormuz does not work with your processor cutting: change it' )
401	ENDIF
402	!
403	END SUBROUTINE cla_init_tam
404
405	SUBROUTINE cla_bab_el_mandeb_tan( cd_td )
406	!!----------------------------------------------------------------------
407	!! * ROUTINE cla_bab_el_mandeb_tan *
408	!!
409	!! ** Purpose : update the now horizontal divergence, the tracer tendancy
410	!! and the after velocity in vicinity of Bab el Mandeb ( Red Sea - Indian ocean).
411	!!
412	!! ** Method : compute the exchanges at each side of the strait :
413	!!
414	!! surf. zio_flow
415	!! (+ balance of emp) /\ \|\\\\\\\\\\\\|
416	!! \|\| \|\\\\\\\\\\\\|
417	!! deep zio_flow \|\| \|\\\\\\\\\\\\|
418	!! \| \|\| \|\| \|\\\\\\\\\\\\|
419	!! 89 \| \|\| \|\| \|\\\\\\\\\\\\|
420	!! \|__\/_v_\|\|__\|____________
421	!! !\\\\\\\\\\\\| surf. zio_flow
422	!! \|\\\\\\\\\\\\|<=== (+ balance of emp)
423	!! \|\\\\\\\\\\\u
424	!! 88 \|\\\\\\\\\\\\|<--- deep zrecirc (upper+deep at 2 different levels)
425	!! \|___________\|__________
426	!! !\\\\\\\\\\\\|
427	!! \|\\\\\\\\\\\\| ---\ deep zrecirc (upper+deep)
428	!! 87 !\\\\\\\\\\\u ===/ + deep zio_flow (all at the same level)
429	!! !\\\\\\\\\\\\|
430	!! !___________\|__________
431	!! 160 161
432	!!
433	!!----------------------------------------------------------------------
434	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
435	! ! ='tra' update the tracers
436	! ! ='spg' update after velocity
437	INTEGER :: ji, jj, jk ! dummy loop indices
438	REAL(wp) :: zemp_red_tl ! temporary scalar
439	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
440	!!---------------------------------------------------------------------
441	!
442	SELECT CASE( cd_td )
443	! ! ---------------- !
444	CASE( 'ini' ) ! initialisation !
445	! ! ---------------- !
446	!
447	zio_flow = 0.4e6 ! imposed in/out flow
448	zrecirc_upp = 0.2e6 ! imposed upper recirculation water
449	zrecirc_bot = 0.5e6 ! imposed bottom recirculation water
450
451	hdiv_161_88(:) = 0.e0 ! (161,88) Gulf of Aden side, north point
452	hdiv_161_87(:) = 0.e0 ! (161,87) Gulf of Aden side, south point
453	hdiv_160_89(:) = 0.e0 ! (160,89) Red sea side
454	hdiv_161_88_tl(:) = 0.e0 ! (161,88) Gulf of Aden side, north point
455	hdiv_161_87_tl(:) = 0.e0 ! (161,87) Gulf of Aden side, south point
456	hdiv_160_89_tl(:) = 0.e0 ! (160,89) Red sea side
457
458	DO jj = mj0(88), mj1(88) !** profile of hdiv at (161,88) (Gulf of Aden side, north point)
459	DO ji = mi0(161), mi1(161) !------------------------------
460	DO jk = 1, 8 ! surface in/out flow (Ind -> Red) (div >0)
461	hdiv_161_88(jk) = + zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
462	END DO
463	! ! recirculation water (Ind -> Red) (div >0)
464	hdiv_161_88(20) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,20) )
465	hdiv_161_88(21) = + ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
466	END DO
467	END DO
468	!
469	DO jj = mj0(87), mj1(87) !** profile of hdiv at (161,88) (Gulf of Aden side, south point)
470	DO ji = mi0(161), mi1(161) !------------------------------
471	! ! deep out flow + recirculation (Red -> Ind) (div <0)
472	hdiv_161_87(21) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
473	END DO
474	END DO
475	!
476	DO jj = mj0(89), mj1(89) !** profile of hdiv at (161,88) (Red sea side)
477	DO ji = mi0(160), mi1(160) !------------------------------
478	DO jk = 1, 8 ! surface inflow (Ind -> Red) (div <0)
479	hdiv_160_89(jk) = - zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
480	END DO
481	! ! deep outflow (Red -> Ind) (div >0)
482	hdiv_160_89(16) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,16) )
483	END DO
484	END DO
485	! ! ---------------- !
486	CASE( 'div' ) ! update hdivn ! (call by divcur module)
487	! ! ---------=====-- !
488	! !** emp on the Red Sea (div >0)
489	zemp_red_tl = 0.e0 !---------------------
490	DO jj = mj0(87), mj1(96) ! sum over the Red sea
491	DO ji = mi0(148), mi1(160)
492	zemp_red_tl = zemp_red_tl + emp_tl(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
493	END DO
494	END DO
495	IF( lk_mpp ) CALL mpp_sum( zemp_red_tl ) ! sum with other processors value
496	zemp_red_tl = zemp_red_tl * 1.e-3 ! convert in m3
497	!
498	! !** Correct hdivn (including emp adjustment)
499	! !-------------------------------------------
500	DO jj = mj0(88), mj1(88) !* profile of hdiv at (161,88) (Gulf of Aden side, north point)
501	DO ji = mi0(161), mi1(161)
502	hdiv_161_88_kt_tl(:) = hdiv_161_88_tl(:)
503	DO jk = 1, 8 ! increase the inflow from the Indian (div >0)
504	hdiv_161_88_kt_tl(jk) = hdiv_161_88_tl(jk) + zemp_red_tl / (8. * e2u(ji,jj) * fse3u(ji,jj,jk) )
505	END DO
506	hdivn_tl(ji,jj,:) = hdivn_tl(ji,jj,:) + hdiv_161_88_kt_tl(:)
507	END DO
508	END DO
509	DO jj = mj0(87), mj1(87) !* profile of divergence at (161,87) (Gulf of Aden side, south point)
510	DO ji = mi0(161), mi1(161)
511	hdivn_tl(ji,jj,:) = hdivn_tl(ji,jj,:) + hdiv_161_87_tl(:)
512	END DO
513	END DO
514	DO jj = mj0(89), mj1(89) !* profile of divergence at (160,89) (Red sea side)
515	DO ji = mi0(160), mi1(160)
516	hdiv_160_89_kt_tl(:) = hdiv_160_89_tl(:)
517	DO jk = 1, 18 ! increase the inflow from the Indian (div <0)
518	hdiv_160_89_kt_tl(jk) = hdiv_160_89_tl(jk) - zemp_red_tl / (10. * e1v(ji,jj) * fse3v(ji,jj,jk) )
519	END DO
520	hdivn_tl(ji, jj,:) = hdivn_tl(ji, jj,:) + hdiv_160_89_kt_tl(:)
521	END DO
522	END DO
523	! ! ---------------- !
524	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
525	! ! --------=======- !
526	!
527	DO jj = mj0(88), mj1(88) !** (161,88) (Gulf of Aden side, north point)
528	DO ji = mi0(161), mi1(161)
529	DO jk = 1, jpkm1 ! surf inflow + reciculation (from Gulf of Aden)
530	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
531	& - hdiv_161_88_kt_tl(jk) * tsn(ji,jj,jk,jp_tem) &
532	& - hdiv_161_88_kt(jk) * tsn_tl(ji,jj,jk,jp_tem)
533	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
534	& - hdiv_161_88_kt_tl(jk) * tsn(ji,jj,jk,jp_sal) &
535	& - hdiv_161_88_kt(jk) * tsn_tl(ji,jj,jk,jp_sal)
536	END DO
537	END DO
538	END DO
539	DO jj = mj0(87), mj1(87) !** (161,87) (Gulf of Aden side, south point)
540	DO ji = mi0(161), mi1(161)
541	jk = 21 ! deep outflow + recirulation (combined flux)
542	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
543	& + hdiv_161_88_tl(20) * tsn(ji ,jj+1,20,jp_tem) & ! upper recirculation from Gulf of Aden
544	& + hdiv_161_88(20) * tsn_tl(ji ,jj+1,20,jp_tem) & ! upper recirculation from Gulf of Aden
545	& + hdiv_161_88_tl(21) * tsn(ji ,jj+1,21,jp_tem) & ! deep recirculation from Gulf of Aden
546	& + hdiv_161_88(21) * tsn_tl(ji ,jj+1,21,jp_tem) & ! deep recirculation from Gulf of Aden
547	& + hdiv_160_89_tl(16) * tsn(ji-1,jj+2,16,jp_tem) & ! deep inflow from Red sea
548	& + hdiv_160_89(16) * tsn_tl(ji-1,jj+2,16,jp_tem) ! deep inflow from Red sea
549	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
550	& + hdiv_161_88_tl(20) * tsn(ji ,jj+1,20,jp_sal) &
551	& + hdiv_161_88(20) * tsn_tl(ji ,jj+1,20,jp_sal) &
552	& + hdiv_161_88_tl(21) * tsn(ji ,jj+1,21,jp_sal) &
553	& + hdiv_161_88(21) * tsn_tl(ji ,jj+1,21,jp_sal) &
554	& + hdiv_160_89_tl(16) * tsn(ji-1,jj+2,16,jp_sal) &
555	& + hdiv_160_89(16) * tsn_tl(ji-1,jj+2,16,jp_sal)
556	END DO
557	END DO
558	DO jj = mj0(89), mj1(89) !** (161,88) (Red sea side)
559	DO ji = mi0(160), mi1(160)
560	DO jk = 1, 14 ! surface inflow (from Gulf of Aden)
561	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
562	& - hdiv_160_89_kt_tl(jk) * tsn(ji+1,jj-1,jk,jp_tem) &
563	& - hdiv_160_89_kt(jk) * tsn_tl(ji+1,jj-1,jk,jp_tem)
564	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
565	& - hdiv_160_89_kt_tl(jk) * tsn(ji+1,jj-1,jk,jp_sal) &
566	& - hdiv_160_89_kt(jk) * tsn_tl(ji+1,jj-1,jk,jp_sal)
567	END DO
568	! ! deep outflow (from Red sea)
569	tsa_tl(ji,jj,16,jp_tem) = tsa_tl(ji,jj,16,jp_tem) &
570	& - hdiv_160_89_tl(16) * tsn(ji,jj,16,jp_tem) &
571	& - hdiv_160_89(16) * tsn_tl(ji,jj,16,jp_tem)
572	tsa_tl(ji,jj,16,jp_sal) = tsa_tl(ji,jj,16,jp_sal) &
573	& - hdiv_160_89_tl(16) * tsn(ji,jj,16,jp_sal) &
574	& - hdiv_160_89(16) * tsn_tl(ji,jj,16,jp_sal)
575	END DO
576	END DO
577	!
578	! ! ---------------- !
579	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
580	! ! --------=======- !
581	! at this stage, (ua,va) are the after velocity, not the tendancy
582	! compute the velocity from the divergence at T-point
583	!
584	DO jj = mj0(88), mj1(88) !** (160,88) (Gulf of Aden side, north point)
585	DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point)
586	ua_tl(ji,jj,:) = - hdiv_161_88_kt_tl(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
587	& * e2u(ji,jj) * fse3u(ji,jj,:)
588	END DO
589	END DO
590	DO jj = mj0(87), mj1(87) !** (160,87) (Gulf of Aden side, south point)
591	DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point)
592	ua_tl(ji,jj,:) = - hdiv_161_87_tl(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
593	& * e2u(ji,jj) * fse3u(ji,jj,:)
594	END DO
595	END DO
596	DO jj = mj0(88), mj1(88) !** profile of divergence at (160,89) (Red sea side)
597	DO ji = mi0(160), mi1(160) ! 88, not 89 as it is a V-point)
598	va_tl(ji,jj,:) = - hdiv_160_89_kt_tl(:) / ( e1t(ji,jj+1) * e2t(ji,jj+1) * fse3t(ji,jj+1,:) ) &
599	& * e1v(ji,jj) * fse3v(ji,jj,:)
600	END DO
601	END DO
602	END SELECT
603	!
604	END SUBROUTINE cla_bab_el_mandeb_tan
605
606	SUBROUTINE cla_bab_el_mandeb_adj( cd_td )
607	!!----------------------------------------------------------------------
608	!! * ROUTINE cla_bab_el_mandeb_adj *
609	!!
610	!! ** Purpose : update the now horizontal divergence, the tracer tendancy
611	!! and the after velocity in vicinity of Bab el Mandeb ( Red Sea - Indian ocean).
612	!!
613	!! ** Method : compute the exchanges at each side of the strait :
614	!!
615	!! surf. zio_flow
616	!! (+ balance of emp) /\ \|\\\\\\\\\\\\|
617	!! \|\| \|\\\\\\\\\\\\|
618	!! deep zio_flow \|\| \|\\\\\\\\\\\\|
619	!! \| \|\| \|\| \|\\\\\\\\\\\\|
620	!! 89 \| \|\| \|\| \|\\\\\\\\\\\\|
621	!! \|__\/_v_\|\|__\|____________
622	!! !\\\\\\\\\\\\| surf. zio_flow
623	!! \|\\\\\\\\\\\\|<=== (+ balance of emp)
624	!! \|\\\\\\\\\\\u
625	!! 88 \|\\\\\\\\\\\\|<--- deep zrecirc (upper+deep at 2 different levels)
626	!! \|___________\|__________
627	!! !\\\\\\\\\\\\|
628	!! \|\\\\\\\\\\\\| ---\ deep zrecirc (upper+deep)
629	!! 87 !\\\\\\\\\\\u ===/ + deep zio_flow (all at the same level)
630	!! !\\\\\\\\\\\\|
631	!! !___________\|__________
632	!! 160 161
633	!!
634	!!----------------------------------------------------------------------
635	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
636	! ! ='tra' update the tracers
637	! ! ='spg' update after velocity
638	INTEGER :: ji, jj, jk ! dummy loop indices
639	REAL(wp) :: zemp_red_ad ! temporary scalar
640	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
641	!!---------------------------------------------------------------------
642	!
643	SELECT CASE( cd_td )
644	! ! ---------------- !
645	CASE( 'ini' ) ! initialisation !
646	! ! ---------------- !
647	!
648	zio_flow = 0.4e6 ! imposed in/out flow
649	zrecirc_upp = 0.2e6 ! imposed upper recirculation water
650	zrecirc_bot = 0.5e6 ! imposed bottom recirculation water
651
652	hdiv_161_88(:) = 0.e0 ! (161,88) Gulf of Aden side, north point
653	hdiv_161_87(:) = 0.e0 ! (161,87) Gulf of Aden side, south point
654	hdiv_160_89(:) = 0.e0 ! (160,89) Red sea side
655	hdiv_161_88_ad(:) = 0.e0 ! (161,88) Gulf of Aden side, north point
656	hdiv_161_87_ad(:) = 0.e0 ! (161,87) Gulf of Aden side, south point
657	hdiv_160_89_ad(:) = 0.e0 ! (160,89) Red sea side
658
659	DO jj = mj0(88), mj1(88) !** profile of hdiv at (161,88) (Gulf of Aden side, north point)
660	DO ji = mi0(161), mi1(161) !------------------------------
661	DO jk = 1, 8 ! surface in/out flow (Ind -> Red) (div >0)
662	hdiv_161_88(jk) = + zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
663	END DO
664	! ! recirculation water (Ind -> Red) (div >0)
665	hdiv_161_88(20) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,20) )
666	hdiv_161_88(21) = + ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
667	END DO
668	END DO
669	!
670	DO jj = mj0(87), mj1(87) !** profile of hdiv at (161,88) (Gulf of Aden side, south point)
671	DO ji = mi0(161), mi1(161) !------------------------------
672	! ! deep out flow + recirculation (Red -> Ind) (div <0)
673	hdiv_161_87(21) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
674	END DO
675	END DO
676	!
677	DO jj = mj0(89), mj1(89) !** profile of hdiv at (161,88) (Red sea side)
678	DO ji = mi0(160), mi1(160) !------------------------------
679	DO jk = 1, 8 ! surface inflow (Ind -> Red) (div <0)
680	hdiv_160_89(jk) = - zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
681	END DO
682	! ! deep outflow (Red -> Ind) (div >0)
683	hdiv_160_89(16) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,16) )
684	END DO
685	END DO
686	! ! ---------------- !
687	CASE( 'div' ) ! update hdivn ! (call by divcur module)
688	! ! ---------=====-- !
689	! !** emp on the Red Sea (div >0)
690	zemp_red_ad = 0.e0 !---------------------
691	DO jj = mj1(89), mj0(89), -1 !* profile of divergence at (160,89) (Red sea side)
692	DO ji = mi1(160), mi0(160), -1
693	hdiv_160_89_kt_ad(:) = hdiv_160_89_kt_ad(:) + hdivn_ad(ji, jj,:)
694	DO jk = 18, 1, -1 ! increase the inflow from the Indian (div <0)
695	zemp_red_ad = zemp_red_ad - hdiv_160_89_ad(jk) / (10. * e1v(ji,jj) * fse3v(ji,jj,jk) )
696	END DO
697	hdiv_160_89_ad(:) = hdiv_160_89_ad(:) + hdiv_160_89_kt_ad(:)
698	END DO
699	END DO
700	DO jj = mj0(87), mj1(87) !* profile of divergence at (161,87) (Gulf of Aden side, south point)
701	DO ji = mi0(161), mi1(161)
702	hdiv_161_87_ad(:) = hdiv_161_87_ad(:) + hdivn_ad(ji,jj,:)
703	END DO
704	END DO
705	! !** Correct hdivn (including emp adjustment)
706	! !-------------------------------------------
707	DO jj = mj1(88), mj0(88), -1 !* profile of hdiv at (161,88) (Gulf of Aden side, north point)
708	DO ji = mi1(161), mi0(161), -1
709	hdiv_161_88_kt_ad(:) = hdiv_161_88_kt_ad(:) + hdivn_ad(ji,jj,:)
710	DO jk = 8, 1, -1 ! increase the inflow from the Indian (div >0)
711	zemp_red_ad = zemp_red_ad + hdiv_161_88_ad(jk) / (8. * e2u(ji,jj) * fse3u(ji,jj,jk) )
712	END DO
713	hdiv_161_88_ad(:) = hdiv_161_88_ad(:) + hdiv_161_88_kt_ad(:)
714	END DO
715	END DO
716	zemp_red_ad = zemp_red_ad * 1.e-3 ! convert in m3
717	IF( lk_mpp ) CALL mpp_sum( zemp_red_ad ) ! sum with other processors value
718
719	DO jj = mj1(96), mj0(87), -1 ! sum over the Red sea
720	DO ji = mi1(160), mi0(148), -1
721	emp_ad(ji,jj) = emp_ad(ji,jj) + zemp_red_ad * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
722	END DO
723	END DO
724	!
725	! ! ---------------- !
726	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
727	! ! --------=======- !
728	!===========================
729	! Direct model recomputation
730	!===========================
731	!
732	DO jj = mj0(88), mj1(88) !** (161,88) (Gulf of Aden side, north point)
733	DO ji = mi0(161), mi1(161)
734	DO jk = 1, jpkm1 ! surf inflow + reciculation (from Gulf of Aden)
735	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - hdiv_161_88_kt(jk) * tsn(ji,jj,jk,jp_tem)
736	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - hdiv_161_88_kt(jk) * tsn(ji,jj,jk,jp_sal)
737	END DO
738	END DO
739	END DO
740	DO jj = mj0(87), mj1(87) !** (161,87) (Gulf of Aden side, south point)
741	DO ji = mi0(161), mi1(161)
742	jk = 21 ! deep outflow + recirulation (combined flux)
743	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + hdiv_161_88(20) * tsn(ji ,jj+1,20,jp_tem) & ! upper recirculation from Gulf of Aden
744	& + hdiv_161_88(21) * tsn(ji ,jj+1,21,jp_tem) & ! deep recirculation from Gulf of Aden
745	& + hdiv_160_89(16) * tsn(ji-1,jj+2,16,jp_tem) ! deep inflow from Red sea
746	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + hdiv_161_88(20) * tsn(ji ,jj+1,20,jp_sal) &
747	& + hdiv_161_88(21) * tsn(ji ,jj+1,21,jp_sal) &
748	& + hdiv_160_89(16) * tsn(ji-1,jj+2,16,jp_sal)
749	END DO
750	END DO
751	DO jj = mj0(89), mj1(89) !** (161,88) (Red sea side)
752	DO ji = mi0(160), mi1(160)
753	DO jk = 1, 14 ! surface inflow (from Gulf of Aden)
754	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - hdiv_160_89_kt(jk) * tsn(ji+1,jj-1,jk,jp_tem)
755	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - hdiv_160_89_kt(jk) * tsn(ji+1,jj-1,jk,jp_sal)
756	END DO
757	! ! deep outflow (from Red sea)
758	tsa(ji,jj,16,jp_tem) = tsa(ji,jj,16,jp_tem) - hdiv_160_89(16) * tsn(ji,jj,16,jp_tem)
759	tsa(ji,jj,16,jp_sal) = tsa(ji,jj,16,jp_sal) - hdiv_160_89(16) * tsn(ji,jj,16,jp_sal)
760	END DO
761	END DO
762	!=============
763	! Adjoint part
764	!=============
765	!
766	DO jj = mj1(89), mj0(89), -1 !** (161,88) (Red sea side)
767	DO ji = mi1(160), mi0(160), -1
768	hdiv_160_89_ad(16) = hdiv_160_89_ad(16) - tsa_ad(ji,jj,16,jp_sal) * tsn(ji,jj,16,jp_sal)
769	tsn_ad(ji,jj,16,jp_sal) = tsn_ad(ji,jj,16,jp_sal) - tsa_ad(ji,jj,16,jp_sal) * tsn(ji,jj,16,jp_sal)
770	hdiv_160_89_ad(16) = hdiv_160_89_ad(16) - tsa_ad(ji,jj,16,jp_tem) * tsn(ji,jj,16,jp_tem)
771	tsn_ad(ji,jj,16,jp_tem) = tsn_ad(ji,jj,16,jp_tem) - tsa_ad(ji,jj,16,jp_tem) * tsn(ji,jj,16,jp_tem)
772	DO jk = 14, 1, -1 ! surface inflow (from Gulf of Aden)
773	hdiv_160_89_kt_ad(jk) = hdiv_160_89_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji+1,jj-1,jk,jp_sal)
774	tsn_ad(ji+1,jj-1,jk,jp_sal) = tsn_ad(ji+1,jj-1,jk,jp_sal) - tsa_ad(ji,jj,jk,jp_sal) * hdiv_160_89_kt(jk)
775	hdiv_160_89_kt_ad(jk) = hdiv_160_89_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji+1,jj-1,jk,jp_tem)
776	tsn_ad(ji+1,jj-1,jk,jp_tem) = tsn_ad(ji+1,jj-1,jk,jp_tem) - tsa_ad(ji,jj,jk,jp_tem) * hdiv_160_89_kt(jk)
777	END DO
778	END DO
779	END DO
780	DO jj = mj1(87), mj0(87), -1 !** (161,87) (Gulf of Aden side, south point)
781	DO ji = mi1(161), mi0(161), -1
782	jk = 21 ! deep outflow + recirulation (combined flux)
783	hdiv_161_88_ad(20) = hdiv_161_88_ad(20) + tsa_ad(ji,jj,jk,jp_sal) * tsn(ji ,jj+1,20,jp_sal)
784	hdiv_161_88_ad(21) = hdiv_161_88_ad(21) + tsa_ad(ji,jj,jk,jp_sal) * tsn(ji ,jj+1,21,jp_sal)
785	hdiv_160_89_ad(16) = hdiv_160_89_ad(16) + tsa_ad(ji,jj,jk,jp_sal) * tsn(ji-1,jj+2,16,jp_sal)
786	tsn_ad(ji ,jj+1,20,jp_sal) = tsn_ad(ji ,jj+1,20,jp_sal) + tsa_ad(ji,jj,jk,jp_sal) * hdiv_161_88_ad(20)
787	tsn_ad(ji ,jj+1,21,jp_sal) = tsn_ad(ji ,jj+1,21,jp_sal) + tsa_ad(ji,jj,jk,jp_sal) * hdiv_161_88_ad(21)
788	tsn_ad(ji-1,jj+2,16,jp_sal) = tsn_ad(ji-1,jj+2,16,jp_sal) + tsa_ad(ji,jj,jk,jp_sal) * hdiv_160_89_ad(16)
789	hdiv_161_88_ad(20) = hdiv_161_88_ad(20) + tsa_ad(ji,jj,jk,jp_tem) * tsn(ji ,jj+1,20,jp_tem)
790	hdiv_161_88_ad(21) = hdiv_161_88_ad(21) + tsa_ad(ji,jj,jk,jp_tem) * tsn(ji ,jj+1,21,jp_tem)
791	hdiv_160_89_ad(16) = hdiv_160_89_ad(16) + tsa_ad(ji,jj,jk,jp_tem) * tsn(ji-1,jj+2,16,jp_tem)
792	tsn_ad(ji ,jj+1,20,jp_tem) = tsn_ad(ji ,jj+1,20,jp_tem) + tsa_ad(ji,jj,jk,jp_tem) * hdiv_161_88_ad(20)
793	tsn_ad(ji ,jj+1,21,jp_tem) = tsn_ad(ji ,jj+1,21,jp_tem) + tsa_ad(ji,jj,jk,jp_tem) * hdiv_161_88_ad(21)
794	tsn_ad(ji-1,jj+2,16,jp_tem) = tsn_ad(ji-1,jj+2,16,jp_tem) + tsa_ad(ji,jj,jk,jp_tem) * hdiv_160_89_ad(16)
795	END DO
796	END DO
797	DO jj = mj1(88), mj0(88), -1 !** (161,88) (Gulf of Aden side, north point)
798	DO ji = mi1(161), mi0(161), -1
799	DO jk = jpkm1, 1, -1 ! surf inflow + reciculation (from Gulf of Aden)
800	hdiv_161_88_kt_ad(jk) = hdiv_161_88_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji,jj,jk,jp_sal)
801	tsn_ad(ji,jj,jk,jp_sal) = tsn_ad(ji,jj,jk,jp_sal) - hdiv_161_88_kt(jk) * tsa_ad(ji,jj,jk,jp_sal)
802	hdiv_161_88_kt_ad(jk) = hdiv_161_88_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji,jj,jk,jp_tem)
803	tsn_ad(ji,jj,jk,jp_tem) = tsn_ad(ji,jj,jk,jp_tem) - hdiv_161_88_kt(jk) * tsa_ad(ji,jj,jk,jp_tem)
804	END DO
805	END DO
806	END DO
807	!
808	! ! ---------------- !
809	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
810	! ! --------=======- !
811	! at this stage, (ua,va) are the after velocity, not the tendancy
812	! compute the velocity from the divergence at T-point
813	!
814	DO jj = mj1(88), mj0(88), -1 !** profile of divergence at (160,89) (Red sea side)
815	DO ji = mi1(160), mi0(160), -1 ! 88, not 89 as it is a V-point)
816	hdiv_160_89_kt_ad(:) = hdiv_160_89_kt_ad(:) - va_ad(ji,jj,:) / ( e1t(ji,jj+1) * e2t(ji,jj+1) * fse3t(ji,jj+1,:) ) &
817	& * e1v(ji,jj) * fse3v(ji,jj,:)
818	va_ad(ji,jj,:) = 0.0_wp
819	END DO
820	END DO
821	DO jj = mj1(87), mj0(87), -1 !** (160,87) (Gulf of Aden side, south point)
822	DO ji = mi1(160), mi0(160), -1 ! 160, not 161 as it is a U-point)
823	hdiv_161_87_ad(:) = hdiv_161_87_ad(:) - ua_ad(ji,jj,:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
824	& * e2u(ji,jj) * fse3u(ji,jj,:)
825	ua_ad(ji,jj,:) = 0.0_wp
826	END DO
827	END DO
828	DO jj = mj1(88), mj0(88), -1 !** (160,88) (Gulf of Aden side, north point)
829	DO ji = mi1(160), mi0(160), -1 ! 160, not 161 as it is a U-point)
830	hdiv_161_88_kt_ad(:) = hdiv_161_88_kt_ad(:) - ua_ad(ji,jj,:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
831	& * e2u(ji,jj) * fse3u(ji,jj,:)
832	ua_ad(ji,jj,:) = 0.0_wp
833	END DO
834	END DO
835	END SELECT
836	!
837	END SUBROUTINE cla_bab_el_mandeb_adj
838
839	SUBROUTINE cla_gibraltar_tan( cd_td )
840	!! -------------------------------------------------------------------
841	!! * ROUTINE cla_gibraltar_tan *
842	!!
843	!! ** Purpose : update the now horizontal divergence, the tracer
844	!! tendancyand the after velocity in vicinity of Gibraltar
845	!! strait ( Persian Gulf - Indian ocean ).
846	!!
847	!! ** Method :
848	!! _______________________
849	!! deep zio_flow /====\|///////\|====> surf. zio_flow
850	!! + deep zrecirc \----\|///////\| (+balance of emp)
851	!! 102 u///////u
852	!! mid. recicul <--\|///////\|<==== deep zio_flow
853	!! _____\|_______\|_____
854	!! surf. zio_flow ====>\|///////\|
855	!! (+balance of emp) \|///////\|
856	!! 101 u///////\|
857	!! mid. recicul -->\|///////\| Caution: zrecirc split into
858	!! deep zrecirc ---->\|///////\| upper & bottom recirculation
859	!! _______\|_______\|_______
860	!! 139 140 141
861	!!
862	!!---------------------------------------------------------------------
863	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
864	! ! ='tra' update the tracers
865	! ! ='spg' update after velocity
866	INTEGER :: ji, jj, jk ! dummy loop indices
867	REAL(wp) :: zemp_med ! temporary scalar
868	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
869	!!---------------------------------------------------------------------
870	!
871	SELECT CASE( cd_td )
872	! ! ---------------- !
873	CASE( 'ini' ) ! initialisation !
874	! ! ---------------- !
875	! !** initialization of the velocity
876	hdiv_139_101(:) = 0.e0 ! 139,101 (Atlantic side, south point)
877	hdiv_139_102(:) = 0.e0 ! 139,102 (Atlantic side, north point)
878	hdiv_141_102(:) = 0.e0 ! 141,102 (Med sea side)
879	hdiv_139_101_tl(:) = 0.e0 ! 139,101 (Atlantic side, south point)
880	hdiv_139_102_tl(:) = 0.e0 ! 139,102 (Atlantic side, north point)
881	hdiv_141_102_tl(:) = 0.e0 ! 141,102 (Med sea side)
882
883	! !** imposed transport
884	zio_flow = 0.8e6 ! inflow surface water
885	zrecirc_mid = 0.7e6 ! middle recirculation water
886	zrecirc_upp = 2.5e6 ! upper recirculation water
887	zrecirc_bot = 3.5e6 ! bottom recirculation water
888	!
889	DO jj = mj0(101), mj1(101) !** profile of hdiv at 139,101 (Atlantic side, south point)
890	DO ji = mi0(139), mi1(139) !-----------------------------
891	DO jk = 1, 14 ! surface in/out flow (Atl -> Med) (div >0)
892	hdiv_139_101(jk) = + zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
893	END DO
894	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div >0)
895	hdiv_139_101(jk) = + zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
896	END DO
897	! ! upper reciculation (Atl 101 -> Atl 101) (div >0)
898	hdiv_139_101(21) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
899	!
900	! ! upper & bottom reciculation (Atl 101 -> Atl 101 & 102) (div >0)
901	hdiv_139_101(22) = ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
902	END DO
903	END DO
904	DO jj = mj0(102), mj1(102) !** profile of hdiv at 139,102 (Atlantic side, north point)
905	DO ji = mi0(139), mi1(139) !-----------------------------
906	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
907	hdiv_139_102(jk) = - zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
908	END DO
909	! ! outflow of Mediterranean sea + deep recirculation (div <0)
910	hdiv_139_102(22) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
911	END DO
912	END DO
913	DO jj = mj0(102), mj1(102) !** velocity profile at 141,102 (Med sea side)
914	DO ji = mi0(141), mi1(141) !------------------------------
915	DO jk = 1, 14 ! surface inflow in the Med (div <0)
916	hdiv_141_102(jk) = - zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
917	END DO
918	! ! deep outflow toward the Atlantic (div >0)
919	hdiv_141_102(21) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
920	END DO
921	END DO
922	! ! ---------------- !
923	CASE( 'div' ) ! update hdivn ! (call by divcur module)
924	! ! ---------=====-- !
925	! !** emp on the Mediterranean Sea (div >0)
926	zemp_med = 0.e0 !-------------------------------
927	DO jj = mj0(96), mj1(110) ! sum over the Med sea
928	DO ji = mi0(141),mi1(181)
929	zemp_med = zemp_med + emp_tl(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
930	END DO
931	END DO
932	DO jj = mj0(96), mj1(96) ! minus 2 points in Red Sea
933	DO ji = mi0(148),mi1(148)
934	zemp_med = zemp_med - emp_tl(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
935	END DO
936	DO ji = mi0(149),mi1(149)
937	zemp_med = zemp_med - emp_tl(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
938	END DO
939	END DO
940	IF( lk_mpp ) CALL mpp_sum( zemp_med ) ! sum with other processors value
941	zemp_med = zemp_med * 1.e-3 ! convert in m3
942	!
943	! !** Correct hdivn (including emp adjustment)
944	! !-------------------------------------------
945	DO jj = mj0(101), mj1(101) !* 139,101 (Atlantic side, south point)
946	DO ji = mi0(139), mi1(139)
947	hdiv_139_101_kt_tl(:) = hdiv_139_101_tl(:)
948	DO jk = 1, 14 ! increase the inflow from the Atlantic (div >0)
949	hdiv_139_101_kt_tl(jk) = hdiv_139_101_tl(jk) + zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
950	END DO
951	hdivn_tl(ji, jj,:) = hdivn_tl(ji, jj,:) + hdiv_139_101_kt_tl(:)
952	END DO
953	END DO
954	DO jj = mj0(102), mj1(102) !* 139,102 (Atlantic side, north point)
955	DO ji = mi0(139), mi1(139)
956	hdivn_tl(ji,jj,:) = hdivn_tl(ji,jj,:) + hdiv_139_102_tl(:)
957	END DO
958	END DO
959	DO jj = mj0(102), mj1(102) !* 141,102 (Med side)
960	DO ji = mi0(141), mi1(141)
961	hdiv_141_102_tl(:) = hdiv_141_102_tl(:)
962	DO jk = 1, 14 ! increase the inflow from the Atlantic (div <0)
963	hdiv_141_102_kt_tl(jk) = hdiv_141_102_tl(jk) - zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
964	END DO
965	hdivn_tl(ji, jj,:) = hdivn_tl(ji, jj,:) + hdiv_141_102_kt_tl(:)
966	END DO
967	END DO
968	! ! ---------------- !
969	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
970	! ! --------=======- !
971	!
972	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point) (div >0)
973	DO ji = mi0(139), mi1(139)
974	DO jk = 1, jpkm1 ! surf inflow + mid. & bottom reciculation (from Atlantic)
975	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
976	& - hdiv_139_101_kt_tl(jk) * tsn(ji,jj,jk,jp_tem) &
977	& - hdiv_139_101_kt(jk) * tsn_tl(ji,jj,jk,jp_tem)
978	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
979	& - hdiv_139_101_kt_tl(jk) * tsn(ji,jj,jk,jp_sal) &
980	& - hdiv_139_101_kt(jk) * tsn_tl(ji,jj,jk,jp_sal)
981	END DO
982	END DO
983	END DO
984	!
985	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point) (div <0)
986	DO ji = mi0(139), mi1(139)
987	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
988	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
989	& - hdiv_139_102_tl(jk) * tsn(ji,jj-1,jk,jp_tem) & ! middle Atlantic recirculation
990	& - hdiv_139_102(jk) * tsn_tl(ji,jj-1,jk,jp_tem) ! middle Atlantic recirculation
991	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
992	& - hdiv_139_102_tl(jk) * tsn(ji,jj-1,jk,jp_sal) &
993	& - hdiv_139_102(jk) * tsn_tl(ji,jj-1,jk,jp_sal)
994	END DO
995	! ! upper & bottom Atl. reciculation (Atl 101 -> Atl 102) - (div <0)
996	! ! deep Med flow (Med 102 -> Atl 102) - (div <0)
997	tsa_tl(ji,jj,22,jp_tem) = tsa_tl(ji,jj,22,jp_tem) &
998	& + hdiv_141_102_tl(21) * tsn(ji+2,jj,21,jp_tem) & ! deep Med flow
999	& + hdiv_141_102(21) * tsn_tl(ji+2,jj,21,jp_tem) & ! deep Med flow
1000	& + hdiv_139_101_tl(21) * tsn(ji,jj-1,21,jp_tem) & ! upper Atlantic recirculation
1001	& + hdiv_139_101(21) * tsn_tl(ji,jj-1,21,jp_tem) & ! upper Atlantic recirculation
1002	& + hdiv_139_101_tl(22) * tsn(ji,jj-1,22,jp_tem) & ! bottom Atlantic recirculation
1003	& + hdiv_139_101(22) * tsn_tl(ji,jj-1,22,jp_tem) ! bottom Atlantic recirculation
1004	tsa_tl(ji,jj,22,jp_sal) = tsa_tl(ji,jj,22,jp_sal) &
1005	& + hdiv_141_102_tl(21) * tsn(ji+2,jj,21,jp_sal) &
1006	& + hdiv_141_102(21) * tsn_tl(ji+2,jj,21,jp_sal) &
1007	& + hdiv_139_101_tl(21) * tsn(ji,jj-1,21,jp_sal) &
1008	& + hdiv_139_101(21) * tsn_tl(ji,jj-1,21,jp_sal) &
1009	& + hdiv_139_101_tl(22) * tsn(ji,jj-1,22,jp_sal) &
1010	& + hdiv_139_101(22) * tsn_tl(ji,jj-1,22,jp_sal)
1011	END DO
1012	END DO
1013	DO jj = mj0(102), mj1(102) !* 141,102 (Med side) (div <0)
1014	DO ji = mi0(141), mi1(141)
1015	DO jk = 1, 14 ! surface flow from Atlantic to Med sea
1016	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
1017	& - hdiv_141_102_kt_tl(jk) * tsn(ji-2,jj-1,jk,jp_tem) &
1018	& - hdiv_141_102_kt(jk) * tsn_tl(ji-2,jj-1,jk,jp_tem)
1019	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
1020	& - hdiv_141_102_kt_tl(jk) * tsn(ji-2,jj-1,jk,jp_sal) &
1021	& - hdiv_141_102_kt(jk) * tsn_tl(ji-2,jj-1,jk,jp_sal)
1022	END DO
1023	! ! deeper flow from Med sea to Atlantic
1024	tsa_tl(ji,jj,21,jp_tem) = tsa_tl(ji,jj,21,jp_tem) &
1025	& - hdiv_141_102_tl(21) * tsn(ji,jj,21,jp_tem) &
1026	& - hdiv_141_102(21) * tsn_tl(ji,jj,21,jp_tem)
1027	tsa_tl(ji,jj,21,jp_sal) = tsa_tl(ji,jj,21,jp_sal) &
1028	& - hdiv_141_102_tl(21) * tsn(ji,jj,21,jp_sal) &
1029	& - hdiv_141_102(21) * tsn_tl(ji,jj,21,jp_sal)
1030	END DO
1031	END DO
1032	! ! ---------------- !
1033	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
1034	! ! --------=======- !
1035	! at this stage, (ua,va) are the after velocity, not the tendancy
1036	! compute the velocity from the divergence at T-point
1037	!
1038	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point)
1039	DO ji = mi0(139), mi1(139) ! div >0 => ua >0, same sign
1040	ua_tl(ji,jj,:) = hdiv_139_101_kt_tl(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
1041	& * e2u(ji,jj) * fse3u(ji,jj,:)
1042	END DO
1043	END DO
1044	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point)
1045	DO ji = mi0(139), mi1(139) ! div <0 => ua <0, same sign
1046	ua_tl(ji,jj,:) = hdiv_139_102_tl(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
1047	& * e2u(ji,jj) * fse3u(ji,jj,:)
1048	END DO
1049	END DO
1050	DO jj = mj0(102), mj1(102) !** 140,102 (Med side) (140 not 141 as it is a U-point)
1051	DO ji = mi0(140), mi1(140) ! div >0 => ua <0, opposite sign
1052	ua_tl(ji,jj,:) = - hdiv_141_102_tl(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
1053	& * e2u(ji,jj) * fse3u(ji,jj,:)
1054	END DO
1055	END DO
1056	!
1057	END SELECT
1058	!
1059	END SUBROUTINE cla_gibraltar_tan
1060
1061	SUBROUTINE cla_gibraltar_adj( cd_td )
1062	!! -------------------------------------------------------------------
1063	!! * ROUTINE cla_gibraltar_adj *
1064	!!
1065	!! ** Purpose : update the now horizontal divergence, the tracer
1066	!! tendancyand the after velocity in vicinity of Gibraltar
1067	!! strait ( Persian Gulf - Indian ocean ).
1068	!!
1069	!! ** Method :
1070	!! _______________________
1071	!! deep zio_flow /====\|///////\|====> surf. zio_flow
1072	!! + deep zrecirc \----\|///////\| (+balance of emp)
1073	!! 102 u///////u
1074	!! mid. recicul <--\|///////\|<==== deep zio_flow
1075	!! _____\|_______\|_____
1076	!! surf. zio_flow ====>\|///////\|
1077	!! (+balance of emp) \|///////\|
1078	!! 101 u///////\|
1079	!! mid. recicul -->\|///////\| Caution: zrecirc split into
1080	!! deep zrecirc ---->\|///////\| upper & bottom recirculation
1081	!! _______\|_______\|_______
1082	!! 139 140 141
1083	!!
1084	!!---------------------------------------------------------------------
1085	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
1086	! ! ='tra' update the tracers
1087	! ! ='spg' update after velocity
1088	INTEGER :: ji, jj, jk ! dummy loop indices
1089	REAL(wp) :: zemp_med ! temporary scalar
1090	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
1091	!!---------------------------------------------------------------------
1092	!
1093	SELECT CASE( cd_td )
1094	! ! ---------------- !
1095	CASE( 'ini' ) ! initialisation !
1096	! ! ---------------- !
1097	! !** initialization of the velocity
1098	hdiv_139_101(:) = 0.e0 ! 139,101 (Atlantic side, south point)
1099	hdiv_139_102(:) = 0.e0 ! 139,102 (Atlantic side, north point)
1100	hdiv_141_102(:) = 0.e0 ! 141,102 (Med sea side)
1101	hdiv_139_101_ad(:) = 0.e0 ! 139,101 (Atlantic side, south point)
1102	hdiv_139_102_ad(:) = 0.e0 ! 139,102 (Atlantic side, north point)
1103	hdiv_141_102_ad(:) = 0.e0 ! 141,102 (Med sea side)
1104
1105	! !** imposed transport
1106	zio_flow = 0.8e6 ! inflow surface water
1107	zrecirc_mid = 0.7e6 ! middle recirculation water
1108	zrecirc_upp = 2.5e6 ! upper recirculation water
1109	zrecirc_bot = 3.5e6 ! bottom recirculation water
1110	!
1111	DO jj = mj0(101), mj1(101) !** profile of hdiv at 139,101 (Atlantic side, south point)
1112	DO ji = mi0(139), mi1(139) !-----------------------------
1113	DO jk = 1, 14 ! surface in/out flow (Atl -> Med) (div >0)
1114	hdiv_139_101(jk) = + zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1115	END DO
1116	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div >0)
1117	hdiv_139_101(jk) = + zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1118	END DO
1119	! ! upper reciculation (Atl 101 -> Atl 101) (div >0)
1120	hdiv_139_101(21) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1121	!
1122	! ! upper & bottom reciculation (Atl 101 -> Atl 101 & 102) (div >0)
1123	hdiv_139_101(22) = ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1124	END DO
1125	END DO
1126	DO jj = mj0(102), mj1(102) !** profile of hdiv at 139,102 (Atlantic side, north point)
1127	DO ji = mi0(139), mi1(139) !-----------------------------
1128	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
1129	hdiv_139_102(jk) = - zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1130	END DO
1131	! ! outflow of Mediterranean sea + deep recirculation (div <0)
1132	hdiv_139_102(22) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1133	END DO
1134	END DO
1135	DO jj = mj0(102), mj1(102) !** velocity profile at 141,102 (Med sea side)
1136	DO ji = mi0(141), mi1(141) !------------------------------
1137	DO jk = 1, 14 ! surface inflow in the Med (div <0)
1138	hdiv_141_102(jk) = - zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1139	END DO
1140	! ! deep outflow toward the Atlantic (div >0)
1141	hdiv_141_102(21) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1142	END DO
1143	END DO
1144	! ! ---------------- !
1145	CASE( 'div' ) ! update hdivn ! (call by divcur module)
1146	! ! ---------=====-- !
1147	! !** Correct hdivn (including emp adjustment)
1148	! !-------------------------------------------
1149	DO jj = mj1(102), mj0(102), -1 !* 141,102 (Med side)
1150	DO ji = mi1(141), mi0(141), -1
1151	hdiv_141_102_kt_ad(:) = hdiv_141_102_kt_ad(:) + hdivn_ad(ji, jj,:)
1152	DO jk = 14, 1, -1 ! increase the inflow from the Atlantic (div <0)
1153	zemp_med = zemp_med - hdiv_141_102_kt_ad(jk) / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1154	hdiv_141_102_ad(jk) = hdiv_141_102_ad(jk) + hdiv_141_102_kt_ad(jk)
1155	END DO
1156	hdiv_141_102_ad(:) = hdiv_141_102_ad(:) + hdiv_141_102_kt_ad(:)
1157	END DO
1158	END DO
1159	DO jj = mj1(102), mj0(102), -1 !* 139,102 (Atlantic side, north point)
1160	DO ji = mi1(139), mi0(139), -1
1161	hdiv_139_102_ad(:) = hdiv_139_102_ad(:) + hdivn_ad(ji,jj,:)
1162	END DO
1163	END DO
1164	DO jj = mj1(101), mj0(101), -1 !* 139,101 (Atlantic side, south point)
1165	DO ji = mi1(139), mi0(139), -1
1166	hdiv_139_101_kt_ad(:) = hdiv_139_101_kt_ad(:) + hdivn_ad(ji, jj,:)
1167	DO jk = 14, 1, -1 ! increase the inflow from the Atlantic (div >0)
1168	zemp_med = zemp_med + hdiv_139_101_kt_ad(jk) / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1169	hdiv_139_101_ad(jk) = hdiv_139_101_ad(jk) + hdiv_139_101_kt_ad(jk)
1170	END DO
1171	hdiv_139_101_kt_ad(:) = hdiv_139_101_ad(:) + hdiv_139_101_kt_ad(:)
1172	END DO
1173	END DO
1174	! !** emp on the Mediterranean Sea (div >0)
1175	zemp_med = zemp_med * 1.e-3 ! convert in m3
1176	IF( lk_mpp ) CALL mpp_sum( zemp_med ) ! sum with other processors value
1177	DO jj = mj1(96), mj0(96), -1 ! minus 2 points in Red Sea
1178	DO ji = mi1(149),mi0(149), -1
1179	emp_ad(ji,jj) = emp_ad(ji,jj) - zemp_med * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
1180	END DO
1181	DO ji = mi1(148),mi0(148), -1
1182	emp_ad(ji,jj) = emp_ad(ji,jj) - zemp_med * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
1183	END DO
1184	END DO
1185	DO jj = mj1(110), mj0(96), -1 ! sum over the Med sea
1186	DO ji = mi1(181), mi0(141), -1
1187	emp_ad(ji,jj) = emp_ad(ji,jj) - zemp_med * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
1188	END DO
1189	END DO
1190	zemp_med = 0.e0 !-------------------------------
1191	! ! ---------------- !
1192	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
1193	! ! --------=======- !
1194	!===========================
1195	! Direct model recomputation
1196	!===========================
1197	!
1198	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point) (div >0)
1199	DO ji = mi0(139), mi1(139)
1200	DO jk = 1, jpkm1 ! surf inflow + mid. & bottom reciculation (from Atlantic)
1201	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - hdiv_139_101_kt(jk) * tsn(ji,jj,jk,jp_tem)
1202	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - hdiv_139_101_kt(jk) * tsn(ji,jj,jk,jp_sal)
1203	END DO
1204	END DO
1205	END DO
1206	!
1207	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point) (div <0)
1208	DO ji = mi0(139), mi1(139)
1209	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
1210	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - hdiv_139_102(jk) * tsn(ji,jj-1,jk,jp_tem) ! middle Atlantic recirculation
1211	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - hdiv_139_102(jk) * tsn(ji,jj-1,jk,jp_sal)
1212	END DO
1213	! ! upper & bottom Atl. reciculation (Atl 101 -> Atl 102) - (div <0)
1214	! ! deep Med flow (Med 102 -> Atl 102) - (div <0)
1215	tsa(ji,jj,22,jp_tem) = tsa(ji,jj,22,jp_tem) + hdiv_141_102(21) * tsn(ji+2,jj,21,jp_tem) & ! deep Med flow
1216	& + hdiv_139_101(21) * tsn(ji,jj-1,21,jp_tem) & ! upper Atlantic recirculation
1217	& + hdiv_139_101(22) * tsn(ji,jj-1,22,jp_tem) ! bottom Atlantic recirculation
1218	tsa(ji,jj,22,jp_sal) = tsa(ji,jj,22,jp_sal) + hdiv_141_102(21) * tsn(ji+2,jj,21,jp_sal) &
1219	& + hdiv_139_101(21) * tsn(ji,jj-1,21,jp_sal) &
1220	& + hdiv_139_101(22) * tsn(ji,jj-1,22,jp_sal)
1221	END DO
1222	END DO
1223	DO jj = mj0(102), mj1(102) !* 141,102 (Med side) (div <0)
1224	DO ji = mi0(141), mi1(141)
1225	DO jk = 1, 14 ! surface flow from Atlantic to Med sea
1226	tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - hdiv_141_102_kt(jk) * tsn(ji-2,jj-1,jk,jp_tem)
1227	tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - hdiv_141_102_kt(jk) * tsn(ji-2,jj-1,jk,jp_sal)
1228	END DO
1229	! ! deeper flow from Med sea to Atlantic
1230	tsa(ji,jj,21,jp_tem) = tsa(ji,jj,21,jp_tem) - hdiv_141_102(21) * tsn(ji,jj,21,jp_tem)
1231	tsa(ji,jj,21,jp_sal) = tsa(ji,jj,21,jp_sal) - hdiv_141_102(21) * tsn(ji,jj,21,jp_sal)
1232	END DO
1233	END DO
1234	!=============
1235	! Adjoint part
1236	!=============
1237	!
1238	DO jj = mj1(102), mj0(102), -1 !* 141,102 (Med side) (div <0)
1239	DO ji = mi1(141), mi0(141), -1
1240	! ! deeper flow from Med sea to Atlantic
1241	hdiv_141_102_ad(21) = hdiv_141_102_ad(21) - tsa_ad(ji,jj,21,jp_sal) * tsn(ji,jj,21,jp_sal)
1242	tsn_ad(ji,jj,21,jp_sal) = tsn_ad(ji,jj,21,jp_sal) - tsa_ad(ji,jj,21,jp_sal) * hdiv_141_102(21)
1243	hdiv_141_102_ad(21) = hdiv_141_102_ad(21) - tsa_ad(ji,jj,21,jp_tem) * tsn(ji,jj,21,jp_tem)
1244	tsn_ad(ji,jj,21,jp_tem) = tsn_ad(ji,jj,21,jp_tem) - tsa_ad(ji,jj,21,jp_tem) * hdiv_141_102(21)
1245	DO jk = 14, 1, -1 ! surface flow from Atlantic to Med sea
1246	hdiv_141_102_kt_ad(jk) = hdiv_141_102_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji-2,jj-1,jk,jp_sal)
1247	tsn_ad(ji-2,jj-1,jk,jp_sal) = tsn_ad(ji-2,jj-1,jk,jp_sal) - tsa_ad(ji,jj,jk,jp_sal) * hdiv_141_102_kt(jk)
1248	hdiv_141_102_kt_ad(jk) = hdiv_141_102_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji-2,jj-1,jk,jp_tem)
1249	tsn_ad(ji-2,jj-1,jk,jp_tem) = tsn_ad(ji-2,jj-1,jk,jp_tem) - tsa_ad(ji,jj,jk,jp_tem) * hdiv_141_102_kt(jk)
1250	END DO
1251	END DO
1252	END DO
1253	!
1254	DO jj = mj1(102), mj0(102), -1 !** 139,102 (Atlantic side, north point) (div <0)
1255	DO ji = mi1(139), mi0(139), -1
1256	! ! upper & bottom Atl. reciculation (Atl 101 -> Atl 102) - (div <0)
1257	! ! deep Med flow (Med 102 -> Atl 102) - (div <0)
1258	hdiv_141_102_ad(21) = hdiv_141_102_ad(21) + tsa_ad(ji,jj,22,jp_sal) * tsn(ji+2,jj,21,jp_sal)
1259	hdiv_139_101_ad(21) = hdiv_139_101_ad(21) + tsa_ad(ji,jj,22,jp_sal) * tsn(ji,jj-1,21,jp_sal)
1260	hdiv_139_101_ad(22) = hdiv_139_101_ad(22) + tsa_ad(ji,jj,22,jp_sal) * tsn(ji,jj-1,22,jp_sal)
1261	tsn_ad(ji+2,jj,21,jp_sal) = tsn_ad(ji+2,jj,21,jp_sal) + tsa_ad(ji,jj,22,jp_sal) * hdiv_141_102(21)
1262	tsn_ad(ji,jj-1,21,jp_sal) = tsn_ad(ji,jj-1,21,jp_sal) + tsa_ad(ji,jj,22,jp_sal) * hdiv_139_101(21)
1263	tsn_ad(ji,jj-1,22,jp_sal) = tsn_ad(ji,jj-1,22,jp_sal) + tsa_ad(ji,jj,22,jp_sal) * hdiv_139_101(22)
1264	hdiv_141_102_ad(21) = hdiv_141_102_ad(21) + tsa_ad(ji,jj,22,jp_tem) * tsn(ji+2,jj,21,jp_tem)
1265	hdiv_139_101_ad(21) = hdiv_139_101_ad(21) + tsa_ad(ji,jj,22,jp_tem) * tsn(ji,jj-1,21,jp_tem)
1266	hdiv_139_101_ad(22) = hdiv_139_101_ad(22) + tsa_ad(ji,jj,22,jp_tem) * tsn(ji,jj-1,22,jp_tem)
1267	tsn_ad(ji+2,jj,21,jp_tem) = tsn_ad(ji+2,jj,21,jp_tem) + tsa_ad(ji,jj,22,jp_tem) * hdiv_141_102(21)
1268	tsn_ad(ji,jj-1,21,jp_tem) = tsn_ad(ji,jj-1,21,jp_tem) + tsa_ad(ji,jj,22,jp_tem) * hdiv_139_101(21)
1269	tsn_ad(ji,jj-1,22,jp_tem) = tsn_ad(ji,jj-1,22,jp_tem) + tsa_ad(ji,jj,22,jp_tem) * hdiv_139_101(22)
1270	DO jk = 20, 15, -1 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
1271	hdiv_139_102_ad(jk) = hdiv_139_102_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji,jj-1,jk,jp_sal)
1272	tsn_ad(ji,jj-1,jk,jp_sal) = tsn_ad(ji,jj-1,jk,jp_sal) - tsa_ad(ji,jj,jk,jp_sal) * hdiv_139_102(jk)
1273	hdiv_139_102_ad(jk) = hdiv_139_102_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji,jj-1,jk,jp_tem)
1274	tsn_ad(ji,jj-1,jk,jp_tem) = tsn_ad(ji,jj-1,jk,jp_tem) - tsa_ad(ji,jj,jk,jp_tem) * hdiv_139_102(jk)
1275	END DO
1276	END DO
1277	END DO
1278	DO jj = mj1(101), mj0(101), -1 !** 139,101 (Atlantic side, south point) (div >0)
1279	DO ji = mi1(139), mi0(139), -1
1280	DO jk = jpkm1, 1, -1 ! surf inflow + mid. & bottom reciculation (from Atlantic)
1281	hdiv_139_101_kt_ad(jk) = hdiv_139_101_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji,jj,jk,jp_sal)
1282	tsn_ad(ji,jj,jk,jp_sal) = tsn_ad(ji,jj,jk,jp_sal) - tsa_ad(ji,jj,jk,jp_sal) * hdiv_139_101_kt(jk)
1283	hdiv_139_101_kt_ad(jk) = hdiv_139_101_kt_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji,jj,jk,jp_tem)
1284	tsn_ad(ji,jj,jk,jp_tem) = tsn_ad(ji,jj,jk,jp_tem) - tsa_ad(ji,jj,jk,jp_tem) * hdiv_139_101_kt(jk)
1285	END DO
1286	END DO
1287	END DO
1288	! ! ---------------- !
1289	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
1290	! ! --------=======- !
1291	! at this stage, (ua,va) are the after velocity, not the tendancy
1292	! compute the velocity from the divergence at T-point
1293	!
1294	DO jj = mj0(102), mj1(102) !** 140,102 (Med side) (140 not 141 as it is a U-point)
1295	DO ji = mi0(140), mi1(140) ! div >0 => ua <0, opposite sign
1296	hdiv_141_102_ad(:) = hdiv_141_102_ad(:) - ua_ad(ji,jj,:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
1297	& * e2u(ji,jj) * fse3u(ji,jj,:)
1298	END DO
1299	END DO
1300	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point)
1301	DO ji = mi0(139), mi1(139) ! div <0 => ua <0, same sign
1302	hdiv_139_102_ad(:) = hdiv_139_102_ad(:) + ua_ad(ji,jj,:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
1303	& * e2u(ji,jj) * fse3u(ji,jj,:)
1304	END DO
1305	END DO
1306	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point)
1307	DO ji = mi0(139), mi1(139) ! div >0 => ua >0, same sign
1308	hdiv_139_101_kt_ad(:) = hdiv_139_101_kt_ad(:) + ua_ad(ji,jj,:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
1309	& * e2u(ji,jj) * fse3u(ji,jj,:)
1310	END DO
1311	END DO
1312	!
1313	END SELECT
1314	!
1315	END SUBROUTINE cla_gibraltar_adj
1316
1317	SUBROUTINE cla_hormuz_tan( cd_td )
1318	!! -------------------------------------------------------------------
1319	!! * ROUTINE div_hormuz_tan *
1320	!!
1321	!! ** Purpose : update the now horizontal divergence, the tracer
1322	!! tendancyand the after velocity in vicinity of Hormuz
1323	!! strait ( Persian Gulf - Indian ocean ).
1324	!!
1325	!! ** Method : Hormuz strait
1326	!! ______________
1327	!! \|/////\|<== surface inflow
1328	!! 94 \|/////\|
1329	!! \|/////\|==> deep outflow
1330	!! \|_____\|_______
1331	!! 171 172
1332	!!---------------------------------------------------------------------
1333	CHARACTER(len=1), INTENT(in) :: cd_td ! ='ini' initialisation
1334	!! ! ='div' update the divergence
1335	!! ! ='tra' update the tracers
1336	!! ! ='spg' update after velocity
1337	!!
1338	INTEGER :: ji, jj, jk ! dummy loop indices
1339	REAL(wp) :: zio_flow ! temporary scalar
1340	!!---------------------------------------------------------------------
1341	!
1342	SELECT CASE( cd_td )
1343	! ! ---------------- !
1344	CASE( 'ini' ) ! initialisation !
1345	! ! ---------------- !
1346	! !** profile of horizontal divergence due to cross-land advection
1347	zio_flow = 1.e6 ! imposed in/out flow
1348	!
1349	hdiv_172_94(:) = 0.e0
1350	hdiv_172_94_tl(:) = 0.e0
1351	!
1352	DO jj = mj0(94), mj1(94) ! in/out flow at (i,j) = (172,94)
1353	DO ji = mi0(172), mi1(172)
1354	DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
1355	hdiv_172_94(jk) = - ( zio_flow / 8.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1356	END DO
1357	DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
1358	hdiv_172_94(jk) = + ( zio_flow / 3.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1359	END DO
1360	END DO
1361	END DO
1362	! !** T & S profile in the Hormuz strait (use in deep outflow)
1363	! Temperature and Salinity
1364	t_171_94_hor(:) = 0.e0 ; s_171_94_hor(:) = 0.e0
1365	t_171_94_hor(16) = 18.4 ; s_171_94_hor(16) = 36.27
1366	t_171_94_hor(17) = 17.8 ; s_171_94_hor(17) = 36.4
1367	t_171_94_hor(18) = 16. ; s_171_94_hor(18) = 36.27
1368	!
1369	! ! ---------------- !
1370	CASE( 'div' ) ! update hdivn ! (call by divcur module)
1371	! ! ---------=====-- !
1372	!
1373	DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side)
1374	DO ji = mi0(172), mi1(172)
1375	hdivn_tl(ji,jj,:) = hdivn_tl(ji,jj,:) + hdiv_172_94_tl(:)
1376	END DO
1377	END DO
1378	! ! ---------------- !
1379	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
1380	! ! --------=======- !
1381	!
1382	DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side)
1383	DO ji = mi0(172), mi1(172)
1384	DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
1385	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) &
1386	& - hdiv_172_94_tl(jk) * tsn(ji,jj,jk,jp_tem) &
1387	& - hdiv_172_94(jk) * tsn_tl(ji,jj,jk,jp_tem)
1388	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) &
1389	& - hdiv_172_94_tl(jk) * tsn(ji,jj,jk,jp_sal) &
1390	& - hdiv_172_94(jk) * tsn_tl(ji,jj,jk,jp_sal)
1391	END DO
1392	DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
1393	tsa_tl(ji,jj,jk,jp_tem) = tsa_tl(ji,jj,jk,jp_tem) - hdiv_172_94_tl(jk) * t_171_94_hor(jk)
1394	tsa_tl(ji,jj,jk,jp_sal) = tsa_tl(ji,jj,jk,jp_sal) - hdiv_172_94_tl(jk) * s_171_94_hor(jk)
1395	END DO
1396	END DO
1397	END DO
1398	! ! ---------------- !
1399	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
1400	! ! --------=======- !
1401	! No barotropic flow through Hormuz strait
1402	! at this stage, (ua,va) are the after velocity, not the tendancy
1403	! compute the velocity from the divergence at T-point
1404	DO jj = mj0(94), mj1(94) !** 171,94 (Indian ocean side) (171 not 172 as it is the western U-point)
1405	DO ji = mi0(171), mi1(171) ! div >0 => ua >0, opposite sign
1406	ua_tl(ji,jj,:) = - hdiv_172_94_tl(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
1407	& * e2u(ji,jj) * fse3u(ji,jj,:)
1408	END DO
1409	END DO
1410	!
1411	END SELECT
1412	!
1413	END SUBROUTINE cla_hormuz_tan
1414
1415
1416	SUBROUTINE cla_hormuz_adj( cd_td )
1417	!! -------------------------------------------------------------------
1418	!! * ROUTINE div_hormuz_adj *
1419	!!
1420	!! ** Purpose : update the now horizontal divergence, the tracer
1421	!! tendancyand the after velocity in vicinity of Hormuz
1422	!! strait ( Persian Gulf - Indian ocean ).
1423	!!
1424	!! ** Method : Hormuz strait
1425	!! ______________
1426	!! \|/////\|<== surface inflow
1427	!! 94 \|/////\|
1428	!! \|/////\|==> deep outflow
1429	!! \|_____\|_______
1430	!! 171 172
1431	!!---------------------------------------------------------------------
1432	CHARACTER(len=1), INTENT(in) :: cd_td ! ='ini' initialisation
1433	!! ! ='div' update the divergence
1434	!! ! ='tra' update the tracers
1435	!! ! ='spg' update after velocity
1436	!!
1437	INTEGER :: ji, jj, jk ! dummy loop indices
1438	REAL(wp) :: zio_flow ! temporary scalar
1439	!!---------------------------------------------------------------------
1440	!
1441	SELECT CASE( cd_td )
1442	! ! ---------------- !
1443	CASE( 'ini' ) ! initialisation !
1444	! ! ---------------- !
1445	! !** profile of horizontal divergence due to cross-land advection
1446	zio_flow = 1.e6 ! imposed in/out flow
1447	!
1448	hdiv_172_94(:) = 0.e0
1449	hdiv_172_94_ad(:) = 0.e0
1450	!
1451	DO jj = mj0(94), mj1(94) ! in/out flow at (i,j) = (172,94)
1452	DO ji = mi0(172), mi1(172)
1453	DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
1454	hdiv_172_94(jk) = - ( zio_flow / 8.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1455	END DO
1456	DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
1457	hdiv_172_94(jk) = + ( zio_flow / 3.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
1458	END DO
1459	END DO
1460	END DO
1461	! !** T & S profile in the Hormuz strait (use in deep outflow)
1462	! Temperature and Salinity
1463	t_171_94_hor(:) = 0.e0 ; s_171_94_hor(:) = 0.e0
1464	t_171_94_hor(16) = 18.4 ; s_171_94_hor(16) = 36.27
1465	t_171_94_hor(17) = 17.8 ; s_171_94_hor(17) = 36.4
1466	t_171_94_hor(18) = 16. ; s_171_94_hor(18) = 36.27
1467	!
1468	! ! ---------------- !
1469	CASE( 'div' ) ! update hdivn ! (call by divcur module)
1470	! ! ---------=====-- !
1471	!
1472	DO jj = mj1(94), mj0(94), -1 !** 172,94 (Indian ocean side)
1473	DO ji = mi1(172), mi0(172), -1
1474	hdiv_172_94_ad(:) = hdiv_172_94_ad(:) + hdivn_ad(ji,jj,:)
1475	END DO
1476	END DO
1477	! ! ---------------- !
1478	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
1479	! ! --------=======- !
1480	!
1481	DO jj = mj1(94), mj0(94), -1 !** 172,94 (Indian ocean side)
1482	DO ji = mi1(172), mi0(172), -1
1483	DO jk = 18, 16, -1 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
1484	hdiv_172_94_ad(jk) = hdiv_172_94_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * s_171_94_hor(jk)
1485	hdiv_172_94_ad(jk) = hdiv_172_94_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * s_171_94_hor(jk)
1486	END DO
1487	DO jk = 8, 1, -1 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
1488	hdiv_172_94_ad(jk) = hdiv_172_94_ad(jk) - tsa_ad(ji,jj,jk,jp_sal) * tsn(ji,jj,jk,jp_sal)
1489	tsn_ad(ji,jj,jk,jp_sal) = tsn_ad(ji,jj,jk,jp_sal) - tsa_ad(ji,jj,jk,jp_sal) * hdiv_172_94(jk)
1490	hdiv_172_94_ad(jk) = hdiv_172_94_ad(jk) - tsa_ad(ji,jj,jk,jp_tem) * tsn(ji,jj,jk,jp_tem)
1491	tsn_ad(ji,jj,jk,jp_tem) = tsn_ad(ji,jj,jk,jp_tem) - tsa_ad(ji,jj,jk,jp_tem) * hdiv_172_94(jk)
1492	END DO
1493	END DO
1494	END DO
1495	! ! ---------------- !
1496	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
1497	! ! --------=======- !
1498	! No barotropic flow through Hormuz strait
1499	! at this stage, (ua,va) are the after velocity, not the tendancy
1500	! compute the velocity from the divergence at T-point
1501	DO jj = mj0(94), mj1(94) !** 171,94 (Indian ocean side) (171 not 172 as it is the western U-point)
1502	DO ji = mi0(171), mi1(171) ! div >0 => ua >0, opposite sign
1503	hdiv_172_94_ad(:) = hdiv_172_94_ad(:) - ua_ad(ji,jj,:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
1504	& * e2u(ji,jj) * fse3u(ji,jj,:)
1505	END DO
1506	END DO
1507	!
1508	END SELECT
1509	!
1510	END SUBROUTINE cla_hormuz_adj
1511
1512	SUBROUTINE cla_div_adj_tst( kumadt )
1513	!!-----------------------------------------------------------------------
1514	!!
1515	!! * ROUTINE cla_divadj_tst *
1516	!!
1517	!! ** Purpose : Test the adjoint routine.
1518	!!
1519	!! ** Method : Verify the scalar product
1520	!!
1521	!! ( L dx )^T W dy = dx^T L^T W dy
1522	!!
1523	!! where L = tangent routine
1524	!! L^T = adjoint routine
1525	!! W = diagonal matrix of scale factors
1526	!! dx = input perturbation (random field)
1527	!! dy = L dx
1528	!!
1529	!!
1530	!! History :
1531	!! ! 08-08 (A. Vidard)
1532	!!-----------------------------------------------------------------------
1533	!! * Modules used
1534
1535	!! * Arguments
1536	INTEGER, INTENT(IN) :: &
1537	& kumadt ! Output unit
1538
1539	!! * Local declarations
1540	INTEGER :: &
1541	& ji, & ! dummy loop indices
1542	& jj, &
1543	& jk, &
1544	& jt
1545
1546	REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
1547	& zemp_tlin, & ! Tangent input
1548	& zemp_tlout, & ! Tangent output
1549	& zemp_adin, & ! adjoint input
1550	& zemp_adout, & ! adjoint output
1551	& zemp
1552
1553	REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
1554	& zhdivn_tlin, & ! Tangent input
1555	& zhdivn_tlout, & ! Tangent output
1556	& zhdivn_adin, & ! adjoint input
1557	& zhdivn_adout, & ! adjoint output
1558	& zhdivn
1559
1560	REAL(KIND=wp) :: &
1561	& zsp1, & ! scalar product involving the tangent routine
1562	& zsp1_1, & ! scalar product components
1563	& zsp1_2, &
1564	& zsp2, & ! scalar product involving the adjoint routine
1565	& zsp2_1, & ! scalar product components
1566	& zsp2_2, &
1567	& zsp2_3, &
1568	& zsp2_4, &
1569	& zsp2_5
1570
1571	CHARACTER(LEN=14) :: cl_name
1572
1573	ALLOCATE( &
1574	& zhdivn_tlin(jpi,jpj,jpk), &
1575	& zhdivn_tlout(jpi,jpj,jpk), &
1576	& zhdivn_adin(jpi,jpj,jpk), &
1577	& zhdivn_adout(jpi,jpj,jpk), &
1578	& zemp_tlin(jpi,jpj), &
1579	& zemp_tlout(jpi,jpj), &
1580	& zemp_adin(jpi,jpj), &
1581	& zemp_adout(jpi,jpj), &
1582	& zhdivn(jpi,jpj,jpk), &
1583	& zemp(jpi,jpj) )
1584
1585	DO jt = 1, 3
1586	!==================================================================
1587	! 1) dx = ( un_tl, vn_tl, hdivn_tl ) and
1588	! dy = ( hdivb_tl, hdivn_tl )
1589	!==================================================================
1590
1591	!--------------------------------------------------------------------
1592	! Reset the tangent and adjoint variables
1593	!--------------------------------------------------------------------
1594
1595	zhdivn_tlin(:,:,:) = 0._wp
1596	zhdivn_tlout(:,:,:) = 0._wp
1597	zhdivn_adin(:,:,:) = 0._wp
1598	zhdivn_adout(:,:,:) = 0._wp
1599	zemp_tlin(:,:) = 0._wp
1600	zemp_tlout(:,:) = 0._wp
1601	zemp_adin(:,:) = 0._wp
1602	zemp_adout(:,:) = 0._wp
1603	zhdivn(:,:,:) = 0._wp
1604	zemp(:,:) = 0._wp
1605
1606	hdivn_tl(:,:,:) = 0._wp
1607	emp_tl(:,:) = 0._wp
1608	hdivn_ad(:,:,:) = 0._wp
1609	emp_ad(:,:) = 0._wp
1610
1611	CALL grid_random( zemp, 'T', 0.0_wp, stdssh )
1612	CALL grid_random( zhdivn, 'T', 0.0_wp, stdu )
1613
1614	DO jj = nldj, nlej
1615	DO ji = nldi, nlei
1616	zemp_tlin(ji,jj) = zemp(ji,jj)
1617	END DO
1618	END DO
1619	DO jk = 1, jpk
1620	DO jj = nldj, nlej
1621	DO ji = nldi, nlei
1622	zhdivn_tlin(ji,jj,jk) = zhdivn(ji,jj,jk)
1623	END DO
1624	END DO
1625	END DO
1626	!--------------------------------------------------------------------
1627	! Call the tangent routine: dy = L dx
1628	!--------------------------------------------------------------------
1629	emp_tl(:,:) = zemp_tlin(:,:)
1630	hdivn_tl(:,:,:) = zhdivn_tlin(:,:,:)
1631
1632	SELECT CASE (jt)
1633	CASE(1)
1634	nbab = 1
1635	ngib = 0
1636	nhor = 0
1637	CALL cla_div_tan( nit000 )
1638	CASE(2)
1639	nbab = 0
1640	ngib = 1
1641	nhor = 0
1642	CALL cla_div_tan( nit000 )
1643	CASE(3)
1644	nbab = 0
1645	ngib = 0
1646	nhor = 1
1647	CALL cla_div_tan( nit000 )
1648	END SELECT
1649
1650	zhdivn_tlout(:,:,:) = hdivn_tl(:,:,:)
1651
1652
1653	!--------------------------------------------------------------------
1654	! Initialize the adjoint variables: dy^* = W dy
1655	!--------------------------------------------------------------------
1656	DO jk = 1, jpk
1657	DO jj = nldj, nlej
1658	DO ji = nldi, nlei
1659	zhdivn_adin(ji,jj,jk) = zhdivn_tlout(ji,jj,jk) &
1660	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
1661	& * tmask(ji,jj,jk)
1662	END DO
1663	END DO
1664	END DO
1665
1666	!--------------------------------------------------------------------
1667	! Compute the scalar product: ( L dx )^T W dy
1668	!--------------------------------------------------------------------
1669
1670	zsp1 = DOT_PRODUCT( zhdivn_tlout, zhdivn_adin )
1671
1672	!--------------------------------------------------------------------
1673	! Call the adjoint routine: dx^* = L^T dy^*
1674	!--------------------------------------------------------------------
1675
1676	hdivn_ad(:,:,:) = zhdivn_adin(:,:,:)
1677
1678	SELECT CASE (jt)
1679	CASE(1)
1680	nbab = 1
1681	ngib = 0
1682	nhor = 0
1683	CALL cla_div_adj( nit000 )
1684	CASE(2)
1685	nbab = 0
1686	ngib = 1
1687	nhor = 0
1688	CALL cla_div_adj( nit000 )
1689	CASE(3)
1690	nbab = 0
1691	ngib = 0
1692	nhor = 1
1693	CALL cla_div_adj( nit000 )
1694	END SELECT
1695
1696	zemp_adout (:,:) = emp_ad (:,:)
1697	zhdivn_adout(:,:,:) = hdivn_ad(:,:,:)
1698
1699	!--------------------------------------------------------------------
1700	! Compute the scalar product: dx^T L^T W dy
1701	!--------------------------------------------------------------------
1702
1703	zsp2_1 = DOT_PRODUCT( zhdivn_tlin, zhdivn_adout )
1704	zsp2_2 = DOT_PRODUCT( zemp_tlin, zemp_adout )
1705	zsp2 = zsp2_1 + zsp2_2
1706
1707	SELECT CASE (jt)
1708	CASE(1)
1709	cl_name = 'cladivadj babm'
1710	CASE(2)
1711	cl_name = 'cladivadj gibr'
1712	CASE(3)
1713	cl_name = 'cladivadj horm'
1714	END SELECT
1715	CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
1716	END DO
1717
1718	DEALLOCATE( &
1719	& zhdivn_tlin, &
1720	& zhdivn_tlout, &
1721	& zhdivn_adin, &
1722	& zhdivn_adout, &
1723	& zemp_tlin, &
1724	& zemp_tlout, &
1725	& zemp_adin, &
1726	& zemp_adout, &
1727	& zhdivn, &
1728	& zemp )
1729
1730	END SUBROUTINE cla_div_adj_tst
1731
1732	SUBROUTINE cla_traadv_adj_tst( kumadt )
1733	!!-----------------------------------------------------------------------
1734	!!
1735	!! * ROUTINE cla_divadj_tst *
1736	!!
1737	!! ** Purpose : Test the adjoint routine.
1738	!!
1739	!! ** Method : Verify the scalar product
1740	!!
1741	!! ( L dx )^T W dy = dx^T L^T W dy
1742	!!
1743	!! where L = tangent routine
1744	!! L^T = adjoint routine
1745	!! W = diagonal matrix of scale factors
1746	!! dx = input perturbation (random field)
1747	!! dy = L dx
1748	!!
1749	!!
1750	!! History :
1751	!! ! 08-08 (A. Vidard)
1752	!!-----------------------------------------------------------------------
1753	!! * Modules used
1754
1755	!! * Arguments
1756	INTEGER, INTENT(IN) :: &
1757	& kumadt ! Output unit
1758
1759	!! * Local declarations
1760	INTEGER :: &
1761	& ji, & ! dummy loop indices
1762	& jj, &
1763	& jk, &
1764	& jt
1765
1766	REAL(KIND=wp) :: &
1767	& zsp1, & ! scalar product involving the tangent routine
1768	& zsp2 ! scalar product involving the adjoint routine
1769	REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
1770	& ztn_tlin , & ! Tangent input
1771	& zsn_tlin , zta_tlin , zsa_tlin , & ! Tangent input
1772	& ztn_adout, & ! Adjoint output
1773	& zsn_adout, zta_adout, zsa_adout, & ! Adjoint output
1774	& zta_tlout, zsa_tlout, & ! Tangent output
1775	& zta_adin , zsa_adin , & ! Adjoint input
1776	& zr ! 3D random field
1777	CHARACTER(LEN=14) ::&
1778	& cl_name
1779	! Allocate memory
1780
1781	ALLOCATE( &
1782	& ztn_tlin( jpi,jpj,jpk), zsn_tlin( jpi,jpj,jpk), zta_tlin( jpi,jpj,jpk), &
1783	& zsa_tlin( jpi,jpj,jpk), zta_tlout(jpi,jpj,jpk), zsa_tlout(jpi,jpj,jpk), &
1784	& zta_adin( jpi,jpj,jpk), zsa_adin( jpi,jpj,jpk), &
1785	& ztn_adout(jpi,jpj,jpk), &
1786	& zsn_adout(jpi,jpj,jpk), zta_adout(jpi,jpj,jpk), zsa_adout(jpi,jpj,jpk), &
1787	& zr( jpi,jpj,jpk) &
1788	& )
1789
1790	DO jt = 1, 3
1791	!==================================================================
1792	! 1) dx = ( un_tl, vn_tl, hdivn_tl ) and
1793	! dy = ( hdivb_tl, hdivn_tl )
1794	!==================================================================
1795
1796	!--------------------------------------------------------------------
1797	! Reset the tangent and adjoint variables
1798	!--------------------------------------------------------------------
1799	ztn_tlin( :,:,:) = 0.0_wp
1800	zsn_tlin( :,:,:) = 0.0_wp
1801	zta_tlin( :,:,:) = 0.0_wp
1802	zsa_tlin( :,:,:) = 0.0_wp
1803	zta_tlout(:,:,:) = 0.0_wp
1804	zsa_tlout(:,:,:) = 0.0_wp
1805	zta_adin( :,:,:) = 0.0_wp
1806	zsa_adin( :,:,:) = 0.0_wp
1807	ztn_adout(:,:,:) = 0.0_wp
1808	zsn_adout(:,:,:) = 0.0_wp
1809	zta_adout(:,:,:) = 0.0_wp
1810	zsa_adout(:,:,:) = 0.0_wp
1811	zr( :,:,:) = 0.0_wp
1812
1813	tsn_ad(:,:,:,jp_tem) = 0.0_wp
1814	tsn_ad(:,:,:,jp_sal) = 0.0_wp
1815
1816	CALL grid_random( zr, 'T', 0.0_wp, stdt )
1817	DO jk = 1, jpk
1818	DO jj = nldj, nlej
1819	DO ji = nldi, nlei
1820	ztn_tlin(ji,jj,jk) = zr(ji,jj,jk)
1821	END DO
1822	END DO
1823	END DO
1824	CALL grid_random( zr, 'T', 0.0_wp, stds )
1825	DO jk = 1, jpk
1826	DO jj = nldj, nlej
1827	DO ji = nldi, nlei
1828	zsn_tlin(ji,jj,jk) = zr(ji,jj,jk)
1829	END DO
1830	END DO
1831	END DO
1832	CALL grid_random( zr, 'T', 0.0_wp, stdt )
1833	DO jk = 1, jpk
1834	DO jj = nldj, nlej
1835	DO ji = nldi, nlei
1836	zta_tlin(ji,jj,jk) = zr(ji,jj,jk)
1837	END DO
1838	END DO
1839	END DO
1840	CALL grid_random( zr, 'T', 0.0_wp, stds )
1841	DO jk = 1, jpk
1842	DO jj = nldj, nlej
1843	DO ji = nldi, nlei
1844	zsa_tlin(ji,jj,jk) = zr(ji,jj,jk)
1845	END DO
1846	END DO
1847	END DO
1848
1849	tsn_tl(:,:,:,jp_tem) = ztn_tlin(:,:,:)
1850	tsn_tl(:,:,:,jp_sal) = zsn_tlin(:,:,:)
1851	tsa_tl(:,:,:,jp_tem) = zta_tlin(:,:,:)
1852	tsa_tl(:,:,:,jp_sal) = zsa_tlin(:,:,:)
1853
1854	!--------------------------------------------------------------------
1855	! Call the tangent routine: dy = L dx
1856	!--------------------------------------------------------------------
1857
1858	SELECT CASE (jt)
1859	CASE(1)
1860	nbab = 1
1861	ngib = 0
1862	nhor = 0
1863	CALL cla_traadv_tan( nit000 )
1864	CASE(2)
1865	nbab = 0
1866	ngib = 1
1867	nhor = 0
1868	CALL cla_traadv_tan( nit000 )
1869	CASE(3)
1870	nbab = 0
1871	ngib = 0
1872	nhor = 1
1873	CALL cla_traadv_tan( nit000 )
1874	END SELECT
1875
1876	zta_tlout(:,:,:) = tsa_tl(:,:,:,jp_tem)
1877	zsa_tlout(:,:,:) = tsa_tl(:,:,:,jp_sal)
1878
1879	!--------------------------------------------------------------------
1880	! Initialize the adjoint variables: dy^* = W dy
1881	!--------------------------------------------------------------------
1882	DO jk = 1, jpk
1883	DO jj = nldj, nlej
1884	DO ji = nldi, nlei
1885	zta_adin(ji,jj,jk) = zta_tlout(ji,jj,jk) &
1886	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
1887	& * tmask(ji,jj,jk) * wesp_t(jk)
1888	zsa_adin(ji,jj,jk) = zsa_tlout(ji,jj,jk) &
1889	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
1890	& * tmask(ji,jj,jk) * wesp_s(jk)
1891	END DO
1892	END DO
1893	END DO
1894	!--------------------------------------------------------------------
1895	! Compute the scalar product: ( L dx )^T W dy
1896	!--------------------------------------------------------------------
1897
1898	zsp1 = DOT_PRODUCT( zta_tlout, zta_adin ) &
1899	& + DOT_PRODUCT( zsa_tlout, zsa_adin )
1900
1901	!--------------------------------------------------------------------
1902	! Call the adjoint routine: dx^* = L^T dy^*
1903	!--------------------------------------------------------------------
1904	tsa_ad(:,:,:,jp_tem) = zta_adin(:,:,:)
1905	tsa_ad(:,:,:,jp_sal) = zsa_adin(:,:,:)
1906
1907	SELECT CASE (jt)
1908	CASE(1)
1909	nbab = 1
1910	ngib = 0
1911	nhor = 0
1912	CALL cla_traadv_adj( nit000 )
1913	CASE(2)
1914	nbab = 0
1915	ngib = 1
1916	nhor = 0
1917	CALL cla_traadv_adj( nit000 )
1918	CASE(3)
1919	nbab = 0
1920	ngib = 0
1921	nhor = 1
1922	CALL cla_traadv_adj( nit000 )
1923	END SELECT
1924
1925	ztn_adout(:,:,:) = tsn_ad(:,:,:,jp_tem)
1926	zsn_adout(:,:,:) = tsn_ad(:,:,:,jp_sal)
1927	zta_adout(:,:,:) = tsa_ad(:,:,:,jp_tem)
1928	zsa_adout(:,:,:) = tsa_ad(:,:,:,jp_sal)
1929
1930	zsp2 = DOT_PRODUCT( ztn_tlin, ztn_adout ) &
1931	& + DOT_PRODUCT( zsn_tlin, zsn_adout ) &
1932	& + DOT_PRODUCT( zta_tlin, zta_adout ) &
1933	& + DOT_PRODUCT( zsa_tlin, zsa_adout )
1934
1935	SELECT CASE (jt)
1936	CASE(1)
1937	cl_name = 'clatraadv babm'
1938	CASE(2)
1939	cl_name = 'clatraadv gibr'
1940	CASE(3)
1941	cl_name = 'clatraadv horm'
1942	END SELECT
1943	CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
1944	END DO
1945
1946	DEALLOCATE( &
1947	& ztn_tlin , zsn_tlin , &
1948	& zta_tlin , zsa_tlin , zta_tlout, zsa_tlout, zta_adin , &
1949	& zsa_adin , ztn_adout, &
1950	& zsn_adout, zta_adout, zsa_adout, zr &
1951	& )
1952
1953	END SUBROUTINE cla_traadv_adj_tst
1954
1955	SUBROUTINE cla_dynspg_adj_tst( kumadt )
1956	!!-----------------------------------------------------------------------
1957	!!
1958	!! * ROUTINE cla_divadj_tst *
1959	!!
1960	!! ** Purpose : Test the adjoint routine.
1961	!!
1962	!! ** Method : Verify the scalar product
1963	!!
1964	!! ( L dx )^T W dy = dx^T L^T W dy
1965	!!
1966	!! where L = tangent routine
1967	!! L^T = adjoint routine
1968	!! W = diagonal matrix of scale factors
1969	!! dx = input perturbation (random field)
1970	!! dy = L dx
1971	!!
1972	!!
1973	!! History :
1974	!! ! 08-08 (A. Vidard)
1975	!!-----------------------------------------------------------------------
1976	!! * Modules used
1977
1978	!! * Arguments
1979	INTEGER, INTENT(IN) :: &
1980	& kumadt ! Output unit
1981
1982	!! * Local declarations
1983	INTEGER :: &
1984	& ji, & ! dummy loop indices
1985	& jj, &
1986	& jk, &
1987	& jt
1988
1989	REAL(KIND=wp) :: &
1990	& zsp1, & ! scalar product involving the tangent routine
1991	& zsp2 ! scalar product involving the adjoint routine
1992	REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: &
1993	& zua_tlin , & ! Tangent input
1994	& zva_tlin , & ! Tangent input
1995	& zua_tlout, & ! Tangent output
1996	& zva_tlout, & ! Tangent output
1997	& zua_adin , & ! Adjoint input
1998	& zva_adin , & ! Adjoint input
1999	& zua_adout, & ! Adjoint output
2000	& zva_adout, & ! Adjoint output
2001	& zr3d ! 3D random field
2002	CHARACTER(LEN=14) :: &
2003	& cl_name
2004	! Allocate memory
2005
2006	ALLOCATE( &
2007	& zua_tlin( jpi,jpj,jpk), &
2008	& zva_tlin( jpi,jpj,jpk), &
2009	& zua_tlout( jpi,jpj,jpk), &
2010	& zva_tlout( jpi,jpj,jpk), &
2011	& zua_adin( jpi,jpj,jpk), &
2012	& zva_adin( jpi,jpj,jpk), &
2013	& zua_adout( jpi,jpj,jpk), &
2014	& zva_adout( jpi,jpj,jpk), &
2015	& zr3d( jpi,jpj,jpk) )
2016
2017	DO jt = 1, 3
2018	!==================================================================
2019	! 1) dx = ( un_tl, vn_tl, hdivn_tl ) and
2020	! dy = ( hdivb_tl, hdivn_tl )
2021	!==================================================================
2022
2023	!--------------------------------------------------------------------
2024	! Reset the tangent and adjoint variables
2025	!--------------------------------------------------------------------
2026	ua_ad( :,:,:) = 0.0_wp
2027	va_ad( :,:,:) = 0.0_wp
2028	!--------------------------------------------------------------------
2029	! Initialize the tangent input with random noise: dx
2030	!--------------------------------------------------------------------
2031
2032	CALL grid_random( zr3d, 'U', 0.0_wp, stdu )
2033	zua_tlin(:,:,:) = zr3d(:,:,:)
2034	CALL grid_random( zr3d, 'V', 0.0_wp, stdv )
2035	zva_tlin(:,:,:) = zr3d(:,:,:)
2036
2037	ua_tl = zua_tlin
2038	va_tl = zva_tlin
2039	!--------------------------------------------------------------------
2040	! Call the tangent routine: dy = L dx
2041	!--------------------------------------------------------------------
2042
2043	SELECT CASE (jt)
2044	CASE(1)
2045	nbab = 1
2046	ngib = 0
2047	nhor = 0
2048	CALL cla_traadv_tan( nit000 )
2049	CASE(2)
2050	nbab = 0
2051	ngib = 1
2052	nhor = 0
2053	CALL cla_traadv_tan( nit000 )
2054	CASE(3)
2055	nbab = 0
2056	ngib = 0
2057	nhor = 1
2058	CALL cla_traadv_tan( nit000 )
2059	END SELECT
2060
2061	zua_tlout = ua_tl
2062	zva_tlout = va_tl
2063
2064	!--------------------------------------------------------------------
2065	! Initialize the adjoint variables: dy^* = W dy
2066	!--------------------------------------------------------------------
2067	DO jk = 1, jpk
2068	DO jj = nldj, nlej
2069	DO ji = nldi, nlei
2070	zua_adin(ji,jj,jk) = zua_tlout(ji,jj,jk) &
2071	& * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) &
2072	& * umask(ji,jj,jk)
2073	zva_adin(ji,jj,jk) = zva_tlout(ji,jj,jk) &
2074	& * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) &
2075	& * vmask(ji,jj,jk)
2076	END DO
2077	END DO
2078	END DO
2079	!--------------------------------------------------------------------
2080	! Compute the scalar product: ( L dx )^T W dy
2081	!--------------------------------------------------------------------
2082
2083	zsp1 = DOT_PRODUCT( zua_tlout, zua_adin ) &
2084	& + DOT_PRODUCT( zva_tlout, zva_adin )
2085
2086	!--------------------------------------------------------------------
2087	! Call the adjoint routine: dx^* = L^T dy^*
2088	!--------------------------------------------------------------------
2089	ua_ad = zua_adin
2090	va_ad = zva_adin
2091
2092	SELECT CASE (jt)
2093	CASE(1)
2094	nbab = 1
2095	ngib = 0
2096	nhor = 0
2097	CALL cla_div_adj( nit000 )
2098	CASE(2)
2099	nbab = 0
2100	ngib = 1
2101	nhor = 0
2102	CALL cla_div_adj( nit000 )
2103	CASE(3)
2104	nbab = 0
2105	ngib = 0
2106	nhor = 1
2107	CALL cla_div_adj( nit000 )
2108	END SELECT
2109
2110	zua_adout = ua_ad
2111	zva_adout = va_ad
2112
2113	zsp2 = DOT_PRODUCT( zua_tlin , zua_adout ) &
2114	& + DOT_PRODUCT( zva_tlin , zva_adout )
2115
2116	SELECT CASE (jt)
2117	CASE(1)
2118	cl_name = 'cladynspg babm'
2119	CASE(2)
2120	cl_name = 'cladynspg gibr'
2121	CASE(3)
2122	cl_name = 'cladynspg horm'
2123	END SELECT
2124	CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
2125	END DO
2126
2127	DEALLOCATE( &
2128	& zua_tlin, &
2129	& zva_tlin, &
2130	& zua_tlout, &
2131	& zva_tlout, &
2132	& zua_adin, &
2133	& zva_adin, &
2134	& zua_adout, &
2135	& zva_adout, &
2136	& zr3d &
2137	& )
2138
2139	END SUBROUTINE cla_dynspg_adj_tst
2140
2141	#else
2142	!!----------------------------------------------------------------------
2143	!! Default key Dummy module
2144	!!----------------------------------------------------------------------
2145	USE lib_mpp, ONLY: ctl_stop
2146	CONTAINS
2147	SUBROUTINE cla_init_tam
2148	CALL ctl_stop( 'cla_init: Cross Land Advection hard coded for ORCA_R2 with 31 levels' )
2149	END SUBROUTINE cla_init_tam
2150	SUBROUTINE cla_div_tan( kt )
2151	WRITE(,) 'cla_div: You should have not see this print! error?', kt
2152	END SUBROUTINE cla_div_tan
2153	SUBROUTINE cla_traadv_tan( kt )
2154	WRITE(,) 'cla_traadv: You should have not see this print! error?', kt
2155	END SUBROUTINE cla_traadv_tan
2156	SUBROUTINE cla_dynspg_tan( kt )
2157	WRITE(,) 'dyn_spg_cla: You should have not see this print! error?', kt
2158	END SUBROUTINE cla_dynspg_tan
2159	SUBROUTINE cla_div_adj( kt )
2160	WRITE(,) 'cla_div: You should have not see this print! error?', kt
2161	END SUBROUTINE cla_div_adj
2162	SUBROUTINE cla_traadv_adj( kt )
2163	WRITE(,) 'cla_traadv: You should have not see this print! error?', kt
2164	END SUBROUTINE cla_traadv_adj
2165	SUBROUTINE cla_dynspg_adj( kt )
2166	WRITE(,) 'dyn_spg_cla: You should have not see this print! error?', kt
2167	END SUBROUTINE cla_dynspg_adj
2168	SUBROUTINE cla_div_adj_tst( kt )
2169	WRITE(,) 'cla_div: You should have not see this print! error?', kt
2170	END SUBROUTINE cla_div_adj_tst
2171	SUBROUTINE cla_traadv_adj_tst( kt )
2172	WRITE(,) 'cla_traadv: You should have not see this print! error?', kt
2173	END SUBROUTINE cla_traadv_adj_tst
2174	SUBROUTINE cla_dynspg_adj_tst( kt )
2175	WRITE(,) 'dyn_spg_cla: You should have not see this print! error?', kt
2176	END SUBROUTINE cla_dynspg_adj_tst
2177	#endif
2178	#endif
2179	END MODULE cla_tam

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

Download in other formats: