Context Navigation

infotrac.F90 @ 3870

Last change on this file since 3870 was 3870, checked in by Laurent Fairhead, 3 years ago
Corrections to r3865 that will hopefully repair the debug compilation TL
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: 37.0 KB

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

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

Download in other formats:

Original Format