Context Navigation

infotrac.F90 @ 4005

Last change on this file since 4005 was 4005, checked in by Ehouarn Millour, 3 years ago
Some cleanup around tracer name string length, and anticipate future needs and make default string max length larger (128 characters instead of 20). EM
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:keywords set to `Id`
File size: 36.8 KB

Line
1	! $Id: infotrac.F90 4005 2021-11-10 09:40:20Z emillour $
2	!
3	MODULE infotrac
4
5	! nqtot : total number of tracers and higher order of moment, water vapor and liquid included
6	INTEGER, SAVE :: nqtot
7	!CR: on ajoute le nombre de traceurs de l eau
8	INTEGER, SAVE :: nqo
9
10	! nbtr : number of tracers not including higher order of moment or water vapor or liquid
11	! number of tracers used in the physics
12	INTEGER, SAVE :: nbtr
13
14	! CRisi: on retranche les isotopes des traceurs habituels
15	! On fait un tableaux d'indices des traceurs qui passeront dans phytrac
16	INTEGER, SAVE :: nqtottr
17	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: itr_indice
18
19	! CRisi: nb traceurs peres= directement advectes par l'air
20	INTEGER, SAVE :: nqperes
21
22	! ThL: nb traceurs INCA
23	INTEGER, SAVE :: nqINCA
24
25	! ThL: nb traceurs CO2
26	INTEGER, SAVE :: nqCO2
27
28	! Name variables
29	INTEGER,PARAMETER :: tname_lenmax=128
30	CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:), SAVE :: tname ! tracer short name for restart and diagnostics
31	CHARACTER(len=tname_lenmax+3), ALLOCATABLE, DIMENSION(:), SAVE :: ttext ! tracer long name for diagnostics
32
33	! iadv : index of trasport schema for each tracer
34	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: iadv
35
36	! niadv : vector keeping the coorspondance between all tracers(nqtot) treated in the
37	! dynamic part of the code and the tracers (nbtr+2) used in the physics part of the code.
38	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: niadv ! equivalent dyn / physique
39
40	! CRisi: tableaux de fils
41	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: nqfils
42	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: nqdesc ! nombres de fils + nombre de tous les petits fils sur toutes les generations
43	INTEGER, SAVE :: nqdesc_tot
44	INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE :: iqfils
45	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: iqpere
46	REAL :: qperemin,masseqmin,ratiomin ! MVals et CRisi
47	PARAMETER (qperemin=1e-16,masseqmin=1e-16,ratiomin=1e-16) ! MVals
48
49	! conv_flg(it)=0 : convection desactivated for tracer number it
50	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: conv_flg
51	! pbl_flg(it)=0 : boundary layer diffusion desactivaded for tracer number it
52	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: pbl_flg
53
54	CHARACTER(len=4),SAVE :: type_trac
55	CHARACTER(len=8),DIMENSION(:),ALLOCATABLE, SAVE :: solsym
56
57	! CRisi: cas particulier des isotopes
58	LOGICAL,SAVE :: ok_isotopes,ok_iso_verif,ok_isotrac,ok_init_iso
59	INTEGER :: niso_possibles
60	PARAMETER ( niso_possibles=5)
61	REAL, DIMENSION (niso_possibles),SAVE :: tnat,alpha_ideal
62	LOGICAL, DIMENSION(niso_possibles),SAVE :: use_iso
63	INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE :: iqiso ! donne indice iq en fn de (ixt,phase)
64	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: iso_num ! donne numero iso entre 1 et niso_possibles en fn de nqtot
65	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: iso_indnum ! donne numero iso entre 1 et niso effectif en fn de nqtot
66	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: zone_num ! donne numero de la zone de tracage en fn de nqtot
67	INTEGER, ALLOCATABLE, DIMENSION(:), SAVE :: phase_num ! donne numero de la zone de tracage en fn de nqtot
68	INTEGER, DIMENSION(niso_possibles), SAVE :: indnum_fn_num ! donne indice entre entre 1 et niso en fonction du numero d isotope entre 1 et niso_possibles
69	INTEGER, ALLOCATABLE, DIMENSION(:,:), SAVE :: index_trac ! numero ixt en fn izone, indnum entre 1 et niso
70	INTEGER,SAVE :: niso,ntraceurs_zone,ntraciso
71
72	#ifdef CPP_StratAer
73	!--CK/OB for stratospheric aerosols
74	INTEGER, SAVE :: nbtr_bin
75	INTEGER, SAVE :: nbtr_sulgas
76	INTEGER, SAVE :: id_OCS_strat
77	INTEGER, SAVE :: id_SO2_strat
78	INTEGER, SAVE :: id_H2SO4_strat
79	INTEGER, SAVE :: id_BIN01_strat
80	INTEGER, SAVE :: id_TEST_strat
81	#endif
82
83	CONTAINS
84
85	SUBROUTINE infotrac_init
86	USE control_mod, ONLY: planet_type, config_inca
87	#ifdef REPROBUS
88	USE CHEM_REP, ONLY : Init_chem_rep_trac
89	#endif
90	IMPLICIT NONE
91	!=======================================================================
92	!
93	! Auteur: P. Le Van /L. Fairhead/F.Hourdin
94	! -------
95	! Modif special traceur F.Forget 05/94
96	! Modif M-A Filiberti 02/02 lecture de traceur.def
97	!
98	! Objet:
99	! ------
100	! GCM LMD nouvelle grille
101	!
102	!=======================================================================
103	! ... modification de l'integration de q ( 26/04/94 ) ....
104	!-----------------------------------------------------------------------
105	! Declarations
106
107	INCLUDE "dimensions.h"
108	INCLUDE "iniprint.h"
109
110	! Local variables
111	INTEGER, ALLOCATABLE, DIMENSION(:) :: hadv ! index of horizontal trasport schema
112	INTEGER, ALLOCATABLE, DIMENSION(:) :: vadv ! index of vertical trasport schema
113
114	INTEGER, ALLOCATABLE, DIMENSION(:) :: hadv_inca ! index of horizontal trasport schema
115	INTEGER, ALLOCATABLE, DIMENSION(:) :: vadv_inca ! index of vertical trasport schema
116
117	INTEGER, ALLOCATABLE, DIMENSION(:) :: conv_flg_inca
118	INTEGER, ALLOCATABLE, DIMENSION(:) :: pbl_flg_inca
119	CHARACTER(len=8), ALLOCATABLE, DIMENSION(:) :: solsym_inca
120
121	CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:) :: tnom_0 ! tracer short name
122	CHARACTER(len=tname_lenmax), ALLOCATABLE, DIMENSION(:) :: tnom_transp ! transporting fluid short name: CRisi
123	CHARACTER(len=3), DIMENSION(30) :: descrq
124	CHARACTER(len=1), DIMENSION(3) :: txts
125	CHARACTER(len=2), DIMENSION(9) :: txtp
126	CHARACTER(len=tname_lenmax) :: str1,str2
127
128	INTEGER :: nqtrue ! number of tracers read from tracer.def, without higer order of moment
129	INTEGER :: iq, new_iq, iiq, jq, ierr,itr
130	INTEGER :: ifils,ipere,generation ! CRisi
131	LOGICAL :: continu,nouveau_traceurdef
132	INTEGER :: IOstatus ! gestion de la retrocompatibilite de traceur.def
133	CHARACTER(len=2*tname_lenmax+1) :: tchaine
134
135	character(len=*),parameter :: modname="infotrac_init"
136
137	!-----------------------------------------------------------------------
138	! Initialization :
139	!
140	txts=(/'x','y','z'/)
141	txtp=(/'x ','y ','z ','xx','xy','xz','yy','yz','zz'/)
142
143	descrq(14)='VLH'
144	descrq(10)='VL1'
145	descrq(11)='VLP'
146	descrq(12)='FH1'
147	descrq(13)='FH2'
148	descrq(16)='PPM'
149	descrq(17)='PPS'
150	descrq(18)='PPP'
151	descrq(20)='SLP'
152	descrq(30)='PRA'
153
154
155	! Coherence test between parameter type_trac, config_inca and preprocessing keys
156	IF (type_trac=='inca') THEN
157	WRITE(lunout,*) 'You have chosen to couple with INCA chemistry model : type_trac=', &
158	type_trac,' config_inca=',config_inca
159	IF (config_inca/='aero' .AND. config_inca/='aeNP' .AND. config_inca/='chem') THEN
160	WRITE(lunout,*) 'Incoherence between type_trac and config_inca. Model stops. Modify run.def'
161	CALL abort_gcm('infotrac_init','Incoherence between type_trac and config_inca',1)
162	ENDIF
163	#ifndef INCA
164	WRITE(lunout,*) 'To run this option you must add cpp key INCA and compile with INCA code'
165	CALL abort_gcm('infotrac_init','You must compile with cpp key INCA',1)
166	#endif
167	ELSE IF (type_trac=='repr') THEN
168	WRITE(lunout,*) 'You have chosen to couple with REPROBUS chemestry model : type_trac=', type_trac
169	#ifndef REPROBUS
170	WRITE(lunout,*) 'To run this option you must add cpp key REPROBUS and compile with REPRPBUS code'
171	CALL abort_gcm('infotrac_init','You must compile with cpp key REPROBUS',1)
172	#endif
173	ELSE IF (type_trac == 'co2i') THEN
174	WRITE(lunout,*) 'You have chosen to run with CO2 cycle: type_trac=', type_trac
175	ELSE IF (type_trac == 'coag') THEN
176	WRITE(lunout,*) 'Tracers are treated for COAGULATION tests : type_trac=', type_trac
177	#ifndef CPP_StratAer
178	WRITE(lunout,*) 'To run this option you must add cpp key StratAer and compile with StratAer code'
179	CALL abort_gcm('infotrac_init','You must compile with cpp key StratAer',1)
180	#endif
181	ELSE IF (type_trac == 'lmdz') THEN
182	WRITE(lunout,*) 'Tracers are treated in LMDZ only : type_trac=', type_trac
183	ELSE IF (type_trac == 'inco') THEN ! ThL
184	WRITE(lunout,*) 'Using jointly INCA and CO2 cycle: type_trac =', type_trac
185	IF (config_inca/='aero' .AND. config_inca/='aeNP' .AND. config_inca/='chem') THEN
186	WRITE(lunout,*) 'Incoherence between type_trac and config_inca. Model stops. Modify run.def'
187	CALL abort_gcm('infotrac_init','Incoherence between type_trac and config_inca',1)
188	ENDIF
189	#ifndef INCA
190	WRITE(lunout,*) 'To run this option you must add cpp key INCA and compilewith INCA code'
191	CALL abort_gcm('infotrac_init','You must compile with cpp key INCA',1)
192	#endif
193	ELSE
194	WRITE(lunout,*) 'type_trac=',type_trac,' not possible. Model stops'
195	CALL abort_gcm('infotrac_init','bad parameter',1)
196	ENDIF
197
198	! Test if config_inca is other then none for run without INCA
199	IF (type_trac/='inca' .AND. type_trac/='inco' .AND. config_inca/='none') THEN
200	WRITE(lunout,*) 'config_inca will now be changed to none as you do not couple with INCA model'
201	config_inca='none'
202	ENDIF
203
204	!-----------------------------------------------------------------------
205	!
206	! 1) Get the true number of tracers + water vapor/liquid
207	! Here true tracers (nqtrue) means declared tracers (only first order)
208	!
209	!-----------------------------------------------------------------------
210	IF (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac == 'coag' .OR. type_trac == 'co2i') THEN
211	IF (type_trac=='co2i') THEN
212	nqCO2 = 1
213	ELSE
214	nqCO2 = 0
215	ENDIF
216	OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
217	IF(ierr.EQ.0) THEN
218	WRITE(lunout,*) trim(modname),': Open traceur.def : ok'
219	READ(90,*) nqtrue
220	write(lunout,*) 'nqtrue=',nqtrue
221	ELSE
222	WRITE(lunout,*) trim(modname),': Failed opening traceur.def'
223	CALL abort_gcm(modname,"file traceur.def not found!",1)
224	ENDIF
225	!jyg<
226	!! if ( planet_type=='earth') then
227	!! ! For Earth, water vapour & liquid tracers are not in the physics
228	!! nbtr=nqtrue-2
229	!! else
230	!! ! Other planets (for now); we have the same number of tracers
231	!! ! in the dynamics than in the physics
232	!! nbtr=nqtrue
233	!! endif
234	!>jyg
235	ELSE ! type_trac=inca or inco
236	IF (type_trac=='inco') THEN
237	nqCO2 = 1
238	ELSE
239	nqCO2 = 0
240	ENDIF
241	!jyg<
242	! The traceur.def file is used to define the number "nqo" of water phases
243	! present in the simulation. Default : nqo = 2.
244	OPEN(90,file='traceur.def',form='formatted',status='old', iostat=ierr)
245	IF(ierr.EQ.0) THEN
246	WRITE(lunout,*) trim(modname),': Open traceur.def : ok'
247	READ(90,*) nqo
248	ELSE
249	WRITE(lunout,*) trim(modname),': Failed opening traceur.def'
250	CALL abort_gcm(modname,"file traceur.def not found!",1)
251	ENDIF
252	IF (nqo /= 2 .AND. nqo /= 3 ) THEN
253	IF (nqo == 4 .AND. type_trac=='inco') THEN ! ThL
254	WRITE(lunout,*) trim(modname),': you are coupling with INCA, and also using CO2i.'
255	nqo = 3 ! A ameliorier... je force 3 traceurs eau... ThL
256	WRITE(lunout,*) trim(modname),': nqo = ',nqo
257	ELSE
258	WRITE(lunout,*) trim(modname),': nqo=',nqo, ' is not allowed. Only 2 or 3 water phases allowed'
259	CALL abort_gcm('infotrac_init','Bad number of water phases',1)
260	ENDIF
261	ENDIF
262	! nbtr has been read from INCA by init_const_lmdz() in gcm.F
263	#ifdef INCA
264	CALL Init_chem_inca_trac(nqINCA)
265	#else
266	nqINCA=0
267	#endif
268	nbtr=nqINCA+nqCO2
269	nqtrue=nbtr+nqo
270	WRITE(lunout,*) trim(modname),': nqo = ',nqo
271	WRITE(lunout,*) trim(modname),': nbtr = ',nbtr
272	WRITE(lunout,*) trim(modname),': nqtrue = ',nqtrue
273	WRITE(lunout,*) trim(modname),': nqCO2 = ',nqCO2
274	WRITE(lunout,*) trim(modname),': nqINCA = ',nqINCA
275	ALLOCATE(hadv_inca(nqINCA), vadv_inca(nqINCA), conv_flg_inca(nqINCA), pbl_flg_inca(nqINCA), solsym_inca(nqINCA))
276	ENDIF ! type_trac 'inca' ou 'inco'
277	!>jyg
278
279	IF ((planet_type=="earth").and.(nqtrue < 2)) THEN
280	WRITE(lunout,*) trim(modname),': nqtrue=',nqtrue, ' is not allowed. 2 tracers is the minimum'
281	CALL abort_gcm('infotrac_init','Not enough tracers',1)
282	ENDIF
283
284	!jyg<
285
286	!
287	! Allocate variables depending on nqtrue
288	!
289	ALLOCATE(tnom_0(nqtrue), hadv(nqtrue), vadv(nqtrue),tnom_transp(nqtrue))
290
291
292	!-----------------------------------------------------------------------
293	! 2) Choix des schemas d'advection pour l'eau et les traceurs
294	!
295	! iadv = 1 schema transport type "humidite specifique LMD"
296	! iadv = 2 schema amont
297	! iadv = 14 schema Van-leer + humidite specifique
298	! Modif F.Codron
299	! iadv = 10 schema Van-leer (retenu pour l'eau vapeur et liquide)
300	! iadv = 11 schema Van-Leer pour hadv et version PPM (Monotone) pour vadv
301	! iadv = 12 schema Frederic Hourdin I
302	! iadv = 13 schema Frederic Hourdin II
303	! iadv = 16 schema PPM Monotone(Collela & Woodward 1984)
304	! iadv = 17 schema PPM Semi Monotone (overshoots autorises)
305	! iadv = 18 schema PPM Positif Defini (overshoots undershoots autorises)
306	! iadv = 20 schema Slopes
307	! iadv = 30 schema Prather
308	!
309	! Dans le tableau q(ij,l,iq) : iq = 1 pour l'eau vapeur
310	! iq = 2 pour l'eau liquide
311	! Et eventuellement iq = 3,nqtot pour les autres traceurs
312	!
313	! iadv(1): choix pour l'eau vap. et iadv(2) : choix pour l'eau liq.
314	!------------------------------------------------------------------------
315	!
316	! Get choice of advection schema from file tracer.def or from INCA
317	!---------------------------------------------------------------------
318	IF (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac == 'coag' .OR. type_trac == 'co2i') THEN
319
320	! Continue to read tracer.def
321	DO iq=1,nqtrue
322
323	write(,) 'infotrac 237: iq=',iq
324	! CRisi: ajout du nom du fluide transporteur
325	! mais rester retro compatible
326	READ(90,'(I2,X,I2,X,A)',IOSTAT=IOstatus) hadv(iq),vadv(iq),tchaine
327	write(lunout,*) 'iq,hadv(iq),vadv(iq)=',iq,hadv(iq),vadv(iq)
328	write(lunout,*) 'tchaine=',trim(tchaine)
329	write(,) 'infotrac 238: IOstatus=',IOstatus
330	if (IOstatus.ne.0) then
331	CALL abort_gcm('infotrac_init','Pb dans la lecture de traceur.def',1)
332	endif
333	! Y-a-t-il 1 ou 2 noms de traceurs? -> On regarde s'il y a un
334	! espace ou pas au milieu de la chaine.
335	continu=.true.
336	nouveau_traceurdef=.false.
337	iiq=1
338	do while (continu)
339	if (tchaine(iiq:iiq).eq.' ') then
340	nouveau_traceurdef=.true.
341	continu=.false.
342	else if (iiq.lt.LEN_TRIM(tchaine)) then
343	iiq=iiq+1
344	else
345	continu=.false.
346	endif
347	enddo
348	write(,) 'iiq,nouveau_traceurdef=',iiq,nouveau_traceurdef
349	if (nouveau_traceurdef) then
350	write(lunout,*) 'C''est la nouvelle version de traceur.def'
351	tnom_0(iq)=tchaine(1:iiq-1)
352	tnom_transp(iq)=tchaine(iiq+1:)
353	else
354	write(lunout,*) 'C''est l''ancienne version de traceur.def'
355	write(lunout,*) 'On suppose que les traceurs sont tous d''air'
356	tnom_0(iq)=tchaine
357	tnom_transp(iq) = 'air'
358	endif
359	write(lunout,*) 'tnom_0(iq)=<',trim(tnom_0(iq)),'>'
360	write(lunout,*) 'tnom_transp(iq)=<',trim(tnom_transp(iq)),'>'
361
362	ENDDO!DO iq=1,nqtrue
363
364	CLOSE(90)
365
366	WRITE(lunout,*) trim(modname),': Valeur de traceur.def :'
367	WRITE(lunout,*) trim(modname),': nombre total de traceurs ',nqtrue
368	DO iq=1,nqtrue
369	WRITE(lunout,*) hadv(iq),vadv(iq),' ',trim(tnom_0(iq)),' ',trim(tnom_transp(iq))
370	END DO
371
372	IF ( planet_type=='earth') THEN
373	!CR: nombre de traceurs de l eau
374	IF (tnom_0(3) == 'H2Oi') THEN
375	nqo=3
376	ELSE
377	nqo=2
378	ENDIF
379	! For Earth, water vapour & liquid tracers are not in the physics
380	nbtr=nqtrue-nqo
381	ELSE
382	! Other planets (for now); we have the same number of tracers
383	! in the dynamics than in the physics
384	nbtr=nqtrue
385	ENDIF
386
387	#ifdef CPP_StratAer
388	IF (type_trac == 'coag') THEN
389	nbtr_bin=0
390	nbtr_sulgas=0
391	DO iq=1,nqtrue
392	IF (tnom_0(iq)(1:3)=='BIN') THEN !check if tracer name contains 'BIN'
393	nbtr_bin=nbtr_bin+1
394	ENDIF
395	IF (tnom_0(iq)(1:3)=='GAS') THEN !check if tracer name contains 'GAS'
396	nbtr_sulgas=nbtr_sulgas+1
397	ENDIF
398	ENDDO
399	print*,'nbtr_bin=',nbtr_bin
400	print*,'nbtr_sulgas=',nbtr_sulgas
401	DO iq=1,nqtrue
402	IF (tnom_0(iq)=='GASOCS') THEN
403	id_OCS_strat=iq-nqo
404	ENDIF
405	IF (tnom_0(iq)=='GASSO2') THEN
406	id_SO2_strat=iq-nqo
407	ENDIF
408	IF (tnom_0(iq)=='GASH2SO4') THEN
409	id_H2SO4_strat=iq-nqo
410	ENDIF
411	IF (tnom_0(iq)=='BIN01') THEN
412	id_BIN01_strat=iq-nqo
413	ENDIF
414	IF (tnom_0(iq)=='GASTEST') THEN
415	id_TEST_strat=iq-nqo
416	ENDIF
417	ENDDO
418	print*,'id_OCS_strat =',id_OCS_strat
419	print*,'id_SO2_strat =',id_SO2_strat
420	print*,'id_H2SO4_strat=',id_H2SO4_strat
421	print*,'id_BIN01_strat=',id_BIN01_strat
422	ENDIF
423	#endif
424
425	ENDIF ! (type_trac == 'lmdz' .OR. type_trac == 'repr' .OR. type_trac = 'coag' .OR. type_trac = 'co2i')
426	!jyg<
427	!
428
429	! Transfert number of tracers to Reprobus
430	IF (type_trac == 'repr') THEN
431	#ifdef REPROBUS
432	CALL Init_chem_rep_trac(nbtr,nqo,tnom_0)
433	#endif
434	ENDIF
435	!
436	! Allocate variables depending on nbtr
437	!
438	ALLOCATE(conv_flg(nbtr), pbl_flg(nbtr), solsym(nbtr))
439	conv_flg(:) = 1 ! convection activated for all tracers
440	pbl_flg(:) = 1 ! boundary layer activated for all tracers
441
442	IF (type_trac == 'inca' .OR. type_trac == 'inco') THEN ! config_inca='aero' ou 'chem'
443	!>jyg
444	! le module de chimie fournit les noms des traceurs
445	! et les schemas d'advection associes. excepte pour ceux lus
446	! dans traceur.def
447
448	DO iq=1,nqo+nqCO2
449
450	write(,) 'infotrac 237: iq=',iq
451	! CRisi: ajout du nom du fluide transporteur
452	! mais rester retro compatible
453	READ(90,'(I2,X,I2,X,A)',IOSTAT=IOstatus) hadv(iq),vadv(iq),tchaine
454	write(lunout,*) 'iq,hadv(iq),vadv(iq)=',iq,hadv(iq),vadv(iq)
455	write(lunout,*) 'tchaine=',trim(tchaine)
456	write(,) 'infotrac 238: IOstatus=',IOstatus
457	if (IOstatus.ne.0) then
458	CALL abort_gcm('infotrac_init','Pb dans la lecture de traceur.def',1)
459	endif
460	! Y-a-t-il 1 ou 2 noms de traceurs? -> On regarde s'il y a un
461	! espace ou pas au milieu de la chaine.
462	continu=.true.
463	nouveau_traceurdef=.false.
464	iiq=1
465
466	do while (continu)
467	if (tchaine(iiq:iiq).eq.' ') then
468	nouveau_traceurdef=.true.
469	continu=.false.
470	else if (iiq.lt.LEN_TRIM(tchaine)) then
471	iiq=iiq+1
472	else
473	continu=.false.
474	endif
475	enddo
476
477	write(,) 'iiq,nouveau_traceurdef=',iiq,nouveau_traceurdef
478
479	if (nouveau_traceurdef) then
480	write(lunout,*) 'C''est la nouvelle version de traceur.def'
481	tnom_0(iq)=tchaine(1:iiq-1)
482	tnom_transp(iq)=tchaine(iiq+1:)
483	else
484	write(lunout,*) 'C''est l''ancienne version de traceur.def'
485	write(lunout,*) 'On suppose que les traceurs sont tous d''air'
486	tnom_0(iq)=tchaine
487	tnom_transp(iq) = 'air'
488	endif
489
490	write(lunout,*) 'tnom_0(iq)=<',trim(tnom_0(iq)),'>'
491	write(lunout,*) 'tnom_transp(iq)=<',trim(tnom_transp(iq)),'>'
492
493	ENDDO !DO iq=1,nqo
494	CLOSE(90)
495
496
497	#ifdef INCA
498	CALL init_transport( &
499	hadv_inca, &
500	vadv_inca, &
501	conv_flg_inca, &
502	pbl_flg_inca, &
503	solsym_inca)
504
505	conv_flg(1+nqCO2:nbtr) = conv_flg_inca
506	pbl_flg(1+nqCO2:nbtr) = pbl_flg_inca
507	solsym(1+nqCO2:nbtr) = solsym_inca
508
509	IF (type_trac == 'inco') THEN
510	conv_flg(1:nqCO2) = 1
511	pbl_flg(1:nqCO2) = 1
512	solsym(1:nqCO2) = 'CO2'
513	ENDIF
514	#endif
515
516	!jyg<
517	DO iq = nqo+nqCO2+1, nqtrue
518	hadv(iq) = hadv_inca(iq-nqo-nqCO2)
519	vadv(iq) = vadv_inca(iq-nqo-nqCO2)
520	tnom_0(iq)=solsym_inca(iq-nqo-nqCO2)
521	tnom_transp(iq) = 'air'
522	END DO
523
524	ENDIF ! (type_trac == 'inca' or 'inco')
525
526	!-----------------------------------------------------------------------
527	!
528	! 3) Verify if advection schema 20 or 30 choosen
529	! Calculate total number of tracers needed: nqtot
530	! Allocate variables depending on total number of tracers
531	!-----------------------------------------------------------------------
532	new_iq=0
533	DO iq=1,nqtrue
534	! Add tracers for certain advection schema
535	IF (hadv(iq)<20 .AND. vadv(iq)<20 ) THEN
536	new_iq=new_iq+1 ! no tracers added
537	ELSE IF (hadv(iq)==20 .AND. vadv(iq)==20 ) THEN
538	new_iq=new_iq+4 ! 3 tracers added
539	ELSE IF (hadv(iq)==30 .AND. vadv(iq)==30 ) THEN
540	new_iq=new_iq+10 ! 9 tracers added
541	ELSE
542	WRITE(lunout,*) trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
543	CALL abort_gcm('infotrac_init','Bad choice of advection schema - 1',1)
544	ENDIF
545	END DO
546
547	IF (new_iq /= nqtrue) THEN
548	! The choice of advection schema imposes more tracers
549	! Assigne total number of tracers
550	nqtot = new_iq
551
552	WRITE(lunout,*) trim(modname),': The choice of advection schema for one or more tracers'
553	WRITE(lunout,*) 'makes it necessary to add tracers'
554	WRITE(lunout,*) trim(modname)//': ',nqtrue,' is the number of true tracers'
555	WRITE(lunout,*) trim(modname)//': ',nqtot, ' is the total number of tracers needed'
556
557	ELSE
558	! The true number of tracers is also the total number
559	nqtot = nqtrue
560	ENDIF
561
562	!
563	! Allocate variables with total number of tracers, nqtot
564	!
565	ALLOCATE(tname(nqtot), ttext(nqtot))
566	ALLOCATE(iadv(nqtot), niadv(nqtot))
567
568	!-----------------------------------------------------------------------
569	!
570	! 4) Determine iadv, long and short name
571	!
572	!-----------------------------------------------------------------------
573	new_iq=0
574	DO iq=1,nqtrue
575	new_iq=new_iq+1
576
577	! Verify choice of advection schema
578	IF (hadv(iq)==vadv(iq)) THEN
579	iadv(new_iq)=hadv(iq)
580	ELSE IF (hadv(iq)==10 .AND. vadv(iq)==16) THEN
581	iadv(new_iq)=11
582	ELSE
583	WRITE(lunout,*)trim(modname),': This choice of advection schema is not available',iq,hadv(iq),vadv(iq)
584
585	CALL abort_gcm('infotrac_init','Bad choice of advection schema - 2',1)
586	ENDIF
587
588	str1=tnom_0(iq)
589	tname(new_iq)= tnom_0(iq)
590	IF (iadv(new_iq)==0) THEN
591	ttext(new_iq)=trim(str1)
592	ELSE
593	ttext(new_iq)=trim(tnom_0(iq))//descrq(iadv(new_iq))
594	ENDIF
595
596	! schemas tenant compte des moments d'ordre superieur
597	str2=ttext(new_iq)
598	IF (iadv(new_iq)==20) THEN
599	DO jq=1,3
600	new_iq=new_iq+1
601	iadv(new_iq)=-20
602	ttext(new_iq)=trim(str2)//txts(jq)
603	tname(new_iq)=trim(str1)//txts(jq)
604	END DO
605	ELSE IF (iadv(new_iq)==30) THEN
606	DO jq=1,9
607	new_iq=new_iq+1
608	iadv(new_iq)=-30
609	ttext(new_iq)=trim(str2)//txtp(jq)
610	tname(new_iq)=trim(str1)//txtp(jq)
611	END DO
612	ENDIF
613	END DO
614
615	!
616	! Find vector keeping the correspodence between true and total tracers
617	!
618	niadv(:)=0
619	iiq=0
620	DO iq=1,nqtot
621	IF(iadv(iq).GE.0) THEN
622	! True tracer
623	iiq=iiq+1
624	niadv(iiq)=iq
625	ENDIF
626	END DO
627
628
629	WRITE(lunout,*) trim(modname),': Information stored in infotrac :'
630	WRITE(lunout,*) trim(modname),': iadv niadv tname ttext :'
631
632	DO iq=1,nqtot
633	WRITE(lunout,*) iadv(iq),niadv(iq), ' ',trim(tname(iq)),' ',trim(ttext(iq))
634	END DO
635
636	!
637	! Test for advection schema.
638	! This version of LMDZ only garantees iadv=10 and iadv=14 (14 only for water vapour) .
639	!
640	DO iq=1,nqtot
641	IF (iadv(iq)/=10 .AND. iadv(iq)/=14 .AND. iadv(iq)/=0) THEN
642	WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
643	CALL abort_gcm('infotrac_init','In this version only iadv=10 and iadv=14 is tested!',1)
644	ELSE IF (iadv(iq)==14 .AND. iq/=1) THEN
645	WRITE(lunout,*)trim(modname),'STOP : The option iadv=',iadv(iq),' is not tested in this version of LMDZ'
646	CALL abort_gcm('infotrac_init','In this version iadv=14 is only permitted for water vapour!',1)
647	ENDIF
648	END DO
649
650
651	! CRisi: quels sont les traceurs fils et les traceurs peres.
652	! initialiser tous les tableaux d'indices lies aux traceurs familiaux
653	! + verifier que tous les peres sont ecrits en premieres positions
654	ALLOCATE(nqfils(nqtot),nqdesc(nqtot))
655	ALLOCATE(iqfils(nqtot,nqtot))
656	ALLOCATE(iqpere(nqtot))
657	nqperes=0
658	nqfils(:)=0
659	nqdesc(:)=0
660	iqfils(:,:)=0
661	iqpere(:)=0
662	nqdesc_tot=0
663	DO iq=1,nqtot
664	if (tnom_transp(iq) == 'air') then
665	! ceci est un traceur pere
666	WRITE(lunout,*) 'Le traceur',iq,', appele ',trim(tnom_0(iq)),', est un pere'
667	nqperes=nqperes+1
668	iqpere(iq)=0
669	else !if (tnom_transp(iq) == 'air') then
670	! ceci est un fils. Qui est son pere?
671	WRITE(lunout,*) 'Le traceur',iq,', appele ',trim(tnom_0(iq)),', est un fils'
672	continu=.true.
673	ipere=1
674	do while (continu)
675	if (tnom_transp(iq) == tnom_0(ipere)) then
676	! Son pere est ipere
677	WRITE(lunout,*) 'Le traceur',iq,'appele ', &
678	& trim(tnom_0(iq)),' est le fils de ',ipere,'appele ',trim(tnom_0(ipere))
679	if (iq.eq.ipere) then
680	CALL abort_gcm('infotrac_init','Un fils est son propre pere',1)
681	endif
682	nqfils(ipere)=nqfils(ipere)+1
683	iqfils(nqfils(ipere),ipere)=iq
684	iqpere(iq)=ipere
685	continu=.false.
686	else !if (tnom_transp(iq) == tnom_0(ipere)) then
687	ipere=ipere+1
688	if (ipere.gt.nqtot) then
689	WRITE(lunout,*) 'Le traceur',iq,'appele ', &
690	& trim(tnom_0(iq)),', est orphelin.'
691	CALL abort_gcm('infotrac_init','Un traceur est orphelin',1)
692	endif !if (ipere.gt.nqtot) then
693	endif !if (tnom_transp(iq) == tnom_0(ipere)) then
694	enddo !do while (continu)
695	endif !if (tnom_transp(iq) == 'air') then
696	enddo !DO iq=1,nqtot
697	WRITE(lunout,*) 'infotrac: nqperes=',nqperes
698	WRITE(lunout,*) 'nqfils=',nqfils
699	WRITE(lunout,*) 'iqpere=',iqpere
700	WRITE(lunout,*) 'iqfils=',iqfils
701
702	! Calculer le nombre de descendants a partir de iqfils et de nbfils
703	DO iq=1,nqtot
704	generation=0
705	continu=.true.
706	ifils=iq
707	do while (continu)
708	ipere=iqpere(ifils)
709	if (ipere.gt.0) then
710	nqdesc(ipere)=nqdesc(ipere)+1
711	nqdesc_tot=nqdesc_tot+1
712	iqfils(nqdesc(ipere),ipere)=iq
713	ifils=ipere
714	generation=generation+1
715	else !if (ipere.gt.0) then
716	continu=.false.
717	endif !if (ipere.gt.0) then
718	enddo !do while (continu)
719	WRITE(lunout,*) 'Le traceur ',iq,', appele ',trim(tnom_0(iq)),' est un traceur de generation: ',generation
720	enddo !DO iq=1,nqtot
721	WRITE(lunout,*) 'infotrac: nqdesc=',nqdesc
722	WRITE(lunout,*) 'iqfils=',iqfils
723	WRITE(lunout,*) 'nqdesc_tot=',nqdesc_tot
724
725	! Interdire autres schemas que 10 pour les traceurs fils, et autres schemas
726	! que 10 et 14 si des peres ont des fils
727	do iq=1,nqtot
728	if (iqpere(iq).gt.0) then
729	! ce traceur a un pere qui n'est pas l'air
730	! Seul le schema 10 est autorise
731	if (iadv(iq)/=10) then
732	WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not implemented for sons'
733	CALL abort_gcm('infotrac_init','Sons should be advected by scheme 10',1)
734	endif
735	! Le traceur pere ne peut etre advecte que par schema 10 ou 14:
736	IF (iadv(iqpere(iq))/=10 .AND. iadv(iqpere(iq))/=14) THEN
737	WRITE(lunout,*)trim(modname),' STOP : The option iadv=',iadv(iq),' is not implemented for fathers'
738	CALL abort_gcm('infotrac_init','Fathers should be advected by scheme 10 ou 14',1)
739	endif !IF (iadv(iqpere(iq))/=10 .AND. iadv(iqpere(iq))/=14) THEN
740	endif !if (iqpere(iq).gt.0) the
741	enddo !do iq=1,nqtot
742
743
744
745	! detecter quels sont les traceurs isotopiques parmi des traceurs
746	call infotrac_isoinit(tnom_0,nqtrue)
747
748	! if (ntraciso.gt.0) then
749	! le 18 sep 2020: on enleve la condition ntraciso.gt.0 car nqtottr doit etre
750	! connu meme si il n'y a pas d'isotopes!
751	write(lunout,*) 'infotrac 702: nbtr,ntraciso=',nbtr,ntraciso
752	! retrancher les traceurs isotopiques de la liste des traceurs qui passent dans
753	! phytrac
754	nbtr=nbtr-nqo*ntraciso
755
756	! faire un tableau d'indice des traceurs qui passeront dans phytrac
757	nqtottr=nqtot-nqo*(1+ntraciso)
758	write(lunout,*) 'infotrac 704: nqtottr,nqtot,nqo=',nqtottr,nqtot,nqo
759	! Rq: nqtottr n'est pas forcement egal a nbtr dans le cas ou new_iq /= nqtrue
760	ALLOCATE (itr_indice(nqtottr))
761	itr_indice(:)=0
762	itr=0
763	do iq=nqo+1, nqtot
764	if (iso_num(iq).eq.0) then
765	itr=itr+1
766	write(,) 'itr=',itr
767	itr_indice(itr)=iq
768	endif !if (iso_num(iq).eq.0) then
769	enddo
770	if (itr.ne.nqtottr) then
771	CALL abort_gcm('infotrac_init','pb dans le calcul de nqtottr',1)
772	endif
773	write(lunout,*) 'itr_indice=',itr_indice
774	! endif !if (ntraciso.gt.0) then
775
776	!-----------------------------------------------------------------------
777	! Finalize :
778	!
779	DEALLOCATE(tnom_0, hadv, vadv,tnom_transp)
780
781	WRITE(lunout,*) 'infotrac init fin'
782
783	END SUBROUTINE infotrac_init
784
785	SUBROUTINE infotrac_isoinit(tnom_0,nqtrue)
786
787	#ifdef CPP_IOIPSL
788	use IOIPSL
789	#else
790	! if not using IOIPSL, we still need to use (a local version of) getin
791	use ioipsl_getincom
792	#endif
793	implicit none
794
795	! inputs
796	INTEGER,INTENT(IN) :: nqtrue
797	CHARACTER(len=*),INTENT(IN) :: tnom_0(nqtrue)
798
799	! locals
800	CHARACTER(len=3), DIMENSION(niso_possibles) :: tnom_iso
801	INTEGER, ALLOCATABLE,DIMENSION(:,:) :: nb_iso,nb_traciso
802	INTEGER, ALLOCATABLE,DIMENSION(:) :: nb_isoind
803	INTEGER :: ntraceurs_zone_prec,iq,phase,ixt,iiso,izone
804	CHARACTER(len=tname_lenmax) :: tnom_trac
805	INCLUDE "iniprint.h"
806
807	tnom_iso=(/'eau','HDO','O18','O17','HTO'/)
808
809	ALLOCATE(nb_iso(niso_possibles,nqo))
810	ALLOCATE(nb_isoind(nqo))
811	ALLOCATE(nb_traciso(niso_possibles,nqo))
812	ALLOCATE(iso_num(nqtot))
813	ALLOCATE(iso_indnum(nqtot))
814	ALLOCATE(zone_num(nqtot))
815	ALLOCATE(phase_num(nqtot))
816
817	iso_num(:)=0
818	iso_indnum(:)=0
819	zone_num(:)=0
820	phase_num(:)=0
821	indnum_fn_num(:)=0
822	use_iso(:)=.false.
823	nb_iso(:,:)=0
824	nb_isoind(:)=0
825	nb_traciso(:,:)=0
826	niso=0
827	ntraceurs_zone=0
828	ntraceurs_zone_prec=0
829	ntraciso=0
830
831	do iq=nqo+1,nqtot
832	! write(lunout,*) 'infotrac 569: iq,tnom_0(iq)=',iq,tnom_0(iq)
833	do phase=1,nqo
834	do ixt= 1,niso_possibles
835	tnom_trac=trim(tnom_0(phase))//'_'
836	tnom_trac=trim(tnom_trac)//trim(tnom_iso(ixt))
837	! write(,) 'phase,ixt,tnom_trac=',phase,ixt,tnom_trac
838	IF (tnom_0(iq) == tnom_trac) then
839	! write(lunout,*) 'Ce traceur est un isotope'
840	nb_iso(ixt,phase)=nb_iso(ixt,phase)+1
841	nb_isoind(phase)=nb_isoind(phase)+1
842	iso_num(iq)=ixt
843	iso_indnum(iq)=nb_isoind(phase)
844	indnum_fn_num(ixt)=iso_indnum(iq)
845	phase_num(iq)=phase
846	! write(lunout,*) 'iso_num(iq)=',iso_num(iq)
847	! write(lunout,*) 'iso_indnum(iq)=',iso_indnum(iq)
848	! write(lunout,*) 'indnum_fn_num(ixt)=',indnum_fn_num(ixt)
849	! write(lunout,*) 'phase_num(iq)=',phase_num(iq)
850	goto 20
851	else if (iqpere(iq).gt.0) then
852	if (tnom_0(iqpere(iq)) == tnom_trac) then
853	! write(lunout,*) 'Ce traceur est le fils d''un isotope'
854	! c'est un traceur d'isotope
855	nb_traciso(ixt,phase)=nb_traciso(ixt,phase)+1
856	iso_num(iq)=ixt
857	iso_indnum(iq)=indnum_fn_num(ixt)
858	zone_num(iq)=nb_traciso(ixt,phase)
859	phase_num(iq)=phase
860	! write(lunout,*) 'iso_num(iq)=',iso_num(iq)
861	! write(lunout,*) 'phase_num(iq)=',phase_num(iq)
862	! write(lunout,*) 'zone_num(iq)=',zone_num(iq)
863	goto 20
864	endif !if (tnom_0(iqpere(iq)) == trim(tnom_0(phase))//trim(tnom_iso(ixt))) then
865	endif !IF (tnom_0(iq) == trim(tnom_0(phase))//trim(tnom_iso(ixt))) then
866	enddo !do ixt= niso_possibles
867	enddo !do phase=1,nqo
868	20 continue
869	enddo !do iq=1,nqtot
870
871	! write(lunout,*) 'iso_num=',iso_num
872	! write(lunout,*) 'iso_indnum=',iso_indnum
873	! write(lunout,*) 'zone_num=',zone_num
874	! write(lunout,*) 'phase_num=',phase_num
875	! write(lunout,*) 'indnum_fn_num=',indnum_fn_num
876
877	do ixt= 1,niso_possibles
878
879	if (nb_iso(ixt,1).eq.1) then
880	! on verifie que toutes les phases ont le meme nombre de
881	! traceurs
882	do phase=2,nqo
883	if (nb_iso(ixt,phase).ne.nb_iso(ixt,1)) then
884	! write(lunout,*) 'ixt,phase,nb_iso=',ixt,phase,nb_iso(ixt,phase)
885	CALL abort_gcm('infotrac_init','Phases must have same number of isotopes',1)
886	endif
887	enddo !do phase=2,nqo
888
889	niso=niso+1
890	use_iso(ixt)=.true.
891	ntraceurs_zone=nb_traciso(ixt,1)
892
893	! on verifie que toutes les phases ont le meme nombre de
894	! traceurs
895	do phase=2,nqo
896	if (nb_traciso(ixt,phase).ne.ntraceurs_zone) then
897	write(lunout,*) 'ixt,phase,nb_traciso=',ixt,phase,nb_traciso(ixt,phase)
898	write(lunout,*) 'ntraceurs_zone=',ntraceurs_zone
899	CALL abort_gcm('infotrac_init','Phases must have same number of tracers',1)
900	endif
901	enddo !do phase=2,nqo
902	! on verifie que tous les isotopes ont le meme nombre de
903	! traceurs
904	if (ntraceurs_zone_prec.gt.0) then
905	if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
906	ntraceurs_zone_prec=ntraceurs_zone
907	else !if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
908	write(,) 'ntraceurs_zone_prec,ntraceurs_zone=',ntraceurs_zone_prec,ntraceurs_zone
909	CALL abort_gcm('infotrac_init', &
910	&'Isotope tracers are not well defined in traceur.def',1)
911	endif !if (ntraceurs_zone.eq.ntraceurs_zone_prec) then
912	endif !if (ntraceurs_zone_prec.gt.0) then
913
914	else if (nb_iso(ixt,1).ne.0) then
915	WRITE(lunout,*) 'nqo,ixt=',nqo,ixt
916	WRITE(lunout,*) 'nb_iso(ixt,1)=',nb_iso(ixt,1)
917	CALL abort_gcm('infotrac_init','Isotopes are not well defined in traceur.def',1)
918	endif !if (nb_iso(ixt,1).eq.1) then
919	enddo ! do ixt= niso_possibles
920
921	! dimensions isotopique:
922	ntraciso=niso*(ntraceurs_zone+1)
923	! WRITE(lunout,*) 'niso=',niso
924	! WRITE(lunout,*) 'ntraceurs_zone,ntraciso=',ntraceurs_zone,ntraciso
925
926	! flags isotopiques:
927	if (niso.gt.0) then
928	ok_isotopes=.true.
929	else
930	ok_isotopes=.false.
931	endif
932	! WRITE(lunout,*) 'ok_isotopes=',ok_isotopes
933
934	if (ok_isotopes) then
935	ok_iso_verif=.false.
936	call getin('ok_iso_verif',ok_iso_verif)
937	ok_init_iso=.false.
938	call getin('ok_init_iso',ok_init_iso)
939	tnat=(/1.0,155.76e-6,2005.2e-6,0.004/100.,0.0/)
940	alpha_ideal=(/1.0,1.01,1.006,1.003,1.0/)
941	endif !if (ok_isotopes) then
942	! WRITE(lunout,*) 'ok_iso_verif=',ok_iso_verif
943	! WRITE(lunout,*) 'ok_init_iso=',ok_init_iso
944
945	if (ntraceurs_zone.gt.0) then
946	ok_isotrac=.true.
947	else
948	ok_isotrac=.false.
949	endif
950	! WRITE(lunout,*) 'ok_isotrac=',ok_isotrac
951
952	! remplissage du tableau iqiso(ntraciso,phase)
953	ALLOCATE(iqiso(ntraciso,nqo))
954	iqiso(:,:)=0
955	do iq=1,nqtot
956	if (iso_num(iq).gt.0) then
957	ixt=iso_indnum(iq)+zone_num(iq)*niso
958	iqiso(ixt,phase_num(iq))=iq
959	endif
960	enddo
961	! WRITE(lunout,*) 'iqiso=',iqiso
962
963	! replissage du tableau index_trac(ntraceurs_zone,niso)
964	ALLOCATE(index_trac(ntraceurs_zone,niso))
965	if (ok_isotrac) then
966	do iiso=1,niso
967	do izone=1,ntraceurs_zone
968	index_trac(izone,iiso)=iiso+izone*niso
969	enddo
970	enddo
971	else !if (ok_isotrac) then
972	index_trac(:,:)=0.0
973	endif !if (ok_isotrac) then
974	! write(lunout,*) 'index_trac=',index_trac
975
976	! Finalize :
977	DEALLOCATE(nb_iso)
978
979	END SUBROUTINE infotrac_isoinit
980
981	END MODULE infotrac

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

Download in other formats:

Original Format