Context Navigation

isotrac_routines_mod.F90 @ 5367

Last change on this file since 5367 was 5158, checked in by abarral, 4 months ago
Add missing klon on strataer_emiss_mod.F90 Correct various missing explicit declarations Replace tabs by spaces (tabs are not part of the fortran charset) Continue cleaning modules Removed unused arguments and variables
Property svn:executable set to ``*
File size: 93.6 KB

Line
1	#ifdef ISO
2	#ifdef ISOTRAC
3	! $Id: $
4
5	MODULE isotrac_routines_mod
6	! on créé ce module pour éviter dépendances circulaires.
7	! isotopes_routines_mod a besoin de isotrac et isotopes_verif
8	! isotopes_verif a besoin de isotopes et isotrac
9	! isotrac n'a besoin que de isotopes
10	USE infotrac_phy, ONLY: ntraciso=>ntiso, niso, index_trac=>itZonIso, ntraceurs_zone=>nzone
11	IMPLICIT NONE
12
13	CONTAINS
14
15
16	SUBROUTINE uncompress_commun_zone(ncas, cas, &
17	xtp_cas,xtp,xtwater_cas,xtwater,xtevap_cas,xtevap, &
18	ncum,izone)
19
20	USE isotopes_mod, ONLY: ridicule,iso_eau
21
22	IMPLICIT NONE
23
24	! decompression des outputs communs à tous les cas dans
25	! appel_stewart
26	! cas des traceurs. Ici, aucun risque de revap franche.
27
28	! inputs
29	INTEGER ncas,ncum
30	INTEGER cas(ncum)
31	REAL xtevap_cas(niso,ncum)
32	REAL xtp_cas(niso,ncum)
33	REAL xtwater_cas(niso,ncum)
34	INTEGER izone
35	! outputs
36	REAL xtwater(ntraciso,ncum)
37	REAL xtp(ntraciso,ncum)
38	REAL xtevap(ntraciso,ncum)
39
40	! locals
41	INTEGER il,ixt,iiso,ixt_revap
42
43
44	DO il=1,ncas
45	DO iiso=1,niso
46	ixt=index_trac(izone,iiso)
47	xtevap(ixt,cas(il))=xtevap_cas(iiso,il)
48	xtp(ixt,cas(il))=xtp_cas(iiso,il)
49	xtwater(ixt,cas(il))=xtwater_cas(iiso,il)
50	enddo !do iiso=1,niso
51	enddo !do il=1,ncas
52
53	END SUBROUTINE uncompress_commun_zone
54
55
56	SUBROUTINE uncompress_commun_zone_revap(ncas, cas, &
57	xtp_cas,xtp,xtwater_cas,xtwater, &
58	xtevap_cas,xtevap, &
59	ncum,izone,Eqi_stewart,fac_ftmr_cas, &
60	#ifdef ISOVERIF
61	Exi_cas,Exi, &
62	#endif
63	xtp_avantevap_cas,liq,hdiag)
64
65	USE isotopes_mod, ONLY: ridicule,iso_eau,iso_HDO,ridicule_evap
66	USE isotrac_mod, ONLY: option_revap,evap_franche,izone_revap, &
67	& ridicule_trac
68	#ifdef ISOVERIF
69	USE isotopes_verif_mod
70	#endif
71
72	IMPLICIT NONE
73
74	! decompression des outputs communs à tous les cas dans
75	! appel_stewart
76	! cas des traceurs: mais ici, risque de révap franche -> on fat
77	! plus attention
78
79	! inputs
80	INTEGER ncas,ncum
81	INTEGER cas(ncum)
82	REAL xtevap_cas(niso,ncum)
83	REAL xtp_cas(niso,ncum)
84	REAL xtwater_cas(niso,ncum)
85	INTEGER izone
86	REAL Eqi_stewart(ncum)
87	REAL xtp_avantevap_cas(niso,ncum)
88	REAL fac_ftmr_cas(ncum)
89	INTEGER liq
90	REAL hdiag(ncas)
91
92	! outputs
93	REAL xtwater(ntraciso,ncum)
94	REAL xtp(ntraciso,ncum)
95	REAL xtevap(ntraciso,ncum)
96
97	! locals
98	INTEGER il,ixt,iiso,ixt_revap
99	REAL xtaddp_tag(niso,ncum)
100	#ifdef ISOVERIF
101	INTEGER ieau,iHDO
102	REAL Exi_cas(niso,ncum)
103	REAL Exi(ntraciso,ncum)
104	!USE isotopes_verif, ONLY:
105	#endif
106
107	! WRITE(,) 'compress_stewart 315 tmp: ',
108	! : 'entrée dans uncompress_commun_zone_revap'
109
110	DO il=1,ncas
111	IF ((option_revap.EQ.1).AND. &
112	((((liq.EQ.1).AND. &
113	(Eqi_stewart(il)*fac_ftmr_cas(il).gt.evap_franche) &
114	.AND.(hdiag(il).lt.0.99)).OR. &
115	((liq.EQ.0).AND. &
116	(Eqi_stewart(il)*fac_ftmr_cas(il).ge.0.0))).OR. &
117	(izone.EQ.izone_revap))) THEN
118	! if ((option_revap.EQ.1).AND.
119	! : (((Eqi_stewart(il)*fac_ftmr_cas(il).gt.evap_franche)
120	! : .OR.((liq.EQ.0).AND.
121	! : (Eqi_stewart(il)*fac_ftmr_cas(il).ge.0.0))).OR.
122	! : (izone.EQ.izone_revap))) THEN
123	! on met la revap dans izone_revap si option_revap=1 et si:
124	! * evap glace (non fractionnante)
125	! * ou evap liq suffisemment forte pour que pas de flux
126	! d'isotopes négatifs.
127	! si option_revap=1 et izone=izone_revap, on met aussi dans izone_revap
128	!#ifdef ISOVERIF
129	! if (il.EQ.1) THEN
130	! WRITE(,) 'compress tmp 341: revap dans izone_revap'
131	! WRITE(,) 'Eqi_stewart(il),fac_ftmr_cas(il),evap_franche=',
132	! : Eqi_stewart(il),fac_ftmr_cas(il),evap_franche
133	! WRITE(,) 'il,xtp_cas(iso_eau,il)=',il,xtp_cas(iso_eau,il)
134	! WRITE(,) 'il,xtp_avantevap_cas(iso_eau,il)=',il,
135	! : xtp_avantevap_cas(iso_eau,il)
136	! WRITE(,) 'xtp(ixt_revap,cas(il))=',
137	! : xtp(index_trac(izone_revap,iso_eau),cas(il))
138	! endif
139	!#endif
140
141	! toute la révap franche va dans izone_revap
142	DO iiso=1,niso
143	ixt=index_trac(izone,iiso)
144	ixt_revap=index_trac(izone_revap,iiso)
145	! le terme d'évap qui était pour la zone izone devient
146	! nul, et à la place on le met dans izone_revap&
147	xtevap(ixt_revap,cas(il))=xtevap(ixt_revap,cas(il)) &
148	+xtevap_cas(iiso,il)
149
150	! ce qui a été ajouté à xtp par rapport à xtp_avantevap
151	! est ajouté à izone_revap, au lieu de izone
152	xtaddp_tag(iiso,il)=xtp_cas(iiso,il) &
153	-xtp_avantevap_cas(iiso,il)
154	xtp(ixt_revap,cas(il))= &
155	xtp(ixt_revap,cas(il)) &
156	+xtaddp_tag(iiso,il)
157	enddo !do iiso=1,niso
158	#ifdef ISOVERIF
159	IF (iso_HDO.gt.0) THEN
160	IF (xtevap_cas(iso_eau,il).gt.ridicule_evap) THEN
161	CALL iso_verif_aberrant_choix( &
162	xtevap_cas(iso_HDO,il),xtevap_cas(iso_eau,il), &
163	ridicule_trac,deltalimtrac*2,'compress 344a')
164	endif
165	ieau=index_trac(izone_revap,iso_eau)
166	iHDO=index_trac(izone_revap,iso_HDO)
167	CALL iso_verif_aberrant_choix(xtevap(iHDO,cas(il)), &
168	xtevap(ieau,cas(il)),ridicule_trac,deltalimtrac*2, &
169	'compress 344b')
170	endif !if (iso_HDO.gt.0) THEN
171	#endif
172
173	! l'évap des autres zones devient nulle
174	IF (izone.NE.izone_revap) THEN
175	DO iiso=1,niso
176	ixt=index_trac(izone,iiso)
177	xtevap(ixt,cas(il))=0.0
178	xtp(ixt,cas(il))=xtp_avantevap_cas(iiso,il)
179	enddo
180	endif
181
182	!#ifdef ISOVERIF
183	! if (il.EQ.1) THEN
184	! WRITE(,) 'compress tmp 341: revap dans izone_revap'
185	! WRITE(,) 'xtp(ixt_revap,cas(il))=',
186	! : xtp(index_trac(izone_revap,iso_eau),cas(il))
187	! WRITE(,) 'xtp(ixt,cas(il))=',
188	! : xtp(index_trac(izone,iso_eau),cas(il))
189	! endif
190	!#endif
191
192	else !if ((Eqi_stewart(il).gt.ridicule_evap*100)
193
194	DO iiso=1,niso
195	ixt=index_trac(izone,iiso)
196	xtevap(ixt,cas(il))=xtevap_cas(iiso,il)
197	xtp(ixt,cas(il))=xtp_cas(iiso,il)
198	enddo !do iiso=1,niso
199	endif !if ((Eqi_stewart(il).gt.ridicule_evap*100)
200	enddo !do il=1,ncas
201
202	DO il=1,ncas
203	DO iiso=1,niso
204	ixt=index_trac(izone,iiso)
205	xtwater(ixt,cas(il))=xtwater_cas(iiso,il)
206	#ifdef ISOVERIF
207	Exi(ixt,cas(il))=Exi_cas(iiso,il)
208	#endif
209	enddo !do iiso=1,niso
210	enddo !do il=1,ncas
211	! WRITE(,) 'compress_stewart 361 tmp: ',
212	! : 'sortie de uncompress_commun_zone_revap'
213
214	END SUBROUTINE uncompress_commun_zone_revap
215
216
217
218	SUBROUTINE compress_cond_facftmr_zone( &
219	ncas, cas, &
220	Eqi_prime_cas,Eqi_prime, &
221	Pqisup_cas,Pqisup, &
222	Pxtisup_cas,Pxtisup, &
223	qp_avantevap_cas,qp_avantevap, &
224	xtp_avantevap_cas,xtp_avantevap, &
225	xtevapsup_cas,xtevap,&
226	water_cas,water,&
227	#ifdef ISOVERIF
228	evap_cas,evap, &
229	#endif
230	nloc,ncum,nd,i,izone)
231
232	USE isotopes_mod, ONLY: iso_eau
233	#ifdef ISOVERIF
234	USE isotopes_verif_mod
235	#endif
236
237	IMPLICIT NONE
238
239	! compression dans le cas condensation_facftmr
240	! inputs
241	INTEGER nd,ncum,nloc
242	INTEGER ncas
243	INTEGER cas(ncum)
244	INTEGER i
245	INTEGER izone
246	real &
247	xtevapsup_cas(niso,ncum),xtevap(ntraciso,ncum), &
248	water_cas(ncum),water(ncum)
249	REAL Eqi_prime_cas(ncum),Eqi_prime(ncum), &
250	Pqisup_cas(ncum),Pqisup(ncum), &
251	Pxtisup_cas(niso,ncum),Pxtisup(ntraciso,ncum), &
252	qp_avantevap_cas(ncum),qp_avantevap(ncum), &
253	xtp_avantevap_cas(niso,ncum),xtp_avantevap(ntraciso,ncum)
254	#ifdef ISOVERIF
255	REAL evap_cas(ncum),evap(ncum)
256	#endif
257	! locals
258	INTEGER il,ixt,ieau,iiso
259
260	ieau=index_trac(izone,iso_eau)
261	DO il=1,ncas
262	IF (qp_avantevap(cas(il)).gt.0.0) THEN
263	Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
264	*(xtp_avantevap(ieau,cas(il))/qp_avantevap(cas(il)))
265	else !if (qp_avantevap_cas(cas(il)).gt.0.0) THEN
266	#ifdef ISOVERIF
267	CALL iso_verif_egalite_choix( &
268	(Eqi_prime(cas(il))),0.0, &
269	'compress_stewart 495',errmax,errmaxrel)
270	#endif
271	Eqi_prime_cas(il)=0.0
272	endif !if (qp_avantevap_cas(cas(il)).gt.0.0) thens
273
274	IF (Pqisup(cas(il))-Eqi_prime(cas(il)).gt.0.0) THEN
275	water_cas(il)=water(cas(il)) &
276	*((Pxtisup(ieau,cas(il))-Eqi_prime_cas(il)) &
277	/(Pqisup(cas(il))-Eqi_prime(cas(il))))
278	else !if (Pqisup(cas(il)).gt.0.0) THEN
279	#ifdef ISOVERIF
280	CALL iso_verif_egalite_choix(water(cas(il)),0.0, &
281	'compress_stewart 507',errmax,errmaxrel)
282	#endif
283	water_cas(il)=0.0
284	endif !if (Pqisup(cas(il)).gt.0.0) THEN
285	Pqisup_cas(il)=Pxtisup(ieau,cas(il))
286	qp_avantevap_cas(il)=xtp_avantevap(ieau,cas(il))
287
288	!#ifdef ISOVERIF
289	! CALL iso_verif_noNaN(water_cas(il),'compress_stewart 518')
290	! evap_cas(il)=evap(cas(il))
291	! : *(xtp_avantevap(ieau,cas(il))/qp_avantevap(cas(il)))
292	!! qp_cas(il)=xtp(ieau,cas(il))
293	!#endif
294	DO iiso=1,niso
295	ixt=index_trac(izone,iiso)
296	Pxtisup_cas(iiso,il)=Pxtisup(ixt,cas(il))
297	xtp_avantevap_cas(iiso,il)=xtp_avantevap(ixt,cas(il))
298	xtevapsup_cas(iiso,il)=xtevap(ixt,cas(il))
299	enddo
300	enddo
301
302	END SUBROUTINE compress_cond_facftmr_zone
303
304	SUBROUTINE compress_cond_nofftmr_zone( &
305	ncas, cas, &
306	Eqi_prime_cas,Eqi_prime, &
307	Pqisup_cas,Pqisup, &
308	Pxtisup_cas,Pxtisup, &
309	water_cas,water, &
310	qp_avantevap_cas,qp_avantevap, &
311	xtp_avantevap_cas,xtp_avantevap, &
312	xt_cas,xt,q_cas,q, &
313	xtevapsup_cas,xtevap, &
314	#ifdef ISOVERIF
315	evap_cas,evap, &
316	#endif
317	nloc,ncum,nd,i,izone)
318
319	USE isotopes_mod, ONLY: iso_eau
320	#ifdef ISOVERIF
321	USE isotopes_verif_mod
322	#endif
323
324	IMPLICIT NONE
325
326	! compression dans le cas condensation_facftmr
327	! inputs
328	INTEGER nloc,nd,ncum
329	INTEGER ncas
330	INTEGER cas(ncum)
331	INTEGER i
332	INTEGER izone
333	real &
334	xt_cas(niso,ncum),q_cas(ncum),xt(ntraciso,ncum),q(ncum), &
335	xtevapsup_cas(niso,ncum),xtevap(ntraciso,ncum), &
336	water_cas(ncum),water(ncum)
337	REAL Eqi_prime_cas(ncum),Eqi_prime(ncum), &
338	Pqisup_cas(ncum),Pqisup(ncum), &
339	Pxtisup_cas(niso,ncum),Pxtisup(ntraciso,ncum), &
340	qp_avantevap_cas(ncum),qp_avantevap(ncum), &
341	xtp_avantevap_cas(niso,ncum),xtp_avantevap(ntraciso,ncum)
342	#ifdef ISOVERIF
343	REAL evap_cas(ncum),evap(ncum)
344	#endif
345	! locals
346	INTEGER il,ixt,ieau,iiso
347
348	ieau=index_trac(izone,iso_eau)
349	DO il=1,ncas
350	IF (qp_avantevap(cas(il)).gt.0) THEN
351	Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
352	*(xtp_avantevap(ieau,cas(il))/qp_avantevap(cas(il)))
353	else
354	Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
355	*(xt(ieau,cas(il))/q(cas(il)))
356	endif
357
358	IF (Pqisup(cas(il))-Eqi_prime(cas(il)).gt.0.0) THEN
359	water_cas(il)=water(cas(il)) &
360	*((Pxtisup(ieau,cas(il))-Eqi_prime_cas(il)) &
361	/(Pqisup(cas(il))-Eqi_prime(cas(il))))
362	else !if (Pqisup(cas(il)).gt.0.0) THEN
363	#ifdef ISOVERIF
364	CALL iso_verif_egalite_choix(water(cas(il)),0.0, &
365	'compress_stewart 654',errmax,errmaxrel)
366	#endif
367	water_cas(il)=0.0
368	endif !if (Pqisup(cas(il)).gt.0.0) THEN
369	Pqisup_cas(il)=Pxtisup(ieau,cas(il))
370	qp_avantevap_cas(il)=xtp_avantevap(ieau,cas(il))
371	q_cas(il)=xt(ieau,cas(il))
372	!#ifdef ISOVERIF
373	! if (qp_avantevap(cas(il)).gt.0.0) THEN
374	! evap_cas(il)=evap(cas(il))
375	! : *(xtp_avantevap(ieau,cas(il))/qp_avantevap(cas(il)))
376	! else
377	! evap_cas(il)=evap(cas(il))
378	! : *(xt(ieau,cas(il))/q(cas(il)))
379	! endif
380	! qp_cas(il)=xtp(ieau,cas(il))
381	!#endif
382	DO iiso=1,niso
383	ixt=index_trac(izone,iiso)
384	Pxtisup_cas(iiso,il)=Pxtisup(ixt,cas(il))
385	xtp_avantevap_cas(iiso,il)=xtp_avantevap(ixt,cas(il))
386	xt_cas(iiso,il)=xt(ixt,cas(il))
387	xtevapsup_cas(iiso,il)=xtevap(ixt,cas(il))
388	enddo
389	enddo
390
391	END SUBROUTINE compress_cond_nofftmr_zone
392
393	SUBROUTINE compress_noevap_zone( &
394	ncas, cas, &
395	Pqisup_cas,Pqisup, &
396	Pxtisup_cas,Pxtisup, &
397	xtp_avantevap_cas,xtp_avantevap, &
398	xtevapsup_cas,xtevap, &
399	water_cas,water, &
400	#ifdef ISOVERIF
401	evap_cas,evap, &
402	#endif
403	nloc,ncum,nd,i,izone)
404
405	USE isotopes_mod, ONLY: ridicule,iso_eau
406	#ifdef ISOVERIF
407	USE isotopes_verif_mod
408	#endif
409	IMPLICIT NONE
410
411	! compression dans le cas condensation_facftmr
412	INTEGER nloc,nd,ncum
413	INTEGER ncas
414	INTEGER cas(ncum)
415	INTEGER i
416	INTEGER izone
417	REAL xtevapsup_cas(niso,ncum),xtevap(ntraciso,ncum), &
418	water_cas(ncum),water(ncum)
419	REAL Pqisup_cas(ncum),Pqisup(ncum), &
420	Pxtisup_cas(niso,ncum),Pxtisup(ntraciso,ncum), &
421	xtp_avantevap_cas(niso,ncum),xtp_avantevap(ntraciso,ncum)
422	#ifdef ISOVERIF
423	REAL evap_cas(ncum),evap(ncum)
424	#endif
425	INTEGER il,ixt,ieau,iiso
426
427	ieau=index_trac(izone,iso_eau)
428	DO il=1,ncas
429	Pqisup_cas(il)=Pxtisup(ieau,cas(il))
430	IF (Pqisup(cas(il)).gt.0.0) THEN
431	water_cas(il)=water(cas(il)) &
432	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
433	else
434	water_cas(il)=0.0
435	#ifdef ISOVERIF
436	CALL iso_verif_egalite_choix(water(cas(il)), &
437	0.0,'compress_stewart 709',errmax,errmaxrel)
438	#endif
439	endif
440	#ifdef ISOVERIF
441	! evap_cas(il)=evap(cas(il))
442	! : *(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
443	! qp_cas(il)=xtp(ieau,cas(il))
444	#endif
445	DO iiso=1,niso
446	ixt=index_trac(izone,iiso)
447	Pxtisup_cas(iiso,il)=Pxtisup(ixt,cas(il))
448	xtp_avantevap_cas(iiso,il)=xtp_avantevap(ixt,cas(il))
449	xtevapsup_cas(iiso,il)=xtevap(ixt,cas(il))
450	enddo
451	enddo
452
453	END SUBROUTINE compress_noevap_zone
454
455	SUBROUTINE compress_evap_liq_zone(iflag_con,ncas, &
456	cas, &
457	Pqisup_cas,Pqisup, &
458	Pxtisup_cas,Pxtisup, &
459	xtp_avantevap_cas,xtp_avantevap, &
460	xtp_avantevaptrac_cas,qp_avantevaptrac_cas, &
461	xtevapsup_cas,xtevap, &
462	water_cas,water, &
463	Eqi_stewart,Pqiinf_stewart,Eqi_prime_cas,&
464	Pqiinf,Eqi_par,Pqiinf_par,Eqi_prime,ptrac, &
465	Eqi,Eqi_cas, &
466	! & qp_cas,
467	#ifdef ISOVERIF
468	evap_cas,evap, &
469	#endif
470	nloc,ncum,nd,izone)
471
472	USE isotopes_mod, ONLY: ridicule,iso_eau
473	#ifdef ISOVERIF
474	USE isotopes_verif_mod
475	#endif
476	IMPLICIT NONE
477
478	! compression dans le cas condensation_facftmr
479	! inputs et outputs
480	INTEGER iflag_con
481	INTEGER izone
482	INTEGER nloc,nd,ncum
483	INTEGER ncas
484	INTEGER cas(ncum)
485	! integer i
486	REAL xtevapsup_cas(niso,ncum),xtevap(ntraciso,ncum)
487	REAL xtp_avantevap_cas(niso,ncum), &
488	xtp_avantevap(ntraciso,ncum)
489	REAL water_cas(ncum),water(ncum)
490	REAL xtp_avantevaptrac_cas(niso,ncum), &
491	qp_avantevaptrac_cas(ncum)
492	REAL ptrac(ncum)
493	! real qp_cas(ncum)
494	#ifdef ISOVERIF
495	REAL evap_cas(ncum),evap(ncum)
496	#endif
497	REAL &
498	Eqi_stewart(ncum),Pqiinf_stewart(ncum),Eqi_prime_cas(ncum), &
499	Pqiinf(ncum),Eqi_par(ncum),Pqiinf_par(ncum), &
500	Eqi_prime(ncum),Pqisup(ncum),Pqisup_cas(ncum), &
501	Pxtisup(ntraciso,ncum),Pxtisup_cas(niso,ncum), &
502	Eqi(ncum),Eqi_cas(ncum)
503	! locals
504	INTEGER il,ixt,iiso,ieau
505
506	! WRITE(,) 'compress 910: xtp_avantevap(iso_eau,cas(1))=',
507	! : xtp_avantevap(iso_eau,cas(1))
508	! WRITE(,) 'compress_evap_liq_zone 510: ncas,ncum=',ncas,ncum
509	ptrac(:)=0. ! CR 31 mars 2023: initialisation de ptrac
510
511	ieau=index_trac(izone,iso_eau)
512	DO il=1,ncas
513	Pqisup_cas(il)=Pxtisup(ieau,cas(il))
514
515	IF (Pqisup(cas(il)).gt.0.0) THEN
516	ptrac(il)=Pxtisup(ieau,cas(il))/Pqisup(cas(il))
517	Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
518	*ptrac(il)
519	else
520	#ifdef ISOVERIF
521	CALL iso_verif_egalite(( &
522	Eqi_prime(cas(il))),0.0, &
523	'compress_stewart 979')
524	#endif
525	Eqi_prime_cas(il)=0.0
526	endif
527
528	IF (Pqisup(cas(il))-Eqi_prime(cas(il)).gt.0.0) THEN
529	water_cas(il)=water(cas(il)) &
530	*((Pxtisup(ieau,cas(il))-Eqi_prime_cas(il)) &
531	/(Pqisup(cas(il))-Eqi_prime(cas(il))))
532	else !if (Pqisup(cas(il)).gt.0.0) THEN
533	#ifdef ISOVERIF
534	CALL iso_verif_egalite_choix(water(cas(il)),0.0, &
535	'compress_stewart 507',errmax,errmaxrel)
536	#endif
537	water_cas(il)=0.0
538	endif !if (Pqisup(cas(il)).gt.0.0) THEN
539	! qp_cas(il)=qp(cas(il))
540	#ifdef ISOVERIF
541	! evap_cas(il)=evap(cas(il))
542	! & *(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
543	! ! ce calcul est faux& normalement,
544	! evap_cas(il)=Eqi_prime_cas(il)//100.0/delP_cas(il)/SIGDg2.0
545	! & -(ieau,cas(il))
546	#endif
547	qp_avantevaptrac_cas(il)=xtp_avantevap(ieau,cas(il))
548	DO iiso=1,niso
549	ixt=index_trac(izone,iiso)
550	Pxtisup_cas(iiso,il)=Pxtisup(ixt,cas(il))
551	xtp_avantevap_cas(iiso,il)=xtp_avantevap(iiso,cas(il))
552	xtp_avantevaptrac_cas(iiso,il)=xtp_avantevap(ixt,cas(il))
553	xtevapsup_cas(iiso,il)=xtevap(ixt,cas(il))
554	enddo
555	enddo !do il=1,ncas
556	! WRITE(,) 'compress_stewart 931: ',
557	! & 'xtp_avantevap_cas(iso_eau,1)=',
558	! & xtp_avantevap_cas(iso_eau,1)
559	! WRITE(,) 'xtp_avantevap(iso_eau,cas(1))=',
560	! & xtp_avantevap(iso_eau,cas(1))
561
562	! calculs des flux de masses à mettre en argument de stewart:
563	! comme l'eau n'est pas bien concervée dans les ddfts, on est
564	! obligé de bidouillé.
565	! 1) soit on considère Pqisup, Eqi, et Pqiinf_par=Pqisup-Eqi
566	! et on suppose que dans la réalité les compositions de
567	! Pqiinf sont les même que Pqiinf_par
568	! 2) soit on considère Pqisup, Eqi_par=Pqisup-Pqiinf, et Pqiinf,
569	! et on suppose que dans la réalité les compositions de
570	! Eqi_prime sont les même que Eqi_par
571	DO il=1,ncas
572	IF (Pqisup(cas(il)).gt.0.0) THEN
573	IF ((water(cas(il)).gt.ridicule/100).AND. &
574	(Pqiinf_par(cas(il)).le.0.0)) THEN
575	! on ne peut pas utiliser la méthode 1, car KE prédit de l'eau
576	! alors que le bilan de masse n'enprédit pas.
577	! Peut-on utiliser la méthode 2?
578	Pqiinf_stewart(il)=Pqiinf(cas(il))*ptrac(il)
579	Eqi_stewart(il)=Eqi_par(cas(il)) *ptrac(il)
580	else !if ((water(il,i).gt.ridicule/100).AND.(Pqiinf_par.le.0.0)) THEN
581	! il n'y a pas d'obstacles à l'utilisation de 1)
582	Pqiinf_stewart(il)=Pqiinf_par(cas(il))*ptrac(il)
583	IF (iflag_con.EQ.30) THEN
584	Eqi_stewart(il)=Eqi_prime(cas(il))*ptrac(il)
585	else
586	Eqi_stewart(il)=Eqi(cas(il))*ptrac(il)
587	endif
588	endif !if ((water(il,i).gt.ridicule/100).AND.(Pqiinf_par.le.0.0)) THEN
589	else ! if (Pqisup(cas(il).gt.0.0) THEN
590	#ifdef ISOVERIF
591	CALL iso_verif_egalite((Pqiinf(cas(il))), &
592	0.0,'compress_stewart 1047a')
593	CALL iso_verif_egalite(( &
594	Eqi_prime(cas(il))),0.0,'compress_stewart 1047b')
595	CALL iso_verif_egalite(( &
596	Pqiinf_par(cas(il))),0.0,'compress_stewart 1047c')
597	CALL iso_verif_egalite((Eqi_par(cas(il))), &
598	0.0,'compress_stewart 1047d')
599	#endif
600	Pqiinf_stewart(il)=0.0
601	Eqi_stewart(il)=0.0
602	endif ! if (Pqisup(cas(il).gt.0.0) THEN
603	enddo !do il=1,ncas_evap_glace
604
605	! petite vérif
606	#ifdef ISOVERIF
607	DO il=1,ncas
608	IF ((abs(water_cas(il)).ge.ridicule/10.) &
609	.AND.(Pqiinf_stewart(il).le.0.0)) THEN
610	WRITE(,) 'compress_stewart 498: evap liq:'
611	WRITE(,) 'water(il,i)=', water_cas(il)
612	WRITE(,) 'Pqiinf=',Pqiinf(cas(il))
613	WRITE(,) 'Pqiinf_par=',Pqiinf_par(cas(il))
614	WRITE(,) 'Pqiinf_stewart=',Pqiinf_stewart(il)
615	stop
616	endif
617	enddo !do il=1,ncas_evap_glace
618	#endif
619
620	END SUBROUTINE compress_evap_liq_zone
621
622	SUBROUTINE compress_evap_glace_zone(iflag_con, &
623	ncas, cas, &
624	water_cas,water, &
625	Pqisup_cas,Pqisup, &
626	Pxtisup_cas,Pxtisup, &
627	xtp_avantevap_cas,xtp_avantevap, &
628	xtp_avantevaptrac_cas,qp_avantevaptrac_cas, &
629	xtevapsup_cas,xtevap, &
630	Eqi_stewart,Pqiinf_stewart,Eqi_prime_cas,Eqi_cas, &
631	Pqiinf,Eqi_par,Pqiinf_par,Eqi_prime,Eqi, &
632	! & qp_cas,
633	#ifdef ISOVERIF
634	evap_cas,evap, &
635	#endif
636	nloc,ncum,nd,i,frac_sublim,izone)
637
638	USE isotopes_mod, ONLY: ridicule,iso_eau
639	#ifdef ISOVERIF
640	USE isotopes_verif_mod
641	#endif
642	IMPLICIT NONE
643
644	! compression dans le cas condensation_facftmr
645	INTEGER iflag_con
646	INTEGER nloc,nd,ncum
647	INTEGER ncas
648	INTEGER cas(ncum)
649	INTEGER i
650	INTEGER izone
651	real &
652	water_cas(ncum),water(ncum), &
653	xtevapsup_cas(niso,ncum),xtevap(ntraciso,ncum)
654	! real qp_cas(ncum)
655	#ifdef ISOVERIF
656	REAL evap_cas(ncum),evap(ncum)
657	#endif
658	real &
659	Pqisup_cas(ncum),Pqisup(ncum), &
660	Pxtisup_cas(niso,ncum),Pxtisup(ntraciso,ncum), &
661	xtp_avantevap_cas(niso,ncum),xtp_avantevap(ntraciso,ncum), &
662	Eqi_stewart(ncum),Pqiinf_stewart(ncum),Eqi_prime_cas(ncum), &
663	Pqiinf(ncum),Eqi_par(ncum),Pqiinf_par(ncum),Eqi_prime(ncum), &
664	Eqi(ncum),Eqi_cas(ncum)
665
666	REAL xtp_avantevaptrac_cas(niso,ncum), &
667	qp_avantevaptrac_cas(ncum)
668	INTEGER frac_sublim
669	! locals
670	INTEGER il,ixt,ieau,iiso
671
672	ieau=index_trac(izone,iso_eau)
673	DO il=1,ncas
674	Pqisup_cas(il)=Pxtisup(ieau,cas(il))
675
676	IF (Pqisup(cas(il)).gt.0.0) THEN
677	Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
678	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
679	Eqi_cas(il)=Eqi(cas(il)) & ! corr bug Camille 15 juin 2024
680	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
681	else
682	#ifdef ISOVERIF
683	CALL iso_verif_egalite(( &
684	Eqi_prime(cas(il))),0.0, &
685	'compress_stewart 979b')
686	#endif
687	Eqi_prime_cas(il)=0.0
688	Eqi_cas(il)=0.0
689	endif
690
691	IF (Pqisup(cas(il))-Eqi_prime(cas(il)).gt.0.0) THEN
692	water_cas(il)=water(cas(il)) &
693	*((Pxtisup(ieau,cas(il))-Eqi_prime_cas(il)) &
694	/(Pqisup(cas(il))-Eqi_prime(cas(il))))
695	else !if (Pqisup(cas(il)).gt.0.0) THEN
696	#ifdef ISOVERIF
697	CALL iso_verif_egalite_choix(water(cas(il)),0.0, &
698	'compress_stewart 507',errmax,errmaxrel)
699	#endif
700	water_cas(il)=0.0
701	endif !if (Pqisup(cas(il)).gt.0.0) THEN
702	qp_avantevaptrac_cas(il)=xtp_avantevap(ieau,cas(il))
703	! qp_cas(il)=xtp(ieau,cas(il))
704	#ifdef ISOVERIF
705	! evap_cas(il)=evap(cas(il))
706	! & *(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
707	! ce calcul est faux& faire plutot:
708	! evap_cas(il)=Eqi_prime_cas(il)//100.0/delP_cas(il)/SIGDg2.0
709	! & -(ieau,cas(il))
710	#endif
711	DO iiso=1,niso
712	ixt=index_trac(izone,iiso)
713	Pxtisup_cas(iiso,il)=Pxtisup(ixt,cas(il))
714	xtp_avantevap_cas(iiso,il)=xtp_avantevap(iiso,cas(il))
715	xtp_avantevaptrac_cas(iiso,il)=xtp_avantevap(ixt,cas(il))
716	xtevapsup_cas(iiso,il)=xtevap(ixt,cas(il))
717	enddo
718	enddo !do il=1,ncas
719
720	! calculs des flux de masses à mettre en argument de stewart:
721	! comme l'eau n'est pas bien concervée dans les ddfts, on est
722	! obligé de bidouillé.
723	! 1) soit on considère Pqisup, Eqi, et Pqiinf_par=Pqisup-Eqi
724	! et on suppose que dans la réalité les compositions de
725	! Pqiinf sont les même que Pqiinf_par
726	! 2) soit on considère Pqisup, Eqi_par=Pqisup-Pqiinf, et Pqiinf,
727	! et on suppose que dans la réalité les compositions de
728	! Eqi_prime sont les même que Eqi_par
729	DO il=1,ncas
730	IF (Pqisup(cas(il)).gt.0.0) THEN
731	IF ((water(cas(il)).gt.ridicule/100).AND. &
732	(Pqiinf_par(cas(il)).le.0.0)) THEN
733	! on ne peut pas utiliser la méthode 1, car KE prédit de l'eau
734	! alors que le bilan de masse n'enprédit pas.
735	! Peut-on utiliser la méthode 2?
736	Pqiinf_stewart(il)=Pqiinf(cas(il)) &
737	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
738	Eqi_stewart(il)=Eqi_par(cas(il)) &
739	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
740	else !if ((water(il,i).gt.ridicule/100).AND.(Pqiinf_par.le.0.0)) THEN
741	! il n'y a pas d'obstacles à l'utilisation de 1)
742	Pqiinf_stewart(il)=Pqiinf_par(cas(il)) &
743	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
744	IF (iflag_con.EQ.30) THEN
745	Eqi_stewart(il)=Eqi_prime(cas(il)) &
746	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
747	else
748	Eqi_stewart(il)=Eqi(cas(il)) &
749	*(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
750	endif
751	endif !if ((water(il,i).gt.ridicule/100).AND.(Pqiinf_par.le.0.0)) THEN
752	else !if (Pqisup(cas(il).gt.0.0) THEN
753	#ifdef ISOVERIF
754	CALL iso_verif_egalite((Pqiinf(cas(il))), &
755	0.0,'compress_stewart 1347a')
756	CALL iso_verif_egalite(( &
757	Eqi_prime(cas(il))),0.0,'compress_stewart 1347b')
758	CALL iso_verif_egalite(( &
759	Pqiinf_par(cas(il))),0.0,'compress_stewart 1347c')
760	CALL iso_verif_egalite((Eqi_par(cas(il))), &
761	0.0,'compress_stewart 1347d')
762	#endif
763	Pqiinf_stewart(il)=0.0
764	Eqi_stewart(il)=0.0
765	endif !if (Pqisup(cas(il).gt.0.0) THEN
766	enddo !do il=1,ncas_evap_glace
767
768	! petite vérif
769	#ifdef ISOVERIF
770	DO il=1,ncas
771	IF ((abs(water_cas(il)).ge.ridicule/10.) &
772	.AND.(Pqiinf_stewart(il).le.0.0)) THEN
773	WRITE(,) 'compress_stewart 498: evap glace:'
774	WRITE(,) 'water(il,i)=', water_cas(il)
775	WRITE(,) 'Pqiinf=',Pqiinf(cas(il))
776	WRITE(,) 'Pqiinf_par=',Pqiinf_par(cas(il))
777	WRITE(,) 'Pqiinf_stewart=',Pqiinf_stewart(il)
778	stop
779	endif
780	enddo !do il=1,ncas_evap_glace
781	#endif
782
783	END SUBROUTINE compress_evap_glace_zone
784
785	SUBROUTINE uncompress_ilp_zone( &
786	ncas,cas, &
787	zxtrfln_cas,zxt_cas,zxtrfl,zxtrfln,zxt,klon, &
788	izone,Eqi,Exi,fac_ftmr, &
789	xtrevap_tag,liq,hdiag)
790
791	USE isotopes_mod, ONLY: ridicule,iso_eau
792	USE isotrac_mod, ONLY: option_revap,evap_franche
793	#ifdef ISOVERIF
794	USE isotopes_verif_mod
795	#endif
796
797	IMPLICIT NONE
798
799	! inputs
800	INTEGER ncas
801	INTEGER cas(ncas)
802	INTEGER klon
803	INTEGER izone
804	REAL zxt_cas(niso,ncas),zxtrfln_cas(niso,ncas)
805	REAL Exi(niso,ncas)
806	REAL Eqi(ncas)
807	REAL fac_ftmr(ncas)
808	INTEGER liq
809	REAL hdiag(ncas)
810
811	! output
812	REAL zxt(ntraciso,klon)
813	REAL zxtrfl(ntraciso,klon),zxtrfln(ntraciso,klon)
814	REAL xtrevap_tag(ntraciso,ncas)
815
816	! locals
817	INTEGER il,ixt,iiso
818
819	! WRITE(,) 'uncompress_stewart 1420 tmp: zxt=',
820	! : zxt(iso_eau:ntraciso:3,cas(1))
821	! WRITE(,) 'Exi,fac_ftmr=',
822	! : Exi(iso_eau,1),fac_ftmr(1)
823	DO il=1,ncas
824	DO iiso=1,niso
825	ixt=index_trac(izone,iiso)
826	zxtrfln(ixt,cas(il))=zxtrfln_cas(iiso,il)
827	zxtrfl(ixt,cas(il))=zxtrfln_cas(iiso,il)
828	enddo !do iiso=1,niso
829	!#ifdef ISOVERIF
830	! if (il.EQ.9) THEN
831	! WRITE(,) 'uncompress 1521'
832	! WRITE(,) 'il,Eqi,fac_ftmr,evap_franche,Exi(2,il)=',
833	! : il,Eqi(il),fac_ftmr(il),evap_franche,Exi(2,il)
834	! endif
835	!#endif
836
837	IF ((option_revap.EQ.1).AND. &
838	(((liq.EQ.1).AND.(Eqi(il)*fac_ftmr(il).gt.evap_franche) &
839	.AND.(hdiag(il).lt.0.99)).OR. &
840	((liq.EQ.0).AND. &
841	(Eqi(il)*fac_ftmr(il).ge.0.0)))) THEN
842	! le flux d'évap va dans un tag particulier
843	! -> zxt est inchangé mais xtrevap_tag(ixt,il) est incrémenté
844	DO iiso=1,niso
845	ixt=index_trac(izone,iiso)
846	xtrevap_tag(ixt,il)=fac_ftmr(il)*Exi(iiso,il)
847	! zxt(ixt,cas(il))=zxt_cas(iiso,il)
848	! : -xtrevap_tag(ixt,il)
849	enddo !do iiso=1,niso
850
851	else !if ((Eqi(il)*fac_ftmr(il).gt.evap_franche).AND.
852	! reequilibration standard
853	DO iiso=1,niso
854	ixt=index_trac(izone,iiso)
855	zxt(ixt,cas(il))=zxt(ixt,cas(il)) &
856	+fac_ftmr(il)*Exi(iiso,il)
857	zxt(ixt,cas(il))=max(0.0,zxt(ixt,cas(il)))
858	xtrevap_tag(ixt,il)=0.0
859	#ifdef ISOVERIF
860	CALL iso_verif_positif_choix(zxt(ixt,cas(il)), &
861	0.0,'compress 1508')
862	#endif
863	enddo ! do iiso=1,niso
864	endif !if ((Eqi(il)*fac_ftmr(il).gt.evap_franche).AND.
865	enddo !do il=1,ncas
866	! WRITE(,) 'compress_stewart 1453 tmp: zxt=',
867	! : zxt(iso_eau:ntraciso:3,cas(1))
868
869	END SUBROUTINE uncompress_ilp_zone
870
871	SUBROUTINE compress_ilp_evap_liq_zone( &
872	ncas,cas, &
873	zxt_cas,zxt, &
874	zxtrfl_cas,zxtrfl_ancien, &
875	zrfln_cas,zrfln, &
876	zrfl_cas,zrfl_ancien, &
877	zqev_diag_cas,zqev_diag, &
878	klon,izone,ptrac)
879
880	USE isotopes_mod, ONLY: ridicule,iso_eau
881	#ifdef ISOVERIF
882	USE isotopes_verif_mod
883	#endif
884	IMPLICIT NONE
885
886	! inputs
887	INTEGER ncas
888	INTEGER cas(ncas)
889	INTEGER klon
890	REAL zxt(ntraciso,klon)
891	REAL zxtrfl(ntraciso,klon)
892	REAL zrfl_ancien(klon)
893	REAL zqev_diag(klon)
894	REAL zrfln(klon)
895	REAL zxtrfl_ancien(ntraciso,klon)
896	INTEGER izone
897
898	! outputs
899	REAL zxt_cas(niso,ncas)
900	REAL zxtrfl_cas(niso,ncas)
901	REAL zqev_diag_cas(ncas)
902	REAL zrfln_cas(ncas)
903	REAL zrfl_cas(ncas)
904	REAL ptrac(ncas)
905
906	! locals
907	INTEGER il,ixt,ieau,iiso
908
909	ieau=index_trac(izone,iso_eau)
910	DO il=1,ncas
911	DO iiso=1,niso
912	ixt=index_trac(izone,iiso)
913	! la compo de la vap à l'extérieure reste la vapeur totale
914	zxt_cas(iiso,il)=zxt(iiso,cas(il))
915	! le flux de pluie est celui le flux de pluie lié à la zone
916	zxtrfl_cas(iiso,il)=zxtrfl_ancien(ixt,cas(il))
917	enddo
918	zrfl_cas(il)=zxtrfl_ancien(ieau,cas(il))
919
920	IF (zrfl_ancien(cas(il)).gt.0.0) THEN
921	! proportion de izone dans l'évap = celle dans la goutte
922	ptrac(il)=zxtrfl_ancien(ieau,cas(il))/zrfl_ancien(cas(il))
923	zrfln_cas(il)=zrfln(cas(il))*ptrac(il)
924	zqev_diag_cas(il)=zqev_diag(cas(il))*ptrac(il)
925	!#ifdef ISOVERIF
926	! if (il.EQ.9) THEN
927	! WRITE(,) 'compress tmp: il, ptrac=',il,ptrac(il)
928	! WRITE(,) 'ieau,zxtrfl_ancien(ieau,cas(il))=',
929	! : ieau,zxtrfl_ancien(ieau,cas(il))
930	! WRITE(,) 'zrfl_ancien(cas(il))=',zrfl_ancien(cas(il))
931	! WRITE(,) 'zrfl_cas(il)=',zrfl_cas(il)
932	! endif
933	!#endif
934
935	else !if (zrfl_ancien(cas(il)).gt.0.0) THEN
936	#ifdef ISOVERIF
937	! WRITE(,) 'il,cas(il),zrfln,zrfl_ancien,zqev_diag=',
938	! : il,cas(il),zrfln(cas(il)),zrfl_ancien(cas(il)),
939	! : zqev_diag(cas(il))
940	CALL iso_verif_egalite(zqev_diag(cas(il)),0.0, &
941	'compress_stewart 1591a')
942	CALL iso_verif_egalite(zrfln(cas(il)),0.0, &
943	'compress_stewart 1591b')
944	#endif
945	zrfln_cas(il)=0.0
946	zqev_diag_cas(il)=0.0
947
948	endif !if (zrfl_ancien(cas(il)).gt.0.0) THEN
949	! les lignes suvantes ne sont pas à recalculer
950	! zt_cas(il)=zt(cas(il))
951	! delP(il)=paprs(cas(il),k)-paprs(cas(il),k+1)
952	enddo !do il=1,ncas
953
954	END SUBROUTINE compress_ilp_evap_liq_zone
955
956
957	SUBROUTINE compress_ilp_evap_glace_zone( &
958	ncas,cas, &
959	zxt_cas,zxt, &
960	zxtrfl_cas,zxtrfl_ancien, &
961	zrfln_cas,zrfln, &
962	zrfl_cas, zrfl_ancien, &
963	zqev_diag_cas,zqev_diag, &
964	klon,izone)
965
966	USE isotopes_mod, ONLY: ridicule,iso_eau
967	#ifdef ISOVERIF
968	USE isotopes_verif_mod
969	#endif
970
971	IMPLICIT NONE
972
973	! inputs
974	INTEGER ncas
975	INTEGER cas(ncas)
976	INTEGER klon
977	REAL zxt(ntraciso,klon)
978	REAL zxtrfl_ancien(ntraciso,klon)
979	REAL zqev_diag(klon)
980	REAL zrfln(klon)
981	INTEGER izone
982	REAL zrfl_ancien(klon)
983
984	! outputs
985	REAL zxt_cas(niso,ncas)
986	REAL zxtrfl_cas(niso,ncas)
987	REAL zqev_diag_cas(ncas)
988	REAL zrfln_cas(ncas)
989	REAL zrfl_cas(ncas)
990
991
992	! locals
993	INTEGER il,ixt,ieau,iiso
994
995	ieau=index_trac(izone,iso_eau)
996	DO il=1,ncas
997	DO iiso=1,niso
998	ixt=index_trac(izone,iiso)
999	zxt_cas(iiso,il)=zxt(iiso,cas(il))
1000	zxtrfl_cas(iiso,il)=zxtrfl_ancien(ixt,cas(il))
1001	enddo
1002	zrfl_cas(il)=zxtrfl_ancien(ieau,cas(il))
1003
1004	IF (zrfl_ancien(cas(il)).gt.0.0) THEN
1005	zrfln_cas(il)=zrfln(cas(il)) &
1006	*(zxtrfl_ancien(ieau,cas(il))/zrfl_ancien(cas(il)))
1007	! car la proportion de traceurs dans zqev_diag et la même
1008	! que dans zrfl_ancien. Comme zrfln=zrfl-zqev_diag*fac_ftmr
1009	! alors cette proportion de traceurs est aussi la même dans
1010	! zrfln
1011	zqev_diag_cas(il)=zqev_diag(cas(il)) &
1012	*zxtrfl_ancien(ieau,cas(il))/zrfl_ancien(cas(il))
1013	else !if (zrfl_ancien(cas(il)).gt.0.0) THEN
1014	#ifdef ISOVERIF
1015	CALL iso_verif_egalite(zqev_diag(cas(il)),0.0, &
1016	'compress_stewart 1791a')
1017	CALL iso_verif_egalite(zrfln(cas(il)),0.0, &
1018	'compress_stewart 1791b')
1019	#endif
1020	zrfln_cas(il)=0.0
1021	zqev_diag_cas(il)=0.0
1022	endif !if (zrfl_ancien(cas(il)).gt.0.0) THEN
1023	enddo
1024
1025	END SUBROUTINE compress_ilp_evap_glace_zone
1026
1027
1028	SUBROUTINE ajoute_revap(ncas,cas, &
1029	klon,izone,zxt,xtrevap_tag)
1030
1031	#ifdef ISOVERIF
1032	USE isotopes_verif_mod
1033	#endif
1034	USE isotrac_mod, ONLY: izone_revap
1035	IMPLICIT NONE
1036
1037	! ajoute xtrevap_tag (evaps des différents traceurs d'isotopes)
1038	! dans la vapeur qui est taggée par la révap des gouttes.
1039
1040	! input
1041	INTEGER ncas
1042	INTEGER cas(ncas)
1043	INTEGER klon
1044	INTEGER izone
1045	REAL xtrevap_tag(ntraciso,ncas)
1046
1047	! inout
1048	REAL zxt(ntraciso,klon)
1049
1050	! local
1051	INTEGER i,ixt,iiso,ixt_revap
1052	!#ifdef ISOVERIF
1053	! integer iso_verif_positif_choix_nostop
1054	!#endif
1055
1056	#ifdef ISOVERIF
1057	DO i=1,ncas
1058	DO ixt=1,ntraciso
1059	CALL iso_verif_positif_choix(zxt(ixt,cas(i)),0.0, &
1060	'ajoute_revap 29')
1061	enddo
1062	enddo !do i=1,ncas
1063	#endif
1064
1065	DO i=1,ncas
1066	DO iiso=1,niso
1067	ixt_revap=index_trac(izone_revap,iiso)
1068	DO izone=1,ntraceurs_zone
1069	ixt=index_trac(izone,iiso)
1070	zxt(ixt_revap,cas(i))=zxt(ixt_revap,cas(i)) &
1071	+xtrevap_tag(ixt,i)
1072	#ifdef ISOVERIF
1073	IF (iso_verif_positif_choix_nostop(zxt(ixt_revap,cas(i)), &
1074	0.0,'ajoute_revap 46').EQ.1) THEN
1075	WRITE(,) 'i,iiso,izone,ixt=',i,iiso,izone,ixt
1076	WRITE(,) 'xtrevap_tag(ixt,i)=',xtrevap_tag(ixt,i)
1077	! stop
1078	endif
1079	#endif
1080	zxt(ixt_revap,cas(i))=max(0.0,zxt(ixt_revap,cas(i)))
1081	enddo !do izone=1,ntraceurs_zone
1082	enddo !do iiso=1,niso
1083	enddo !do i=1,ncas_evap_liq
1084
1085	#ifdef ISOVERIF
1086	DO i=1,ncas
1087	DO ixt=1,ntraciso
1088	CALL iso_verif_positif_choix(zxt(ixt,cas(i)),0.0, &
1089	'ajoute_revap 40')
1090	enddo
1091	enddo !do i=1,ncas
1092	#endif
1093
1094
1095	END SUBROUTINE ajoute_revap
1096
1097
1098	function is_in_bassin(lat,lon,bassin)
1099	USE isotrac_mod, ONLY: use_bassin_atlantic,use_bassin_medit, &
1100	& use_bassin_indian,use_bassin_austral,use_bassin_pacific, &
1101	& use_bassin_MerArabie,use_bassin_BengalGolf,use_bassin_SouthIndian, &
1102	& use_bassin_tropics,use_bassin_midlats,use_bassin_HighLats, &
1103	& bassin_atlantic,bassin_medit, &
1104	& bassin_indian,bassin_austral,bassin_pacific, &
1105	& bassin_MerArabie,bassin_BengalGolf,bassin_SouthIndian, &
1106	& bassin_tropics,bassin_midlats,bassin_HighLats
1107	IMPLICIT NONE
1108	! répond true si lat,lon se trouve dans le bassin numéroté bassin
1109
1110	! input
1111	INTEGER bassin
1112	REAL lat,lon
1113
1114	! outut
1115	LOGICAL is_in_bassin
1116
1117	! locals
1118	!LOGICAL is_in_rectangle
1119	!LOGICAL is_in_triangle
1120
1121	is_in_bassin=.FALSE.
1122	#ifdef ISOVERIF
1123	WRITE(,) 'is_in_basin 84: entree,bassin=',bassin
1124	#endif
1125	IF (use_bassin_atlantic .AND. bassin==bassin_atlantic) THEN
1126	#ifdef ISOVERIF
1127	WRITE(,) 'bassin Atlantique?'
1128	#endif
1129	IF (is_in_rectangle(lon,lat,-67.0,28.0,20.0,-45.0)) THEN
1130	! boite sud
1131	is_in_bassin=.TRUE.
1132	RETURN
1133	endif
1134	IF (is_in_rectangle(lon,lat,-100.0,40.0,-5.3,28.0)) THEN
1135	! ouest gibraltar
1136	is_in_bassin=.TRUE.
1137	RETURN
1138	endif
1139	IF (is_in_rectangle(lon,lat,-100.0,48.0,0.0,40.0)) THEN
1140	! Ouest France
1141	is_in_bassin=.TRUE.
1142	RETURN
1143	endif
1144	IF (is_in_rectangle(lon,lat,-90.0,80.0,10.0,46.0)) THEN
1145	! Atlantic Nord
1146	is_in_bassin=.TRUE.
1147	RETURN
1148	endif
1149	IF (is_in_triangle(lon,lat, &
1150	-62.0,0.0,-62.0,30.0,-112.0,30.0)) THEN
1151	! golfe du Mexique
1152	is_in_bassin=.TRUE.
1153	RETURN
1154	endif
1155
1156	ELSE IF (use_bassin_medit .AND. bassin==bassin_medit) THEN
1157	#ifdef ISOVERIF
1158	WRITE(,) 'bassin Medit?'
1159	#endif
1160	IF (is_in_rectangle(lon,lat,0.0,48.0,45.0,29.0)) THEN
1161	is_in_bassin=.TRUE.
1162	RETURN
1163	endif
1164	IF (is_in_rectangle(lon,lat,-5.3,42.0,45.0,29.0)) THEN
1165	is_in_bassin=.TRUE.
1166	RETURN
1167	endif
1168
1169	ELSE IF (use_bassin_indian .AND. bassin==bassin_indian) THEN
1170	#ifdef ISOVERIF
1171	WRITE(,) 'bassin indian?'
1172	#endif
1173	IF (is_in_rectangle(lon,lat,20.0,30.0,110.0,-45.0)) THEN
1174	is_in_bassin=.TRUE.
1175	RETURN
1176	endif
1177	IF (is_in_triangle(lon,lat, &
1178	90.0,30.0,90.0,-45.0,150.0,-45.0)) THEN
1179	! Ouest Australie
1180	is_in_bassin=.TRUE.
1181	RETURN
1182	endif
1183
1184	ELSE IF (use_bassin_SouthIndian .AND. bassin==bassin_SouthIndian) THEN
1185	#ifdef ISOVERIF
1186	WRITE(,) 'bassin indian hemisphere Sud?'
1187	#endif
1188	IF (is_in_rectangle(lon,lat,20.0,0.0,120.0,-45.0)) THEN
1189	is_in_bassin=.TRUE.
1190	RETURN
1191	endif
1192
1193	ELSE IF (use_bassin_MerArabie .AND. bassin==bassin_MerArabie) THEN
1194	#ifdef ISOVERIF
1195	WRITE(,) 'bassin Mer d''Arabie?'
1196	#endif
1197	IF (is_in_rectangle(lon,lat,20.0,30.0,76.0,0.0)) THEN
1198	is_in_bassin=.TRUE.
1199	RETURN
1200	endif
1201
1202	ELSE IF (use_bassin_BengalGolf .AND. bassin==bassin_BengalGolf) THEN
1203	#ifdef ISOVERIF
1204	WRITE(,) 'bassin Golfe du Bengale?'
1205	#endif
1206	IF (is_in_rectangle(lon,lat,76.0,30.0,110.0,0.0)) THEN
1207	is_in_bassin=.TRUE.
1208	RETURN
1209	endif
1210
1211	ELSE IF (use_bassin_pacific .AND. bassin==bassin_pacific) THEN
1212	#ifdef ISOVERIF
1213	WRITE(,) 'bassin Pacific?'
1214	#endif
1215	IF (is_in_rectangle(lon,lat,-180.0,80.0,-100.0,-45.0)) THEN
1216	! pacifique Est
1217	is_in_bassin=.TRUE.
1218	RETURN
1219	endif
1220	IF (is_in_rectangle(lon,lat,110.0,80.0,180.0,28.0)) THEN
1221	! Pacifique Nord Ouest
1222	is_in_bassin=.TRUE.
1223	RETURN
1224	endif
1225	IF (is_in_rectangle(lon,lat,120.0,80.0,180.0,-45.0)) THEN
1226	! Pacifique central Sud
1227	is_in_bassin=.TRUE.
1228	RETURN
1229	endif
1230	IF (is_in_triangle(lon,lat, &
1231	90.0,28.0,150.0,-45.0,150.0,28.0)) THEN
1232	! Pacifque Sud Ouest
1233	is_in_bassin=.TRUE.
1234	RETURN
1235	endif
1236	IF (is_in_triangle(lon,lat, &
1237	-62.0,0.0,-112.0,30.0,-112.0,0.0)) THEN
1238	! Ouest Amérique centrale
1239	is_in_bassin=.TRUE.
1240	RETURN
1241	endif
1242	IF (is_in_rectangle(lon,lat,-180.0,0.0,-67.0,-45.0)) THEN
1243	! Ouest Chili
1244	is_in_bassin=.TRUE.
1245	RETURN
1246	endif
1247
1248	ELSE IF (use_bassin_austral .AND. bassin==bassin_austral) THEN
1249	#ifdef ISOVERIF
1250	WRITE(,) 'bassin austral?'
1251	#endif
1252	IF (lat.lt.-45.0+0.2) THEN
1253	is_in_bassin=.TRUE.
1254	return
1255	endif
1256
1257	ELSE IF (use_bassin_HighLats .AND. bassin==bassin_HighLats) THEN
1258	#ifdef ISOVERIF
1259	WRITE(,) 'bassin hautes lats?'
1260	#endif
1261	IF (abs(lat).gt.35.0) THEN
1262	is_in_bassin=.TRUE.
1263	return
1264	endif
1265
1266	ELSE IF (use_bassin_tropics .AND. bassin==bassin_tropics) THEN
1267	#ifdef ISOVERIF
1268	WRITE(,) 'bassin tropics?'
1269	#endif
1270	IF (abs(lat).lt.15.0) THEN
1271	is_in_bassin=.TRUE.
1272	return
1273	endif
1274
1275	ELSE IF (use_bassin_midlats .AND. bassin==bassin_midlats) THEN
1276	#ifdef ISOVERIF
1277	WRITE(,) 'bassin mid lats?'
1278	#endif
1279	IF ((abs(lat).ge.15.0).AND. &
1280	(abs(lat).le.35.0)) THEN
1281	is_in_bassin=.TRUE.
1282	return
1283	endif
1284
1285	else
1286	WRITE(,) 'iso_traceurs_routines 59: bassin inconnu'
1287	WRITE(,) 'bassin_atlantic=' ,bassin_atlantic
1288	WRITE(,) 'bassin_medit=' ,bassin_medit
1289	WRITE(,) 'bassin_indian=' ,bassin_indian
1290	WRITE(,) 'bassin_austral=' ,bassin_austral
1291	WRITE(,) 'bassin_MerArabie=' ,bassin_MerArabie
1292	WRITE(,) 'bassin_BengalGolf=' ,bassin_BengalGolf
1293	WRITE(,) 'bassin_SouthIndian=' ,bassin_SouthIndian
1294	WRITE(,) 'use_bassin_atlantic=' ,use_bassin_atlantic
1295	WRITE(,) 'use_bassin_medit=' ,use_bassin_medit
1296	WRITE(,) 'use_bassin_indian=' ,use_bassin_indian
1297	WRITE(,) 'use_bassin_austral=' ,use_bassin_austral
1298	WRITE(,) 'use_bassin_MerArabie=' ,use_bassin_MerArabie
1299	WRITE(,) 'use_bassin_BengalGolf=' ,use_bassin_BengalGolf
1300	WRITE(,) 'use_bassin_SouthIndian=' ,use_bassin_SouthIndian
1301	stop
1302	endif
1303
1304
1305	END FUNCTION is_in_bassin
1306
1307	SUBROUTINE find_bassin(lat,lon,bassin)
1308	USE isotrac_mod, ONLY: izone_poubelle,ntraceurs_zone=>ntiso,option_traceurs, &
1309	& bassin_map
1310	#ifdef ISOVERIF
1311	USE isotopes_verif_mod
1312	#endif
1313	IMPLICIT NONE
1314
1315	! inputs
1316	REAL lat,lon
1317	! output
1318	INTEGER bassin
1319	!locals
1320	LOGICAL continu
1321	!LOGICAL is_in_bassin
1322
1323	continu=.TRUE.
1324	bassin=1
1325
1326	#ifdef ISOVERIF
1327	WRITE(,) ''
1328	WRITE(,) 'find bassin lat,lon=',lat,lon
1329	#endif
1330	DO while (continu)
1331	!#ifdef ISOVERIF
1332	!! WRITE(,) 'find_bassin 169: lat,lon,bassin=',lat,lon,bassin
1333	!#endif
1334	IF (is_in_bassin(lat,lon,bassin)) THEN
1335	continu=.FALSE.
1336	#ifdef ISOVERIF
1337	WRITE(,) 'find_bassin 173: trouve: bassin=',bassin
1338	#endif
1339	else
1340	!#ifdef ISOVERIF
1341	! WRITE(,) 'find_bassin 175: pas trouve: bassin=',bassin
1342	!#endif
1343	bassin=bassin+1
1344	endif
1345	IF (bassin.EQ.izone_poubelle) THEN
1346	continu=.FALSE.
1347	bassin=izone_poubelle
1348	!#ifdef ISOVERIF
1349	! WRITE(,) 'find_bassin 179: poubelle: bassin=',bassin
1350	!#endif
1351	endif
1352	enddo
1353
1354	! normalement, le bassin est soit un bassin oce, soit un résidu
1355	! donc bassin<=ntraceurs_zone-1
1356	#ifdef ISOVERIF
1357	CALL iso_verif_positif(float(ntraceurs_zone-1-bassin), &
1358	'find_bassin 195')
1359	#endif
1360
1361
1362	END SUBROUTINE find_bassin
1363
1364	SUBROUTINE initialise_bassins_boites(presnivs)
1365	USE dimphy, ONLY: klev
1366	USE lmdz_geometry, ONLY: longitude_deg, latitude_deg
1367	USE isotrac_mod, ONLY: bassin_map,option_traceurs,boite_map
1368	#ifdef ISOVERIF
1369	USE isotopes_verif_mod
1370	#endif
1371	IMPLICIT NONE
1372	REAL presnivs(klev)
1373
1374	IF (option_traceurs.EQ.3) THEN
1375	! initialisation de bassin_map
1376	CALL bassin_map_init(latitude_deg,longitude_deg,bassin_map)
1377	ELSE IF (option_traceurs.EQ.20) THEN
1378	! initialisation de bassin_map selon < ou > 35° lat
1379	WRITE(,) 'physiq 1681: init de la map pour tag 20'
1380	CALL bassin_map_init_opt20(latitude_deg,bassin_map)
1381	ELSE IF (option_traceurs.EQ.5) THEN
1382	CALL boite_AMMA_init(latitude_deg,longitude_deg,presnivs,boite_map)
1383	ELSE IF (option_traceurs.EQ.21) THEN
1384	CALL boite_UT_extra_init(latitude_deg,longitude_deg,presnivs,boite_map)
1385	endif
1386
1387
1388	END SUBROUTINE initialise_bassins_boites
1389
1390	SUBROUTINE bassin_map_init(lat,lon,bassin_map)
1391	USE dimphy, ONLY: klon
1392	#ifdef ISOVERIF
1393	USE isotopes_verif_mod
1394	#endif
1395
1396	IMPLICIT NONE
1397
1398	! inputs
1399	REAL lat(klon),lon(klon)
1400	! output
1401	INTEGER bassin_map(klon)
1402	! locals
1403	INTEGER i
1404
1405	DO i=1,klon
1406	CALL find_bassin(lat(i),lon(i),bassin_map(i))
1407	#ifdef ISOVERIF
1408	WRITE(,) 'init 233: i,lat,lon,bassin=',i,lat(i),lon(i), &
1409	bassin_map(i)
1410	#endif
1411	enddo
1412
1413
1414	END SUBROUTINE bassin_map_init
1415
1416	function is_in_rectangle(x,y,x1,y1,x2,y2)
1417	IMPLICIT NONE
1418	! inputs
1419	REAL x,y
1420	! point en haut à gauche
1421	REAL x1,y1
1422	! point en bas à droite
1423	REAL x2,y2
1424	! output
1425	LOGICAL is_in_rectangle
1426
1427	!#ifdef ISOVERIF
1428	! WRITE(,) 'is_in_rectange 237: x,y=',x,y
1429	! WRITE(,) 'x1,y1,x2,y2=',x1,y1,x2,y2
1430	!#endif
1431	IF ((x-x2.lt.0.1).AND.(x-x1.gt.-0.1).AND. &
1432	(y-y1.lt.0.1).AND.(y-y2.gt.-0.1)) THEN
1433	is_in_rectangle=.TRUE.
1434	else
1435	is_in_rectangle=.FALSE.
1436	endif
1437	!#ifdef ISOVERIF
1438	! WRITE(,) 'is_in_rectangle=',is_in_rectangle
1439	!#endif
1440
1441	END FUNCTION is_in_rectangle
1442
1443	function is_in_triangle(x,y,x1,y1,x2,y2,x3,y3)
1444	IMPLICIT NONE
1445	! inputs
1446	REAL x,y
1447	! points dans le sens trigo
1448	! à gauche
1449	REAL x1,y1
1450	! en bas
1451	REAL x2,y2
1452	! à droite
1453	REAL x3,y3
1454	! output
1455	LOGICAL is_in_triangle
1456	! locals
1457	REAL det1
1458	REAL det2
1459	REAL det3
1460	!#ifdef ISOVERIF
1461	! WRITE(,) 'is_in_triange 271: x,y=',x,y
1462	! WRITE(,) 'x1,y1,x2,y2,x3,y3=',x1,y1,x2,y2,x3,y3
1463	!#endif
1464	det1=(x1-x)(y2-y)-(y1-y)(x2-x)
1465	det2=(x2-x)(y3-y)-(y2-y)(x3-x)
1466	det3=(x3-x)(y1-y)-(y3-y)(x1-x)
1467	is_in_triangle=.FALSE.
1468	IF ((det1det2.gt.0.0).AND.(det2det3.gt.0.0)) THEN
1469	is_in_triangle=.TRUE.
1470	else
1471	is_in_triangle=.FALSE.
1472	endif
1473	!#ifdef ISOVERIF
1474	! WRITE(,) 'det1,det2,det3,is_in_triangle',
1475	! : det1,det2,det3,is_in_triangle
1476	!#endif
1477
1478	END FUNCTION is_in_triangle
1479
1480
1481	SUBROUTINE isotrac_recolorise_tmin(xt,t)
1482	USE dimphy, ONLY: klon, klev
1483	USE isotrac_mod, ONLY: zone_temp,nzone_temp
1484	#ifdef ISOVERIF
1485	USE isotopes_verif_mod
1486	#endif
1487	IMPLICIT NONE
1488
1489	! inout
1490	REAL xt(ntraciso,klon,klev)
1491	! input
1492	REAL t(klon,klev)
1493	! locals
1494	INTEGER izone_temp
1495	INTEGER ixt,ixt_recoit
1496	INTEGER k,i,izone,iiso
1497	!integer find_index
1498
1499
1500	DO k=1,klev
1501	DO i=1,klon
1502	!#ifdef ISOVERIF
1503	! if (i.EQ.1) THEN
1504	! WRITE(,) 'recolorise 396: i,k,t=',i,k,t(i,k)
1505	! WRITE(,) 'zone_temp=',zone_temp
1506	! endif
1507	!#endif
1508	! trouver la zone de cette température
1509	izone_temp=find_index(t(i,k),nzone_temp,zone_temp)
1510
1511	!#ifdef ISOVERIF
1512	! if (i.EQ.1) THEN
1513	! WRITE(,) 'recolorise 414: izone_temp=',izone_temp
1514	! endif
1515	!#endif
1516	DO izone=1,nzone_temp-1
1517	! tous les tags de zone < nzone_temp se trouvant à des
1518	! températures plus basses sont convertis
1519	!#ifdef ISOVERIF
1520	! if (i.EQ.1) THEN
1521	! WRITE(,) 'recolorise 405: izone,xt_eau=',
1522	! : izone,xt(index_trac(izone,iso_eau),i,k)
1523	! endif
1524	!#endif
1525	IF (izone.lt.izone_temp) THEN
1526	DO iiso=1,niso
1527	ixt=index_trac(izone,iiso) ! emmetteur
1528	ixt_recoit=index_trac(izone_temp,iiso) ! recepteur
1529	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k)+xt(ixt,i,k)
1530	xt(ixt,i,k)=0.0
1531	enddo !do iiso=1,niso
1532	endif !if (izone.lt.izone_temp) THEN
1533	!#ifdef ISOVERIF
1534	! if (i.EQ.1) THEN
1535	! WRITE(,) 'recolorise 419: xt_eau,xt_recoit=',
1536	! : xt(index_trac(izone,iso_eau),i,k),
1537	! : xt(index_trac(izone_temp,iso_eau),i,k)
1538	! endif
1539	!#endif
1540	enddo !do izone=1,zone_pot(k)-1
1541
1542	! conversion de l'évap en surf et de la revap des gouttes
1543	DO izone=nzone_temp+1,ntraceurs_zone
1544	DO iiso=1,niso
1545	ixt=index_trac(izone,iiso) ! emmetteur
1546	ixt_recoit=index_trac(izone_temp,iiso) ! recepteur
1547	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k)+xt(ixt,i,k)
1548	xt(ixt,i,k)=0.0
1549	enddo !do iiso=1,niso
1550	enddo !do izone=nzone_temp+1,ntraceurs_zone
1551	enddo !do i=1,klon
1552	enddo !do k=1,klev
1553
1554	#ifdef ISOVERIF
1555	DO k=1,klev
1556	DO i=1,klon
1557	CALL iso_verif_traceur(xt(1,i,k),'recolorise 403')
1558	enddo !do i=1,klon
1559	enddo !do k=1,klev
1560	#endif
1561
1562
1563	END SUBROUTINE isotrac_recolorise_tmin
1564
1565	SUBROUTINE isotrac_recolorise_tmin_sfrev(xt,t)
1566	USE dimphy, ONLY: klon,klev
1567	USE isotrac_mod, ONLY: nzone_temp,zone_temp
1568	#ifdef ISOVERIF
1569	USE isotopes_verif_mod
1570	#endif
1571	IMPLICIT NONE
1572
1573	! recolorise selon la température minimum, mais les tags de
1574	! revap sont laissés en revap
1575
1576	! inout
1577	REAL xt(ntraciso,klon,klev)
1578	! input
1579	REAL t(klon,klev)
1580	! locals
1581	INTEGER izone_temp
1582	INTEGER ixt,ixt_recoit
1583	INTEGER k,i,izone,iiso
1584	!integer find_index
1585
1586	DO k=1,klev
1587	DO i=1,klon
1588	! trouver la zone de cette température
1589	izone_temp=find_index(t(i,k),nzone_temp,zone_temp)
1590
1591	DO izone=1,nzone_temp-1
1592	! tous les tags de zone < nzone_temp se trouvant à des
1593	! températures plus basses sont convertis
1594	! sauf la revap
1595	! le tag de la revap est nzone_temp+1=ntraceurs_zone
1596	IF (izone.lt.izone_temp) THEN
1597	DO iiso=1,niso
1598	ixt=index_trac(izone,iiso) ! emmetteur
1599	ixt_recoit=index_trac(izone_temp,iiso) ! recepteur
1600	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k)+xt(ixt,i,k)
1601	xt(ixt,i,k)=0.0
1602	enddo !do iiso=1,niso
1603	endif !if (izone.lt.izone_temp) THEN
1604	enddo !do izone=1,zone_pot(k)-1
1605
1606	enddo !do i=1,klon
1607	enddo !do k=1,klev
1608
1609	#ifdef ISOVERIF
1610	DO k=1,klev
1611	DO i=1,klon
1612	CALL iso_verif_traceur(xt(1,i,k),'recolorise 594')
1613	enddo !do i=1,klon
1614	enddo !do k=1,klev
1615	#endif
1616
1617
1618	END SUBROUTINE isotrac_recolorise_tmin_sfrev
1619
1620
1621	SUBROUTINE isotrac_recolorise_saturation(xt,rh,lat,pres)
1622	USE dimphy, ONLY: klon,klev
1623	#ifdef ISOVERIF
1624	USE isotopes_verif_mod
1625	#endif
1626	IMPLICIT NONE
1627
1628	! recolorise selon la température minimum, mais les tags de
1629	! revap sont laissés en revap
1630
1631	! inout
1632	REAL xt(ntraciso,klon,klev)
1633	! input
1634	REAL rh(klon,klev)
1635	REAL lat(klon)
1636	REAL pres(klev)
1637	! locals
1638	INTEGER izone_recoit
1639	INTEGER ixt,ixt_recoit
1640	INTEGER k,i,izone,iiso
1641	LOGICAL continu
1642	REAL rh_seuil
1643	parameter (rh_seuil=0.90)
1644	!integer index_zone_latpres
1645
1646	#ifdef ISOVERIF
1647	DO k=1,klev
1648	DO i=1,klon
1649	CALL iso_verif_traceur(xt(1,i,k),'recolorise 612')
1650	enddo !do i=1,klon
1651	enddo !do k=1,klev
1652	#endif
1653
1654	! on ne sature pas les 2 premières couches: on les laisse se
1655	! recharger en evap de surface
1656	DO k=3,klev
1657	DO i=1,klon
1658	IF (rh(i,k).gt.rh_seuil) THEN
1659	izone_recoit=index_zone_latpres(lat(i),pres(k))
1660	DO izone=1,ntraceurs_zone
1661	IF (izone.NE.izone_recoit) THEN
1662	DO iiso=1,niso
1663	ixt=index_trac(izone,iiso) ! emmetteur
1664	ixt_recoit=index_trac(izone_recoit,iiso) ! recepteur
1665	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k)+xt(ixt,i,k)
1666	xt(ixt,i,k)=0.0
1667	enddo !do iiso=1,niso
1668	endif
1669	enddo !do izone=1,ntraceurs_zone
1670	endif !if (rh(i,k).gt.rh_seuil) THEN
1671	enddo !do i=1,klon
1672	enddo !do k=1,klev
1673
1674	#ifdef ISOVERIF
1675	DO k=1,klev
1676	DO i=1,klon
1677	CALL iso_verif_traceur(xt(1,i,k),'recolorise 637')
1678	enddo !do i=1,klon
1679	enddo !do k=1,klev
1680	#endif
1681
1682
1683	END SUBROUTINE isotrac_recolorise_saturation
1684
1685	SUBROUTINE isotrac_recolorise_boite(xt,boite_map)
1686	USE dimphy, ONLY: klon,klev
1687	#ifdef ISOVERIF
1688	USE isotopes_verif_mod
1689	#endif
1690	IMPLICIT NONE
1691
1692	! SUBROUTINE écrite à la base pour tagguer 3 boites AMMA.
1693	! Mais ça peut être générique, selon comment est initialisée boite_map
1694
1695	! inout
1696	REAL xt(ntraciso,klon,klev)
1697	! input
1698	INTEGER boite_map(klon,klev)
1699	! locals
1700	INTEGER i,k
1701	INTEGER izone_recoit,izone,iiso
1702	INTEGER ixt,ixt_recoit
1703
1704	DO k=1,klev
1705	DO i=1,klon
1706	izone_recoit=boite_map(i,k)
1707	IF (izone_recoit.gt.0) THEN
1708	! on est dans une boite connue
1709	! toutes les molécules sont converties en cete couleur
1710	DO izone=1,ntraceurs_zone
1711	IF (izone.NE.izone_recoit) THEN
1712	! on met les traceurs izone dans izone_recoit
1713	DO iiso=1,niso
1714	ixt=index_trac(izone,iiso)
1715	ixt_recoit=index_trac(izone_recoit,iiso)
1716	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k) &
1717	+xt(ixt,i,k)
1718	xt(ixt,i,k)=0.0
1719	enddo
1720	endif !if (izone.NE.izone_recoit) THEN
1721	enddo !do izone=2,ntraceurs_zone
1722	endif !if (izone_recoit.gt.0) THEN
1723	enddo !do i=1,klon
1724	enddo !do k=1,klev
1725
1726	#ifdef ISOVERIF
1727	DO k=1,klev
1728	DO i=1,klon
1729	CALL iso_verif_traceur(xt(1,i,k),'recolorise 514')
1730	enddo !do i=1,klon
1731	enddo !do k=1,klev
1732	#endif
1733
1734
1735	END SUBROUTINE isotrac_recolorise_boite
1736
1737	SUBROUTINE isotrac_recolorise_extra(xt,rlat)
1738	USE dimphy, ONLY: klon,klev
1739	usE isotrac_mod, ONLY: lim_tag20,izone_trop,izone_extra
1740	#ifdef ISOVERIF
1741	USE isotopes_verif_mod
1742	#endif
1743	IMPLICIT NONE
1744
1745	! SUBROUTINE écrite pour l'option de taggage 20
1746	! permet de retagguer la vapeur tropicale en vapeur
1747	! extratropicale dès qu'elle atteint 35° de latitude
1748
1749	! inout
1750	REAL xt(ntraciso,klon,klev)
1751	! input
1752	REAL rlat(klon)
1753	! locals
1754	INTEGER i,k
1755	INTEGER iiso,ixt,ixt_recoit
1756
1757	! WRITE(,) 'iso_traceurs_routines 723: lim_tag20=',lim_tag20
1758	DO k=1,klev
1759	DO i=1,klon
1760	IF (abs(rlat(i)).gt.lim_tag20) THEN
1761	! on met les traceurs izone_trop dans izone_extra
1762	DO iiso=1,niso
1763	ixt=index_trac(izone_trop,iiso)
1764	ixt_recoit=index_trac(izone_extra,iiso)
1765	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k) &
1766	+xt(ixt,i,k)
1767	xt(ixt,i,k)=0.0
1768	enddo
1769	endif ! if (abs(rlat(i)).lt.35.0) THEN
1770	enddo !do i=1,klon
1771	enddo !do k=1,klev
1772
1773	#ifdef ISOVERIF
1774	DO k=1,klev
1775	DO i=1,klon
1776	CALL iso_verif_traceur(xt(1,i,k),'recolorise 741')
1777	enddo !do i=1,klon
1778	enddo !do k=1,klev
1779	#endif
1780
1781
1782	END SUBROUTINE isotrac_recolorise_extra
1783
1784	SUBROUTINE isotrac_recolorise_conv(xt,rlat,presnivs,rain_con)
1785	USE dimphy, ONLY: klon,klev
1786	USE isotrac_mod, ONLY: lim_precip_tag22, &
1787	& izone_conv_BT,izone_conv_UT
1788	#ifdef ISOVERIF
1789	USE isotopes_verif_mod
1790	#endif
1791	IMPLICIT NONE
1792
1793	! SUBROUTINE écrite pour l'option de taggage 20
1794	! permet de retagguer la vapeur tropicale en vapeur
1795	! extratropicale dès qu'elle atteint 35° de latitude
1796
1797	! inout
1798	REAL xt(ntraciso,klon,klev)
1799	! input
1800	REAL rlat(klon)
1801	REAL presnivs(klev)
1802	REAL rain_con(klon)
1803	! locals
1804	INTEGER i,k
1805	INTEGER iiso,ixt,ixt_recoit,izone
1806
1807	! WRITE(,) 'iso_traceurs_routines 723: lim_tag20=',lim_tag20
1808	! WRITE(,) 'presnivs=',presnivs
1809	! stop
1810	DO k=1,klev
1811	DO i=1,klon
1812	#ifdef ISOVERIF
1813	IF ((abs(rlat(i)).lt.30.0).AND.(k.EQ.1)) THEN
1814	endif
1815	#endif
1816	IF ((abs(rlat(i)).lt.30.0).AND. &
1817	(rain_con(i)*86400.gt.lim_precip_tag22)) THEN
1818	! on met les traceurs izone_trop dans izone_conv fn z
1819	DO iiso=1,niso
1820	IF (presnivs(k).gt.650.0*100.0) THEN
1821	ixt_recoit=index_trac(izone_conv_BT,iiso)
1822	else
1823	ixt_recoit=index_trac(izone_conv_UT,iiso)
1824	endif
1825	DO izone=1,ntraceurs_zone
1826	ixt=index_trac(izone,iiso)
1827	IF (ixt.NE.ixt_recoit) THEN
1828	xt(ixt_recoit,i,k)=xt(ixt_recoit,i,k) &
1829	+xt(ixt,i,k)
1830	xt(ixt,i,k)=0.0
1831	endif !if (ixt.NE.ixt_recoit) THEN
1832	enddo !do izone=1,ntraceurs_zone
1833	enddo !do iiso=1,niso
1834	! WRITE(,) 'k,presnivs,ixt,ixt_recoit=',k,presnivs(k),
1835	! : ixt,ixt_recoit
1836	! WRITE(,) 'xt(:,i,k)=',xt(:,i,k)
1837	endif ! if (abs(rlat(i)).lt.35.0) THEN
1838	enddo !do i=1,klon
1839	enddo !do k=1,klev
1840
1841	#ifdef ISOVERIF
1842	DO k=1,klev
1843	DO i=1,klon
1844	CALL iso_verif_traceur(xt(1,i,k),'recolorise 741')
1845	enddo !do i=1,klon
1846	enddo !do k=1,klev
1847	#endif
1848
1849
1850	END SUBROUTINE isotrac_recolorise_conv
1851
1852
1853	SUBROUTINE boite_AMMA_init(lat,lon,presnivs,boite_map)
1854	USE dimphy, ONLY: klon,klev
1855	#ifdef ISOVERIF
1856	USE isotopes_verif_mod
1857	#endif
1858	USE isotrac_mod, ONLY: izone_aej,izone_mousson,izone_harmattan
1859	IMPLICIT NONE
1860
1861	REAL lat(klon),lon(klon)
1862	REAL presnivs(klev)
1863	! output
1864	INTEGER boite_map(klon,klev)
1865	! locals
1866	INTEGER i,k
1867
1868	! WRITE(,) 'izone_aej,izone_mousson,izone_harmattan=',
1869	! : izone_aej,izone_mousson,izone_harmattan
1870	DO k=1,klev
1871	DO i=1,klon
1872	boite_map(i,k)=0.0
1873	! WRITE(,) 'i,k,lat,lon,pres=',
1874	! : i,k,lat(i),lon(i),presnivs(k)
1875	IF ((presnivs(k).le.700.0*100.0).AND. &
1876	(presnivs(k).gt.400.0*100.0).AND. &
1877	(lat(i).gt.8.0).AND. &
1878	(lat(i).lt.20.0).AND. &
1879	(lon(i).gt.10.0).AND. &
1880	(lon(i).lt.30.0)) THEN
1881	boite_map(i,k)=izone_aej
1882	! WRITE(,) ' -> zone AEJ'
1883	ELSE IF ((presnivs(k).ge.850.0*100.0).AND. &
1884	(lat(i).gt.-5.0).AND. &
1885	(lat(i).le.8.0).AND. &
1886	(lon(i).gt.-40.0).AND. &
1887	(lon(i).lt.15.0)) THEN
1888	boite_map(i,k)=izone_mousson
1889	! WRITE(,) ' -> zone flux de mousson'
1890	ELSE IF ((presnivs(k).gt.700.0*100.0).AND. &
1891	(lat(i).ge.20.0).AND. &
1892	(lat(i).lt.30.0).AND. &
1893	(lon(i).gt.-10.0).AND. &
1894	(lon(i).lt.40.0)) THEN
1895	boite_map(i,k)=izone_harmattan
1896	! WRITE(,) ' -> zone Harmattan'
1897	endif
1898	! WRITE(,) ' ** boite_map=',boite_map(i,k)
1899	enddo
1900	enddo
1901
1902
1903	END SUBROUTINE boite_AMMA_init
1904
1905
1906	SUBROUTINE boite_UT_extra_init(lat,lon,presnivs,boite_map)
1907	USE dimphy, ONLY: klon,klev
1908	USE isotrac_mod, ONLY: izone_extra,izone_trop
1909	#ifdef ISOVERIF
1910	USE isotopes_verif_mod
1911	#endif
1912	IMPLICIT NONE
1913
1914	REAL lat(klon),lon(klon)
1915	REAL presnivs(klev)
1916	! output
1917	INTEGER boite_map(klon,klev)
1918	! locals
1919	INTEGER i,k
1920
1921	! WRITE(,) 'izone_trop,izone_extra=',
1922	! : izone_trop,izone_extra
1923	DO k=1,klev
1924	DO i=1,klon
1925	boite_map(i,k)=0.0
1926	! WRITE(,) 'i,k,lat,lon,pres=',
1927	! : i,k,lat(i),lon(i),presnivs(k)
1928	IF ((presnivs(k).le.500.0*100.0) &
1929	.AND.(abs(lat(i)).lt.15.0)) THEN
1930	boite_map(i,k)=izone_trop
1931	! WRITE(,) ' -> zone trop'
1932	ELSE IF (abs(lat(i)).gt.35.0) THEN
1933	boite_map(i,k)=izone_extra
1934	! WRITE(,) ' -> zone extratropiques'
1935	endif
1936	! WRITE(,) ' ** boite_map=',boite_map(i,k)
1937	enddo
1938	enddo
1939
1940
1941	END SUBROUTINE boite_UT_extra_init
1942
1943
1944	function index_zone_lat(lat)
1945	USE isotrac_mod, ONLY: lattag_min,dlattag,nzone_lat
1946	IMPLICIT NONE
1947
1948	! inputs
1949	REAL lat,pres
1950	! output
1951	INTEGER index_zone_lat
1952
1953	IF (lat.lt.lattag_min) THEN
1954	index_zone_lat=1
1955	else
1956	index_zone_lat=int((lat-lattag_min)/dlattag)+2
1957	index_zone_lat=min(index_zone_lat,nzone_lat)
1958	endif
1959
1960
1961	END FUNCTION index_zone_lat
1962
1963
1964	function index_zone_pres(pres)
1965	USE isotrac_mod, ONLY: nzone_pres,zone_pres
1966	IMPLICIT NONE
1967
1968	! inputs
1969	REAL lat,pres
1970	! output
1971	INTEGER index_zone_pres
1972	!integer find_index
1973
1974
1975	index_zone_pres=find_index(pres,nzone_pres,zone_pres)
1976	WRITE(,) 'iso_traceurs_routines 802: pres,index_zone_pres=', &
1977	pres,index_zone_pres
1978	WRITE(,) 'zone_pres=',zone_pres(1:nzone_pres-1)
1979
1980
1981	END FUNCTION index_zone_pres
1982
1983	function find_index(pres,nzone_pres,zone_pres)
1984	IMPLICIT NONE
1985
1986	! inputs
1987	REAL pres
1988	INTEGER nzone_pres
1989	REAL zone_pres(nzone_pres)
1990	! output
1991	INTEGER find_index
1992	LOGICAL continu
1993
1994	IF (nzone_pres.gt.1) THEN
1995	IF (pres.ge.zone_pres(1)) THEN
1996	find_index=1
1997	ELSE IF (pres.lt.zone_pres(nzone_pres-1)) THEN
1998	find_index=nzone_pres
1999	else !if (t(i,k).ge.zone_temp1) THEN
2000	continu=.TRUE.
2001	find_index=2
2002	DO while (continu)
2003	IF (pres.ge.zone_pres(find_index)) THEN
2004	continu=.FALSE.
2005	! c'est izone_temp, zone trouvée
2006	else
2007	find_index=find_index+1
2008	endif
2009	enddo !do while (continu)
2010	endif !if (t(i,k).ge.zone_temp1) THEN
2011	else !if (nzone_pres.gt.1) THEN
2012	find_index=1
2013	endif !if (nzone_pres.gt.1) THEN
2014	END FUNCTION find_index
2015
2016	function index_zone_latpres(lat,pres)
2017	USE isotrac_mod, ONLY: nzone_lat
2018	IMPLICIT NONE
2019
2020	! inputs
2021	REAL lat,pres
2022	! output
2023	INTEGER index_zone_latpres
2024	! locals
2025	INTEGER index_lat
2026	INTEGER index_pres
2027	!integer index_zone_lat
2028	!integer index_zone_pres
2029
2030	index_lat=index_zone_lat(lat)
2031	index_pres=index_zone_pres(pres)
2032	index_zone_latpres=index_lat+(index_pres-1)*nzone_lat
2033
2034
2035	END FUNCTION index_zone_latpres
2036
2037	SUBROUTINE iso_recolorise_condensation(qt,cond, &
2038	xt,zxtcond,tcond,ep,xtres, &
2039	seuil_in)
2040	USE dimphy, ONLY: klon,klev
2041	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau
2042	USE isotrac_mod, ONLY: option_seuil_tag_tmin,izone_cond, &
2043	& nzone_temp,zone_temp
2044	#ifdef ISOVERIF
2045	USE isotopes_verif_mod
2046	#endif
2047	IMPLICIT NONE
2048
2049	! on recolorise la vapeur résiduelle selon la température de condensation
2050	! on supose qu'une vapeur xt,q condense en cond,zxtcond, à une
2051	! température tcond. A ce stade, la vapeur initiale n'est pas
2052	! retranchée de son condensat. On calcule les tags dans la vepru
2053	! résiduelle xtres qu'on aurait si on retranchait un fraction ep du
2054	! condensat
2055
2056	! inputs
2057	REAL qt
2058	REAL cond
2059	REAL tcond
2060	REAL ep
2061	REAL xt(ntraciso)
2062	REAL zxtcond(ntraciso)
2063	REAL seuil_in
2064	! outputs
2065	REAL xtres(ntraciso)
2066	! locals
2067	INTEGER izone_temp,izone
2068	INTEGER ixt,ixt_recoit
2069	INTEGER iiso
2070	!integer find_index
2071	REAL fcond, qmicro
2072	! real f
2073
2074	IF ((cond.gt.0.0).AND.(qt.gt.0.0)) THEN
2075	izone_temp=find_index(tcond,nzone_temp,zone_temp)
2076	! WRITE(,) 'pgm 901 tmp: izone_temp=',izone_temp
2077	#ifdef ISOVERIF
2078	DO ixt=1,ntraciso
2079	CALL iso_verif_positif(xt(ixt)-zxtcond(ixt), &
2080	'iso_trac 898')
2081	enddo !do ixt=1,ntraciso
2082	CALL iso_verif_traceur_justmass(xt, &
2083	'iso_trac_routines 906')
2084	CALL iso_verif_traceur_justmass(zxtcond, &
2085	'iso_trac_routines 908')
2086	#endif
2087	! bidouille
2088	IF (bidouille_anti_divergence) THEN
2089	CALL iso_verif_traceur_jbidouille(xt)
2090	CALL iso_verif_traceur_jbidouille(zxtcond)
2091	endif
2092
2093	DO ixt=1,niso
2094	xtres(ixt)=xt(ixt)-ep*zxtcond(ixt)
2095	enddo
2096	DO ixt=1+niso,ntraciso
2097	xtres(ixt)=0.0
2098	enddo
2099	! WRITE(,) 'iso_trac_routines tmp 916: xtres=',xtres
2100
2101	#ifdef ISOVERIF
2102	DO ixt=1,ntraciso
2103	CALL iso_verif_positif(xtres(ixt), &
2104	'iso_trac_routines 921')
2105	enddo
2106	#endif
2107
2108	! cas de izone sfc et izone precip et izone cond et izone< izone_temp
2109
2110	! WRITE(,) 'iso_trac 940: cond/qt,seuil_in,izone_temp=',
2111	! : cond/qt,seuil_in,izone_temp
2112
2113	IF (option_seuil_tag_tmin.EQ.2) THEN
2114	qmicro=0.0
2115	DO izone=nzone_temp+1,ntraceurs_zone
2116	ixt= index_trac(izone,iso_eau)
2117	qmicro=qmicro+xt(ixt)
2118	enddo !do izone=nzone_temp+1,ntraceurs_zone
2119	IF (qt-qmicro.gt.0.0) THEN
2120	fcond=(cond-qmicro)/(qt-qmicro)
2121	else
2122	fcond=0.0
2123	endif
2124	else
2125	fcond=cond/qt
2126	endif
2127
2128	IF (fcond.gt.seuil_in) THEN
2129	! on les transfert à izone_temp
2130	DO izone=1,ntraceurs_zone
2131	IF ((izone.gt.nzone_temp).OR.(izone.lt.izone_temp)) THEN
2132	! ieau=index_trac(izone,iso_eau)
2133	DO iiso=1,niso
2134	ixt= index_trac(izone,iiso)
2135	ixt_recoit=index_trac(izone_temp,iiso) ! recepteur
2136	xtres(ixt_recoit)=xtres(ixt_recoit) &
2137	+(xt(ixt)-zxtcond(ixt))
2138	xtres(ixt)=0.0
2139	! WRITE(,) 'iso_trac 920: izone,ixt,',
2140	! : 'ixt_recoit=',
2141	! : izone,ixt,ixt_recoit
2142	! WRITE(,) 'isotrac 924: xt=',xt
2143	! WRITE(,) 'isotrac 925: zxtcond=',zxtcond
2144	enddo !do iiso=1,niso
2145	! WRITE(,) 'iso_trac tmp 944: izone,xtres=',izone,xtres
2146	endif !if (izone.NE.izone_cond) THEN
2147	enddo !do izone=nzones_temp+1,ntraceurs_zone
2148
2149	else !if (cond/qt.gt.seuil_in) THEN
2150	! on les laisse sur place
2151	DO izone=1,ntraceurs_zone
2152	IF ((izone.gt.nzone_temp).OR.(izone.lt.izone_temp)) THEN
2153	DO iiso=1,niso
2154	ixt= index_trac(izone,iiso)
2155	xtres(ixt)=(xt(ixt)-zxtcond(ixt))
2156	enddo !do iiso=1,niso
2157	endif !if (izone.NE.izone_cond) THEN
2158	enddo !do izone=nzones_temp+1,ntraceurs_zone
2159
2160	endif !if (cond/qt.gt.seuil_in) THEN
2161	! izone_temp est conservé, on lui enlève juste son
2162	! condesat
2163	DO iiso=1,niso
2164	ixt_recoit=index_trac(izone_temp,iiso) ! recepteur
2165	xtres(ixt_recoit)=xtres(ixt_recoit) &
2166	+(xt(ixt_recoit)-zxtcond(ixt_recoit))
2167	enddo !do iiso=1,niso
2168
2169	#ifdef ISOVERIF
2170	DO ixt=1,ntraciso
2171	CALL iso_verif_positif(xtres(ixt), &
2172	'iso_trac_routines 940')
2173	enddo
2174	#endif
2175
2176	! cas des zones > izone temp
2177	! on conserve le condensat résiduel
2178	DO izone=izone_temp+1,nzone_temp
2179	DO iiso=1,niso
2180	ixt= index_trac(izone,iiso)
2181	xtres(ixt)=xt(ixt)-zxtcond(ixt)
2182	! WRITE(,) 'iso_trac 931: izone,ixt,ixt_recoit=',
2183	! : izone,ixt,ixt_recoit
2184	! WRITE(,) 'isotrac 934: xt=',xt
2185	! WRITE(,) 'isotrac 935: zxtice=',zxtice
2186	enddo !do iiso=1,niso
2187	! WRITE(,) 'iso_trac tmp 965: izone,xtres=',izone,xtres
2188	enddo !do izone=izone_temp+1,nzones_temp
2189
2190	! on rajoute le condensat qui ne precipite pas
2191	IF (ep.lt.1.0) THEN
2192	DO iiso=1,niso
2193	ixt= index_trac(izone_cond,iiso)
2194	xtres(ixt)=xtres(ixt)+(1.0-ep)*zxtcond(iiso)
2195	! WRITE(,) 'iso_trac 940: izone,ixt,ixt_recoit=',
2196	! : izone,ixt,ixt_recoit
2197	! WRITE(,) 'isotrac 1014: xt=',xt
2198	! WRITE(,) 'isotrac 945: zxtice=',zxtice
2199	enddo !do iiso=1,niso
2200	endif !if (ep.lt.0.0) THEN
2201	else
2202	! si cond=0 ou qt=0, tot reste pareil
2203	DO ixt=1,ntraciso
2204	xtres(ixt)=xt(ixt)
2205	enddo !do ixt=1,ntraciso
2206
2207	endif ! if (qt.gt.0.0) THEN
2208	#ifdef ISOVERIF
2209	IF (iso_verif_traceur_jm_nostop(xtres, &
2210	'iso_trac_routines 166').EQ.1) THEN
2211	WRITE(,) 'isotrac 1024: xt=',xt
2212	WRITE(,) 'zxtcond=',zxtcond
2213	WRITE(,) 'xtres=',xtres
2214	WRITE(,) 'ep=',ep
2215	stop
2216	endif
2217	#endif
2218	#ifdef ISOVERIF
2219	DO ixt=1,ntraciso
2220	CALL iso_verif_positif(xtres(ixt),'iso_trac_routines 953')
2221	enddo
2222	IF (nzone_temp.ge.5) THEN
2223	IF (iso_verif_tag17_q_deltaD_chns(xtres, &
2224	'iso_trac_routines 1025').EQ.1) THEN
2225	WRITE(,) 'xt=',xt
2226	WRITE(,) 'zxtcond=',zxtcond
2227	WRITE(,) 'xtres=',xtres
2228	WRITE(,) 'ep=',ep
2229	WRITE(,) 'tcond=',tcond
2230	WRITE(,) 'izone_temp=',izone_temp
2231	stop
2232	endif
2233	endif
2234	! WRITE(,) 'isotrac 1048: sortie de iso_recolorise_condensation'
2235	#endif
2236
2237
2238	END SUBROUTINE iso_recolorise_condensation
2239
2240	SUBROUTINE bassin_map_init_opt20(lat,bassin_map)
2241	USE dimphy, ONLY: klon
2242	USE isotrac_mod, ONLY: izone_cont,izone_trop,lim_tag20
2243	#ifdef ISOVERIF
2244	USE isotopes_verif_mod
2245	#endif
2246	IMPLICIT NONE
2247
2248	! inputs
2249	REAL lat(klon)
2250	! output
2251	INTEGER bassin_map(klon)
2252	! locals
2253	INTEGER i
2254
2255	WRITE(,) 'iso_traceurs_routines 1142: lim_tag20=',lim_tag20
2256	DO i=1,klon
2257	IF (abs(lat(i)).gt.lim_tag20) THEN
2258	bassin_map(i)=izone_cont
2259	else
2260	bassin_map(i)=izone_trop
2261	endif
2262	enddo !do i=1,klon
2263
2264
2265	END SUBROUTINE bassin_map_init_opt20
2266
2267	SUBROUTINE isotrac_recolorise_general(xt_seri,t_seri,zx_rh,presnivs)
2268	USE lmdz_geometry, ONLY: latitude_deg
2269	USE dimphy, ONLY: klon,klev
2270	USE isotrac_mod, ONLY: option_traceurs,boite_map
2271	IMPLICIT NONE
2272
2273	! inputs
2274	REAL, DIMENSION(ntraciso,klon,klev), INTENT(IN) :: xt_seri
2275	REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri
2276	REAL, DIMENSION(klon,klev), INTENT(IN) :: zx_rh
2277	REAL, DIMENSION(klev), INTENT(IN) :: presnivs
2278
2279	IF (option_traceurs.EQ.4) THEN
2280	CALL isotrac_recolorise_tmin(xt_seri,t_seri)
2281	elseif ((option_traceurs.EQ.5).OR. &
2282	(option_traceurs.EQ.21)) THEN
2283	CALL isotrac_recolorise_boite(xt_seri,boite_map)
2284	elseif (option_traceurs.EQ.13) THEN
2285	CALL isotrac_recolorise_tmin_sfrev(xt_seri,t_seri)
2286	elseif (option_traceurs.EQ.14) THEN
2287	CALL isotrac_recolorise_saturation(xt_seri,zx_rh,latitude_deg,presnivs)
2288	elseif (option_traceurs.EQ.20) THEN
2289	CALL isotrac_recolorise_extra(xt_seri,latitude_deg)
2290	endif !if (option_traceurs.EQ.4) THEN
2291	END SUBROUTINE isotrac_recolorise_general
2292
2293
2294
2295
2296	SUBROUTINE iso_verif_traceur_jbid_vect(x,n,m)
2297	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2298	!use isotrac_mod, ONLY: ntraceurs_zone=>nzone
2299	USE infotrac_phy, ONLY: ntraceurs_zone=>nzone
2300	IMPLICIT NONE
2301
2302	! version vectrisée de iso_verif_traceur_jbidouille
2303	! inputs
2304	INTEGER n,m
2305	REAL x(ntraciso,n,m)
2306
2307	! locals
2308	INTEGER iiso,izone,ixt,i,j
2309	REAL xtractot(n,m)
2310
2311	IF (bidouille_anti_divergence) THEN
2312	DO iiso=1,niso
2313
2314	DO j=1,m
2315	DO i=1,n
2316	xtractot(i,j)=0.0
2317	enddo !do j=1,m
2318	enddo !do j=1,m
2319
2320	DO izone=1,ntraceurs_zone
2321	ixt=index_trac(izone,iiso)
2322	DO j=1,m
2323	DO i=1,n
2324	xtractot(i,j)=xtractot(i,j)+x(ixt,i,j)
2325	enddo !do j=1,m
2326	enddo !do j=1,m
2327	enddo !do izone=1,ntraceurs_zone
2328
2329	! on réajuste pour que les traceurs fasses bien la somme
2330	! des traceurs
2331	DO izone=1,ntraceurs_zone
2332	ixt=index_trac(izone,iiso)
2333	DO j=1,m
2334	DO i=1,n
2335	! if (abs(xtractot(i,j)).gt.ridicule*10) THEN
2336	IF (abs(xtractot(i,j)).gt.ridicule) THEN
2337	! modif le 19 fev 2011
2338	x(ixt,i,j)=x(ixt,i,j)/xtractot(i,j)*x(iiso,i,j)
2339	endif !if (abs(xtractot(i,j)).gt.ridicule*10) THEN
2340	enddo !do i=1,n
2341	enddo !do j=1,m
2342	enddo !do izone=1,ntraceurs_zone
2343
2344	! ! ajout le 19 fev 2011
2345	! ! on rend plutot les vérifs plus strictes
2346	! ixt=index_trac(izone_poubelle,iiso)
2347	! do j=1,m
2348	! do i=1,n
2349	! if ((abs(xtractot(i,j)).lt.1e-18).AND.
2350	! : (x(iiso,i,j).gt.ridicule)) THEN
2351	! x(ixt,i,j)=x(iiso,i,j)
2352	! endif !if (abs(xtractot(i,j)).gt.ridicule*10) THEN
2353	! enddo ! do i=1,n
2354	! enddo !do j=1,m
2355
2356	enddo !do iiso=1,ntraceurs_iso
2357	endif !if (bidouille_anti_divergence) THEN
2358	END SUBROUTINE iso_verif_traceur_jbid_vect
2359
2360	SUBROUTINE iso_verif_traceur_jbidouille(x)
2361	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2362	IMPLICIT NONE
2363
2364	! on réajuste aussi les valeurs des traceurs pour la
2365	! conservation de la masse, dans le cas bidouille
2366
2367	! inputs
2368	REAL x(ntraciso)
2369
2370	! locals
2371	INTEGER iiso,izone,ixt
2372	REAL xtractot
2373
2374	IF (bidouille_anti_divergence) THEN
2375	DO iiso=1,niso
2376
2377	xtractot=0.0
2378	DO izone=1,ntraceurs_zone
2379	ixt=index_trac(izone,iiso)
2380	xtractot=xtractot+x(ixt)
2381	enddo !do izone=1,ntraceurs_zone
2382
2383	! on réajuste pour que les traceurs fasses bien la somme
2384	! des traceurs
2385	IF (abs(xtractot).gt.ridicule*10) THEN
2386	DO izone=1,ntraceurs_zone
2387	ixt=index_trac(izone,iiso)
2388	x(ixt)=x(ixt)/xtractot*x(iiso)
2389	enddo !do izone=1,ntraceurs_zone
2390	endif
2391
2392	enddo !do iiso=1,ntraceurs_iso
2393	endif !if (bidouille_anti_divergence) THEN
2394	END SUBROUTINE iso_verif_traceur_jbidouille
2395
2396
2397	SUBROUTINE iso_verif_traceur_jbid_pos(x)
2398	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2399	!#ifdef ISOVERIF
2400	! use isotopes_verif_mod, ONLY: iso_verif_traceur_pbidouille
2401	!#endif
2402	IMPLICIT NONE
2403
2404	! on réajuste les valeurs des traceurs pour qu'il n'y ai pas de
2405	! valeurs négatives. Si valeurs négatives -> on pompe les autres
2406	! traceurs
2407	! attention: fait la même chose pour tous les isos -> peut
2408	! induire des fractionnements.
2409	! Pour ne pas induire des fractionnements, prendre
2410	! iso_verif_traceur_jbid_pos2
2411	! avantage de cette subroutine: conserve la masse en isotopes
2412	! légers, ce qui nest pas le cas de pos2
2413
2414	! inputs
2415	REAL x(ntraciso)
2416
2417	! locals
2418	INTEGER iiso,izone,ixt
2419	REAL xtractot,xtractotprec
2420
2421	IF (bidouille_anti_divergence) THEN
2422	! WRITE(,) 'pgm 532 tmp: x=',x
2423	DO iiso=1,niso
2424
2425	xtractot=0.0
2426	xtractotprec=0.0
2427	DO izone=1,ntraceurs_zone
2428	ixt=index_trac(izone,iiso)
2429	xtractotprec=xtractotprec+x(ixt)
2430	x(ixt)=max(x(ixt),0.0)
2431	xtractot=xtractot+x(ixt)
2432	enddo !do izone=1,ntraceurs_zone
2433	! WRITE(,) 'iiso,xtractotprec,xtractot=',
2434	! : iiso,xtractotprec,xtractot
2435
2436	IF (xtractot.gt.xtractotprec) THEN
2437	! on réajuste pour que les traceurs fasses bien la somme
2438	! des traceurs
2439	IF (abs(xtractot).gt.ridicule) THEN
2440	DO izone=1,ntraceurs_zone
2441	ixt=index_trac(izone,iiso)
2442	x(ixt)=x(ixt)*xtractotprec/xtractot
2443	enddo !do izone=1,ntraceurs_zone
2444	! on modifie aussi l'isotope de base si lui aussi était
2445	! négatif
2446	! x(iiso)=xtractot
2447	else !if (abs(xtractot).gt.ridicule) THEN
2448	! normallement, valeurs restantes très faibles
2449	! on ne fait rien.
2450	! on met juste un max
2451	x(iiso)=max(x(iiso),0.0)
2452	DO izone=1,ntraceurs_zone
2453	ixt=index_trac(izone,iiso)
2454	x(ixt)=max(x(ixt),0.0)
2455	enddo !do izone=1,ntraceurs_zone
2456	endif !if (abs(xtractot).gt.ridicule) THEN
2457	endif !if (xtractot.gt.xtractotprec) THEN
2458	enddo !do iiso=1,ntraceurs_iso
2459	#ifdef ISOVERIF
2460	CALL iso_verif_traceur_pbidouille(x,'iso_verif_trac 558')
2461	#else
2462	CALL iso_verif_traceur_jbidouille(x)
2463	#endif
2464	endif !if (bidouille_anti_divergence) THEN
2465	END SUBROUTINE iso_verif_traceur_jbid_pos
2466
2467	SUBROUTINE iso_verif_traceur_jbid_pos_vect(n,m,x)
2468	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2469	#ifdef ISOVERIF
2470	USE isotopes_verif_mod
2471	#endif
2472	IMPLICIT NONE
2473
2474	! inputs
2475	INTEGER n,m
2476	REAL x(ntraciso,n,m)
2477
2478	! locals
2479	INTEGER iiso,izone,ixt
2480	REAL xtractot(n,m),xtractotprec(n,m)
2481	INTEGER i,j
2482
2483	IF (bidouille_anti_divergence) THEN
2484	! WRITE(,) 'pgm 532 tmp: x=',x
2485
2486	DO iiso=1,niso
2487	DO j=1,m
2488	DO i=1,n
2489	xtractot(i,j)=0.0
2490	xtractotprec(i,j)=0.0
2491	enddo !do j=1,m
2492	enddo !do i=1,n
2493	DO izone=1,ntraceurs_zone
2494	ixt=index_trac(izone,iiso)
2495
2496	DO j=1,m
2497	DO i=1,n
2498	xtractotprec(i,j)=xtractotprec(i,j)+x(ixt,i,j)
2499	x(ixt,i,j)=max(x(ixt,i,j),0.0)
2500	xtractot(i,j)=xtractot(i,j)+x(ixt,i,j)
2501	enddo !do i=1,n
2502	enddo !do j=1,m
2503
2504	enddo !do izone=1,ntraceurs_zone
2505	! WRITE(,) 'iiso,xtractotprec,xtractot=',
2506	! : iiso,xtractotprec,xtractot
2507
2508	DO j=1,m
2509	DO i=1,n
2510	IF (xtractot(i,j).gt.xtractotprec(i,j)) THEN
2511	! on réajuste pour que les traceurs fasses bien la somme
2512	! des traceurs
2513	IF (abs(xtractot(i,j)).gt.ridicule) THEN
2514	DO izone=1,ntraceurs_zone
2515	ixt=index_trac(izone,iiso)
2516	x(ixt,i,j)=x(ixt,i,j)*xtractotprec(i,j)/xtractot(i,j)
2517	enddo !do izone=1,ntraceurs_zone
2518	! on modifie aussi l'isotope de base si lui aussi était
2519	! négatif
2520	! x(iiso)=xtractot
2521	else !if (abs(xtractot).gt.ridicule) THEN
2522	! normallement, valeurs restantes très faibles
2523	! on ne fait rien.
2524	! on met juste un max
2525	x(iiso,i,j)=max(x(iiso,i,j),0.0)
2526	DO izone=1,ntraceurs_zone
2527	ixt=index_trac(izone,iiso)
2528	x(ixt,i,j)=max(x(ixt,i,j),0.0)
2529	enddo !do izone=1,ntraceurs_zone
2530	endif !if (abs(xtractot).gt.ridicule) THEN
2531	endif !if (xtractot.gt.xtractotprec) THEN
2532	enddo !do i=1,n
2533	enddo !do j=1,m
2534
2535
2536	enddo !do iiso=1,ntraceurs_iso
2537	#ifdef ISOVERIF
2538	CALL iso_verif_traceur_pbid_vect(x,n,m,'iso_verif_trac 558')
2539	#else
2540	CALL iso_verif_traceur_jbid_vect(x,n,m)
2541	#endif
2542	endif !if (bidouille_anti_divergence) THEN
2543	END SUBROUTINE iso_verif_traceur_jbid_pos_vect
2544
2545	SUBROUTINE iso_verif_traceur_jbid_pos2(x,q)
2546	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2547	#ifdef ISOVERIF
2548	USE isotopes_verif_mod
2549	#endif
2550	IMPLICIT NONE
2551
2552	! même but que iso_verif_traceur_jbid_pos, mais n'induit
2553	! pas de fractionnement.
2554	! on regarde si xteau est positif. S'il ne l'est pas, on pompe
2555	! dans les autres tags pour le mettre à 0. On conserve la compo
2556	! iso.
2557	! Pb: ne conserve pas la masse d'isotopes légers.
2558
2559	! inputs
2560	REAL x(ntraciso),q
2561
2562	! locals
2563	INTEGER iiso,izone,ixt,ieau
2564	REAL dqtmp,factmp
2565
2566	IF (bidouille_anti_divergence) THEN
2567	! WRITE(,) 'iso_verif_trac 578 tmp: q,xt=',
2568	! : q,x(1:ntraciso)
2569	IF (q.gt.0.0) THEN
2570	dqtmp=0.0
2571	DO izone=1,ntraceurs_zone
2572	ieau=index_trac(izone,iso_eau)
2573	IF (x(ieau).lt.0.0) THEN
2574	! WRITE(,) 'local_x<0 pour izone=',izone
2575	dqtmp=dqtmp-x(ieau)
2576	DO iiso=1,niso
2577	ixt=index_trac(izone,iiso)
2578	x(ixt) =0.0
2579	enddo !do iiso=1,niso
2580	endif !if (local_xt(ieau,i,k).lt.0.0) THEN
2581	enddo !do izone=1,ntraceurs_zone
2582	! WRITE(,) 'dqtmp=',dqtmp
2583	IF (dqtmp.gt.0.0) THEN
2584	! WRITE(,) 'iso_verif_trac 593 warning: q,dqtmp,xt=',
2585	! : q,dqtmp,x(1:ntraciso)
2586	! on redistribue la négativité des traceurs dans les
2587	! traceurs positifs
2588	! factmp=(1.0-dqtmp/(local_q(i,k)+dqtmp))
2589	! correction janv 2010
2590	factmp=(q/(q+dqtmp))
2591	! WRITE(,) 'factmp=',factmp
2592	DO izone=1,ntraceurs_zone
2593	ieau=index_trac(izone,iso_eau)
2594	IF (x(ieau).gt.0.0) THEN
2595	DO iiso=1,niso
2596	ixt=index_trac(izone,iiso)
2597	x(ixt)=x(ixt)*factmp
2598	enddo !do iiso=1,niso
2599	endif !if (local_xt(ieau,i,k).gt.0.0) THEN
2600	enddo ! do izone=1,ntraceurs_zone
2601	! WRITE(,) 'apres bidouille: xt=',x(1:ntraciso)
2602	endif !if (dqtmp.gt.0.0) THEN
2603	#ifdef ISOVERIF
2604	CALL iso_verif_traceur(x,'iso_verif_traceurs 612')
2605	#endif
2606	endif !if (local_q(i,k).lt.0.0) THEN
2607	#ifdef ISOVERIF
2608	CALL iso_verif_traceur_pbidouille(x,'iso_verif_trac 625')
2609	#endif
2610	endif ! if (bidouille_anti_divergence) THEN
2611	END SUBROUTINE iso_verif_traceur_jbid_pos2
2612
2613	SUBROUTINE iso_verif_traceur_jbid_vect1D(x,n)
2614	USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
2615	IMPLICIT NONE
2616
2617	! version vectrisée de iso_verif_traceur_jbidouille
2618
2619	! inputs
2620	INTEGER n
2621	REAL x(ntraciso,n)
2622
2623	! locals
2624	INTEGER iiso,izone,ixt,i
2625	REAL xtractot
2626
2627	IF (bidouille_anti_divergence) THEN
2628	DO i=1,n
2629	DO iiso=1,niso
2630
2631	xtractot=0.0
2632	DO izone=1,ntraceurs_zone
2633	ixt=index_trac(izone,iiso)
2634	xtractot=xtractot+x(ixt,i)
2635	enddo !do izone=1,ntraceurs_zone
2636
2637	! on réajuste pour que les traceurs fasses bien la somme
2638	! des traceurs
2639	IF (abs(xtractot).gt.ridicule*10) THEN
2640	DO izone=1,ntraceurs_zone
2641	ixt=index_trac(izone,iiso)
2642	x(ixt,i)=x(ixt,i)/xtractot*x(iiso,i)
2643	enddo !do izone=1,ntraceurs_zone
2644	endif
2645
2646	enddo !do iiso=1,ntraceurs_iso
2647	enddo !do i=1,n
2648	endif !if (bidouille_anti_divergence) THEN
2649	END SUBROUTINE iso_verif_traceur_jbid_vect1D
2650
2651	! on met ces routines ici pour éviter dépendances circulaires
2652	#ifdef ISOVERIF
2653
2654	SUBROUTINE iso_verif_traceur_pbidouille(x,err_msg)
2655	USE isotopes_verif_mod
2656	IMPLICIT NONE
2657	! vérifier des choses sur les traceurs
2658	! * toutes les zones donne t l'istope total
2659	! * pas de deltaD aberrant
2660	! on réajuste aussi les valeurs des traceurs pour la
2661	! conservation de la masse, dans le cas bidouille
2662
2663	! on prend les valeurs pas défaut pour
2664	! errmax,errmaxrel,ridicule_trac,deltalimtrac
2665
2666	! inputs
2667	REAL x(ntraciso)
2668	CHARACTER() err_msg ! message d''erreur à afficher
2669
2670	! local
2671	!integer iso_verif_traceur_pbid_ns
2672
2673	IF (iso_verif_traceur_pbid_ns(x,err_msg).EQ.1) THEN
2674	stop
2675	endif
2676
2677	END SUBROUTINE iso_verif_traceur_pbidouille
2678
2679	function iso_verif_traceur_pbid_ns(x,err_msg)
2680	USE isotopes_mod, ONLY: iso_HDO,bidouille_anti_divergence
2681	USE isotrac_mod, ONLY: ridicule_trac
2682	USE isotopes_verif_mod
2683	IMPLICIT NONE
2684	! vérifier des choses sur les traceurs
2685	! * toutes les zones donne t l'istope total
2686	! * pas de deltaD aberrant
2687	! on réajuste aussi les valeurs des traceurs pour la
2688	! conservation de la masse, dans le cas bidouille
2689
2690	! on prend les valeurs pas défaut pour
2691	! errmax,errmaxrel,ridicule_trac,deltalimtrac
2692
2693	! inputs
2694	REAL x(ntraciso)
2695	CHARACTER() err_msg ! message d''erreur à afficher
2696
2697	! output
2698	INTEGER iso_verif_traceur_pbid_ns
2699
2700	! locals
2701	!integer iso_verif_traceur_noNaN_nostop
2702	!integer iso_verif_tracm_choix_nostop
2703	!integer iso_verif_tracdD_choix_nostop
2704	INTEGER iiso,izone,ixt
2705	REAL xtractot
2706
2707	! verif noNaN
2708	iso_verif_traceur_pbid_ns=0
2709
2710	IF (iso_verif_traceur_noNaN_nostop(x,err_msg).EQ.1) THEN
2711	! stop
2712	iso_verif_traceur_pbid_ns=1
2713	endif
2714
2715	! verif masse
2716	IF (iso_verif_tracm_choix_nostop(x,err_msg, &
2717	errmax10,errmaxrel50).EQ.1) THEN
2718	! on est plus laxiste car ça vient en général après une
2719	! bidouille pour iso_eau normal
2720	! stop
2721	iso_verif_traceur_pbid_ns=1
2722	endif
2723
2724	IF (bidouille_anti_divergence) THEN
2725	! on réajuste pour que les traceurs fasses bien la somme
2726	! des traceurs
2727	CALL iso_verif_traceur_jbidouille(x)
2728	endif !if (bidouille_anti_divergence) THEN
2729	! verif deltaD
2730	IF (iso_HDO.gt.0) THEN
2731	IF (iso_verif_tracdD_choix_nostop(x,err_msg, &
2732	ridicule_trac,deltalimtrac).EQ.1) THEN
2733	! stop
2734	iso_verif_traceur_pbid_ns=1
2735	endif
2736	endif !if (iso_HDO.gt.0) THEN
2737	END FUNCTION iso_verif_traceur_pbid_ns
2738
2739	SUBROUTINE iso_verif_traceur_pbid_vect(x,n,m,err_msg)
2740	USE isotopes_mod, ONLY: iso_HDO,bidouille_anti_divergence
2741	USE isotopes_verif_mod
2742	IMPLICIT NONE
2743
2744	! inputs
2745	INTEGER n,m
2746	REAL x(ntraciso,n,m)
2747	CHARACTER() err_msg ! message d''erreur à afficher
2748
2749	! locals
2750	INTEGER iiso,izone,ixt
2751	REAL xtractot
2752
2753	! verif noNaN
2754	CALL iso_verif_traceur_noNaN_vect(x,n,m,err_msg)
2755
2756	! verif masse
2757	CALL iso_verif_trac_masse_vect(x,n,m,err_msg,errmax*10, &
2758	errmaxrel*50)
2759
2760	IF (bidouille_anti_divergence) THEN
2761	! on réajuste pour que les traceurs fasses bien la somme
2762	! des traceurs
2763	CALL iso_verif_traceur_jbid_vect(x,n,m)
2764	endif !if (bidouille_anti_divergence) THEN
2765	! verif deltaD
2766	IF (iso_HDO.gt.0) THEN
2767	CALL iso_verif_tracdd_vect(x,n,m,err_msg)
2768	endif
2769
2770	END SUBROUTINE iso_verif_traceur_pbid_vect
2771	#endif
2772
2773	END MODULE isotrac_routines_mod
2774	#endif
2775	#endif

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

Download in other formats:

Original Format