Context Navigation

initracer.F @ 2616

Last change on this file since 2616 was 2616, checked in by romain.vande, 3 years ago
LMDZ_MARS RV : Open_MP; Put all the "save" variables as "!$OMP THREADPRIVATE" in phymars. The code can now be tested, see README for more info
File size: 33.5 KB

Line
1	SUBROUTINE initracer(ngrid,nq,qsurf)
2
3	use tracer_mod
4	use comcstfi_h, only: pi
5	use dust_param_mod, only: doubleq, submicron, dustbin
6	IMPLICIT NONE
7	c=======================================================================
8	c subject:
9	c --------
10	c Initialization related to tracer
11	c (transported dust, water, chemical species, ice...)
12	c
13	c Name of the tracer
14	c
15	c Test of dimension :
16	c Initialize tracer related data in tracer_mod, using tracer names provided
17	c by the dynamics in "infotrac"
18	c
19	c
20	c author: F.Forget
21	c ------
22	c Modifs: Franck Montmessin, Sebastien Lebonnois (june 2003)
23	c Ehouarn Millour (oct. 2008) identify tracers by their names
24	c=======================================================================
25
26
27	include "callkeys.h"
28
29	integer,intent(in) :: ngrid ! number of atmospheric columns
30	integer,intent(in) :: nq ! number of tracers
31	real,intent(out) :: qsurf(ngrid,nq) ! tracer on surface (e.g. kg.m-2)
32
33	integer iq,ig,count
34	real r0_lift , reff_lift, nueff_lift
35	real r0_storm,reff_storm
36	c Ratio of small over large dust particles (used when both
37	c doubleq and the submicron mode are active); In Montmessin
38	c et al. (2002), a value of 25 has been deduced;
39	real, parameter :: popratio = 25.
40	character(len=30) :: txt ! to store some text
41
42	c-----------------------------------------------------------------------
43	c radius(nq) ! aerosol particle radius (m)
44	c rho_q(nq) ! tracer densities (kg.m-3)
45	c alpha_lift(nq) ! saltation vertical flux/horiz flux ratio (m-1)
46	c alpha_devil(nq) ! lifting coeeficient by dust devil
47	c rho_dust ! Mars dust density
48	c rho_ice ! Water ice density
49	c nuice_ref ! Effective variance nueff of the
50	c ! water-ice size distributions
51	c doubleq ! if method with mass (iq=1) and number(iq=2) mixing ratio
52	c varian ! Characteristic variance of log-normal distribution
53	c-----------------------------------------------------------------------
54
55
56	c------------------------------------------------------------
57	c NAME and molar mass of the tracer
58	c------------------------------------------------------------
59
60	! Identify tracers by their names: (and set corresponding values of mmol)
61	! 0. initialize tracer indexes to zero:
62	igcm_dustbin(1:nq)=0
63	igcm_co2_ice=0
64	igcm_ccnco2_mass=0
65	igcm_ccnco2_number=0
66	igcm_ccnco2_meteor_mass=0
67	igcm_ccnco2_meteor_number=0
68	igcm_ccnco2_h2o_mass_ice=0
69	igcm_ccnco2_h2o_mass_ccn=0
70	igcm_ccnco2_h2o_number=0
71	igcm_dust_mass=0
72	igcm_dust_number=0
73	igcm_ccn_mass=0
74	igcm_ccn_number=0
75	igcm_dust_submicron=0
76	igcm_h2o_vap=0
77	igcm_h2o_ice=0
78	igcm_hdo_vap=0
79	igcm_hdo_ice=0
80	igcm_stormdust_mass=0
81	igcm_stormdust_number=0
82	igcm_topdust_mass=0
83	igcm_topdust_number=0
84	igcm_co2=0
85	igcm_co=0
86	igcm_o=0
87	igcm_o1d=0
88	igcm_o2=0
89	igcm_o3=0
90	igcm_h=0
91	igcm_d=0
92	igcm_hd=0
93	igcm_h2=0
94	igcm_od=0
95	igcm_do2=0
96	igcm_hdo2=0
97	igcm_oh=0
98	igcm_ho2=0
99	igcm_h2o2=0
100	igcm_ch4=0
101	igcm_n2=0
102	igcm_ar=0
103	igcm_ar_n2=0
104	igcm_n=0
105	igcm_no=0
106	igcm_no2=0
107	igcm_n2d=0
108	igcm_he=0
109	igcm_co2plus=0
110	igcm_oplus=0
111	igcm_o2plus=0
112	igcm_coplus=0
113	igcm_cplus=0
114	igcm_nplus=0
115	igcm_noplus=0
116	igcm_n2plus=0
117	igcm_hplus=0
118	igcm_hco2plus=0
119	igcm_hcoplus=0
120	igcm_h2oplus=0
121	igcm_h3oplus=0
122	igcm_ohplus=0
123	igcm_elec=0
124
125	! 1. find dust tracers
126	count=0
127	if (dustbin.gt.0) then
128	do iq=1,nq
129	txt=" "
130	write(txt,'(a4,i2.2)')'dust',count+1
131	if (noms(iq).eq.txt) then
132	count=count+1
133	igcm_dustbin(count)=iq
134	mmol(iq)=100.
135	endif
136	enddo !do iq=1,nq
137	endif ! of if (dustbin.gt.0)
138	if (doubleq) then
139	do iq=1,nq
140	if (noms(iq).eq."dust_mass") then
141	igcm_dust_mass=iq
142	count=count+1
143	endif
144	if (noms(iq).eq."dust_number") then
145	igcm_dust_number=iq
146	count=count+1
147	endif
148	enddo
149	endif ! of if (doubleq)
150	if (microphys) then
151	do iq=1,nq
152	if (noms(iq).eq."ccn_mass") then
153	igcm_ccn_mass=iq
154	count=count+1
155	endif
156	if (noms(iq).eq."ccn_number") then
157	igcm_ccn_number=iq
158	count=count+1
159	endif
160	enddo
161	endif ! of if (microphys)
162	if (submicron) then
163	do iq=1,nq
164	if (noms(iq).eq."dust_submicron") then
165	igcm_dust_submicron=iq
166	mmol(iq)=100.
167	count=count+1
168	endif
169	enddo
170	endif ! of if (submicron)
171	if (rdstorm) then
172	do iq=1,nq
173	if (noms(iq).eq."stormdust_mass") then
174	igcm_stormdust_mass=iq
175	count=count+1
176	endif
177	if (noms(iq).eq."stormdust_number") then
178	igcm_stormdust_number=iq
179	count=count+1
180	endif
181	enddo
182	endif ! of if (rdstorm)
183	if (slpwind) then
184	do iq=1,nq
185	if (noms(iq).eq."topdust_mass") then
186	igcm_topdust_mass=iq
187	count=count+1
188	endif
189	if (noms(iq).eq."topdust_number") then
190	igcm_topdust_number=iq
191	count=count+1
192	endif
193	enddo
194	endif ! of if (slpwind)
195	! 2. find chemistry and water tracers
196	do iq=1,nq
197	if (noms(iq).eq."co2") then
198	igcm_co2=iq
199	mmol(igcm_co2)=44.
200	count=count+1
201	endif
202	if (noms(iq).eq."co") then
203	igcm_co=iq
204	mmol(igcm_co)=28.
205	count=count+1
206	endif
207	if (noms(iq).eq."o") then
208	igcm_o=iq
209	mmol(igcm_o)=16.
210	count=count+1
211	endif
212	if (noms(iq).eq."o1d") then
213	igcm_o1d=iq
214	mmol(igcm_o1d)=16.
215	count=count+1
216	endif
217	if (noms(iq).eq."o2") then
218	igcm_o2=iq
219	mmol(igcm_o2)=32.
220	count=count+1
221	endif
222	if (noms(iq).eq."o3") then
223	igcm_o3=iq
224	mmol(igcm_o3)=48.
225	count=count+1
226	endif
227	if (noms(iq).eq."h") then
228	igcm_h=iq
229	mmol(igcm_h)=1.
230	count=count+1
231	endif
232	if (noms(iq).eq."h2") then
233	igcm_h2=iq
234	mmol(igcm_h2)=2.
235	count=count+1
236	endif
237	if (noms(iq).eq."oh") then
238	igcm_oh=iq
239	mmol(igcm_oh)=17.
240	count=count+1
241	endif
242	if (noms(iq).eq."ho2") then
243	igcm_ho2=iq
244	mmol(igcm_ho2)=33.
245	count=count+1
246	endif
247	if (noms(iq).eq."h2o2") then
248	igcm_h2o2=iq
249	mmol(igcm_h2o2)=34.
250	count=count+1
251	endif
252	if (noms(iq).eq."n2") then
253	igcm_n2=iq
254	mmol(igcm_n2)=28.
255	count=count+1
256	endif
257	if (noms(iq).eq."ch4") then
258	igcm_ch4=iq
259	mmol(igcm_ch4)=16.
260	count=count+1
261	endif
262	if (noms(iq).eq."ar") then
263	igcm_ar=iq
264	mmol(igcm_ar)=40.
265	count=count+1
266	endif
267	if (noms(iq).eq."n") then
268	igcm_n=iq
269	mmol(igcm_n)=14.
270	count=count+1
271	endif
272	if (noms(iq).eq."no") then
273	igcm_no=iq
274	mmol(igcm_no)=30.
275	count=count+1
276	endif
277	if (noms(iq).eq."no2") then
278	igcm_no2=iq
279	mmol(igcm_no2)=46.
280	count=count+1
281	endif
282	if (noms(iq).eq."n2d") then
283	igcm_n2d=iq
284	mmol(igcm_n2d)=28.
285	count=count+1
286	endif
287	if (noms(iq).eq."he") then
288	igcm_he=iq
289	mmol(igcm_he)=4.
290	count=count+1
291	endif
292	if (noms(iq).eq."co2plus") then
293	igcm_co2plus=iq
294	mmol(igcm_co2plus)=44.
295	count=count+1
296	endif
297	if (noms(iq).eq."oplus") then
298	igcm_oplus=iq
299	mmol(igcm_oplus)=16.
300	count=count+1
301	endif
302	if (noms(iq).eq."o2plus") then
303	igcm_o2plus=iq
304	mmol(igcm_o2plus)=32.
305	count=count+1
306	endif
307	if (noms(iq).eq."coplus") then
308	igcm_coplus=iq
309	mmol(igcm_coplus)=28.
310	count=count+1
311	endif
312	if (noms(iq).eq."cplus") then
313	igcm_cplus=iq
314	mmol(igcm_cplus)=12.
315	count=count+1
316	endif
317	if (noms(iq).eq."nplus") then
318	igcm_nplus=iq
319	mmol(igcm_nplus)=14.
320	count=count+1
321	endif
322	if (noms(iq).eq."noplus") then
323	igcm_noplus=iq
324	mmol(igcm_noplus)=30.
325	count=count+1
326	endif
327	if (noms(iq).eq."n2plus") then
328	igcm_n2plus=iq
329	mmol(igcm_n2plus)=28.
330	count=count+1
331	endif
332	if (noms(iq).eq."hplus") then
333	igcm_hplus=iq
334	mmol(igcm_hplus)=1.
335	count=count+1
336	endif
337	if (noms(iq).eq."hco2plus") then
338	igcm_hco2plus=iq
339	mmol(igcm_hco2plus)=45.
340	count=count+1
341	endif
342	if (noms(iq).eq."hcoplus") then
343	igcm_hcoplus=iq
344	mmol(igcm_hcoplus)=29.
345	count=count+1
346	endif
347	if (noms(iq).eq."h2oplus") then
348	igcm_h2oplus=iq
349	mmol(igcm_h2oplus)=18.
350	count=count+1
351	endif
352	if (noms(iq).eq."h3oplus") then
353	igcm_h3oplus=iq
354	mmol(igcm_h3oplus)=19.
355	count=count+1
356	endif
357	if (noms(iq).eq."ohplus") then
358	igcm_ohplus=iq
359	mmol(igcm_ohplus)=17.
360	count=count+1
361	endif
362	if (noms(iq).eq."elec") then
363	igcm_elec=iq
364	mmol(igcm_elec)=1./1822.89
365	count=count+1
366	endif
367	if (noms(iq).eq."h2o_vap") then
368	igcm_h2o_vap=iq
369	mmol(igcm_h2o_vap)=18.
370	count=count+1
371	endif
372	if (noms(iq).eq."hdo_vap") then
373	igcm_hdo_vap=iq
374	mmol(igcm_hdo_vap)=19.
375	count=count+1
376	endif
377	if (noms(iq).eq."od") then
378	igcm_od=iq
379	mmol(igcm_od)=18.
380	count=count+1
381	endif
382	if (noms(iq).eq."d") then
383	igcm_d=iq
384	mmol(igcm_d)=2.
385	count=count+1
386	endif
387	if (noms(iq).eq."hd") then
388	igcm_hd=iq
389	mmol(igcm_hd)=3.
390	count=count+1
391	endif
392	if (noms(iq).eq."do2") then
393	igcm_do2=iq
394	mmol(igcm_do2)=34.
395	count=count+1
396	endif
397	if (noms(iq).eq."hdo2") then
398	igcm_hdo2=iq
399	mmol(igcm_hdo2)=35.
400	count=count+1
401	endif
402	if (noms(iq).eq."co2_ice") then
403	igcm_co2_ice=iq
404	mmol(igcm_co2_ice)=44.
405	count=count+1
406	endif
407	if (noms(iq).eq."h2o_ice") then
408	igcm_h2o_ice=iq
409	mmol(igcm_h2o_ice)=18.
410	count=count+1
411	endif
412	if (noms(iq).eq."hdo_ice") then
413	igcm_hdo_ice=iq
414	mmol(igcm_hdo_ice)=19.
415	count=count+1
416	endif
417	! Other stuff: e.g. for simulations using co2 + neutral gaz
418	if (noms(iq).eq."Ar_N2") then
419	igcm_ar_n2=iq
420	mmol(igcm_ar_n2)=30.
421	count=count+1
422	endif
423	if (co2clouds) then
424	if (noms(iq).eq."ccnco2_mass") then
425	igcm_ccnco2_mass=iq
426	count=count+1
427	endif
428	if (noms(iq).eq."ccnco2_number") then
429	igcm_ccnco2_number=iq
430	count=count+1
431	endif
432	if (meteo_flux) then
433	if (noms(iq).eq."ccnco2_meteor_mass") then
434	igcm_ccnco2_meteor_mass=iq
435	count=count+1
436	endif
437	if (noms(iq).eq."ccnco2_meteor_number") then
438	igcm_ccnco2_meteor_number=iq
439	count=count+1
440	endif
441	end if
442	if (co2useh2o) then
443	if (noms(iq).eq."ccnco2_h2o_number") then
444	igcm_ccnco2_h2o_number=iq
445	count=count+1
446	endif
447	if (noms(iq).eq."ccnco2_h2o_mass_ice") then
448	igcm_ccnco2_h2o_mass_ice=iq
449	count=count+1
450	endif
451	if (noms(iq).eq."ccnco2_h2o_mass_ccn") then
452	igcm_ccnco2_h2o_mass_ccn=iq
453	count=count+1
454	endif
455	end if
456	endif
457	enddo ! of do iq=1,nq
458
459	! check that we identified all tracers:
460	if (count.ne.nq) then
461	write(,) "initracer: found only ",count," tracers"
462	write(,) " expected ",nq
463	do iq=1,count
464	write(,)' ',iq,' ',trim(noms(iq))
465	enddo
466	call abort_physic("initracer","tracer mismatch",1)
467	else
468	write(,) "initracer: found all expected tracers, namely:"
469	do iq=1,nq
470	write(,)' ',iq,' ',trim(noms(iq))
471	enddo
472	endif
473
474	! if water cycle but iceparty=.false., there will nevertheless be
475	! water ice at the surface (iceparty is not used anymore, but this
476	! part is still relevant, as we want to stay compatible with the
477	! older versions).
478	if (water.and.(igcm_h2o_ice.eq.0)) then
479	igcm_h2o_ice=igcm_h2o_vap ! so that qsurf(i_h2o_ice) is identified
480	! even though there is no q(i_h2o_ice)
481	else
482	! surface ice qsurf(i_h2o_ice) was loaded twice by phyetat0,
483	! as qsurf(i_h2o_vap) & as qsurf(i_h2o_ice), so to be clean:
484	if (igcm_h2o_vap.ne.0) then
485	qsurf(1:ngrid,igcm_h2o_vap)=0
486	endif
487	endif
488
489	! Additional test required for HDO
490	! We need to compute some things for H2O before HDO
491	if (hdo) then
492	if (igcm_h2o_vap.gt.igcm_hdo_vap) then
493	call abort_physic("initracer",
494	& "Tracer H2O must be initialized before HDO",1)
495	else if ((nqfils(igcm_h2o_ice).lt.1)
496	& .or. (nqfils(igcm_h2o_vap).lt.1)) then
497	write(,) "HDO must be transported as a son",
498	& " of H2O: complete traceur.def"
499	call abort_physic("initracer","adapt your tracer.def",1)
500	else if ((igcm_hdo_ice.lt.nq-2)
501	& .or. (igcm_hdo_vap.lt.nq-2)) then
502	write(,) "The isotopes (HDO) must be placed at",
503	& " the end of the list in traceur.def"
504	call abort_physic("initracer","adapt your tracer.def",1)
505	endif
506	endif
507
508	c------------------------------------------------------------
509	c Initialize tracers .... (in tracer_mod)
510	c------------------------------------------------------------
511	! start by setting everything to (default) zero
512	rho_q(1:nq)=0 ! tracer density (kg.m-3)
513	radius(1:nq)=0. ! tracer particle radius (m)
514	alpha_lift(1:nq) =0. ! tracer saltation vertical flux/horiz flux ratio (m-1)
515	alpha_devil(1:nq)=0. ! tracer lifting coefficient by dust devils
516
517
518	! some reference values
519	rho_dust=2500. ! Mars dust density (kg.m-3)
520	rho_ice=920. ! Water ice density (kg.m-3)
521	rho_ice_co2=1650.
522	!Mangan et al., Icarus 2017 :CO2 density = 1.72391-2.53×10−4T – 2.87×10−6T^2
523	nuice_ref=0.1 ! Effective variance nueff of the
524	! water-ice size distribution
525	!!!nuice_sed=0.45 ! Sedimentation effective variance
526	! of the water-ice size distribution
527
528	nuiceco2_ref=0.2 !C.M. Effective variance nueff of the
529	! co2-ice size distribution
530	if (doubleq) then
531	c "doubleq" technique
532	c -------------------
533	c (transport of mass and number mixing ratio)
534	c iq=1: Q mass mixing ratio, iq=2: N number mixing ratio
535
536	if( (nq.lt.2).or.(water.and.(nq.lt.4))
537	* .or.(hdo.and.(nq.lt.6) )) then
538	write(,)'initracer: nq is too low : nq=', nq
539	write(,)'water= ',water,' doubleq= ',doubleq
540	end if
541
542	nueff_lift = 0.5
543	varian=sqrt(log(1.+nueff_lift))
544
545	rho_q(igcm_dust_mass)=rho_dust
546	rho_q(igcm_dust_number)=rho_dust
547
548	c Intermediate calcul for computing geometric mean radius r0
549	c as a function of mass and number mixing ratio Q and N
550	c (r0 = (r3n_q * Q/ N)^(1/3))
551	r3n_q = exp(-4.5varian2)(3./4.)/(pi*rho_dust)
552
553	c Intermediate calcul for computing effective radius reff
554	c from geometric mean radius r0
555	c (reff = ref_r0 * r0)
556	ref_r0 = exp(2.5varian*2)
557
558	c lifted dust :
559	c '''''''''''
560	reff_lift = 3.0e-6 !3.e-6 !Effective radius of lifted dust (m)
561	alpha_devil(igcm_dust_mass)=9.e-9 ! dust devil lift mass coeff
562	c alpha_lift(igcm_dust_mass)=3.0e-15 ! Lifted mass coeff
563
564	!! default lifting settings
565	!! -- GCM: alpha_lift not zero because large-scale lifting by default
566	!! -- MESOSCALE: alpha_lift zero because no lifting at all in mesoscale by default
567	#ifdef MESOSCALE
568	alpha_lift(igcm_dust_mass)=0.0
569	#else
570	alpha_lift(igcm_dust_mass)=1.e-6 !1.e-6 !Lifted mass coeff
571	IF (dustinjection.ge.1) THEN
572	reff_lift = 3.0e-6 ! Effective radius of lifted dust (m)
573	alpha_lift(igcm_dust_mass)=(4/3.)reff_liftrho_dust
574	& /2.4
575	ENDIF
576	#endif
577
578	r0_lift = reff_lift/ref_r0
579	alpha_devil(igcm_dust_number)=r3n_q*
580	& alpha_devil(igcm_dust_mass)/r0_lift**3
581	alpha_lift(igcm_dust_number)=r3n_q*
582	& alpha_lift(igcm_dust_mass)/r0_lift**3
583
584	radius(igcm_dust_mass) = reff_lift
585	radius(igcm_dust_number) = reff_lift
586
587	write(,) "initracer: doubleq_param reff_lift:", reff_lift
588	write(,) "initracer: doubleq_param nueff_lift:", nueff_lift
589	write(,) "initracer: doubleq_param alpha_lift:",
590	& alpha_lift(igcm_dust_mass)
591	!c ----------------------------------------------------------------------
592	!c rocket dust storm scheme
593	!c lifting tracer stormdust using same distribution than
594	!c normal dust
595	if (rdstorm) then
596	reff_storm=3.e-6 ! reff_lift !3.e-6
597	r0_storm=reff_storm/ref_r0
598	rho_q(igcm_stormdust_mass)=rho_dust
599	rho_q(igcm_stormdust_number)=rho_dust
600
601	alpha_devil(igcm_stormdust_mass)=9.e-9 ! dust devil lift mass coeff
602	alpha_lift(igcm_stormdust_mass)=4./3./2.4reff_stormrho_dust
603
604	write(,) 'alpha_lift(rds):',alpha_lift(igcm_stormdust_mass)
605
606	alpha_devil(igcm_stormdust_number)=r3n_q*
607	& alpha_devil(igcm_stormdust_mass)/r0_storm**3
608	alpha_lift(igcm_stormdust_number)=r3n_q*
609	& alpha_lift(igcm_stormdust_mass)/r0_storm**3
610
611	radius(igcm_stormdust_mass) = reff_storm
612	radius(igcm_stormdust_number) = reff_storm
613	end if !(rdstorm)
614	!c ----------------------------------------------------------------------
615	!c slope wind scheme
616	!c you need a radius value for topdust to active its sedimentation
617	!c we take the same value as for the normal dust
618	if (slpwind) then
619	rho_q(igcm_topdust_mass)=rho_dust
620	rho_q(igcm_topdust_number)=rho_dust
621	radius(igcm_topdust_mass) = 3.e-6
622	radius(igcm_topdust_number) = 3.e-6
623	end if !(slpwind)
624	!c ----------------------------------------------------------------------
625
626	else
627
628	! initialize varian, which may be used (e.g. by surfacearea)
629	! even with conrath dust
630	nueff_lift = 0.5
631	varian=sqrt(log(1.+nueff_lift))
632
633	if (dustbin.gt.1) then
634	print*,'initracer: STOP!',
635	$ ' properties of dust need to be set in initracer !!!'
636	call abort_physic("initracer","dustbin properties issue",1)
637
638	else if (dustbin.eq.1) then
639
640	c This will be used for 1 dust particle size:
641	c ------------------------------------------
642	radius(igcm_dustbin(1))=3.e-6
643	alpha_lift(igcm_dustbin(1))=0.0e-6
644	alpha_devil(igcm_dustbin(1))=7.65e-9
645	rho_q(igcm_dustbin(1))=rho_dust
646
647	endif
648	end if ! (doubleq)
649
650
651	c Scavenging of dust particles by H2O clouds:
652	c ------------------------------------------
653	c Initialize the two tracers used for the CCNs
654	if (water.AND.doubleq.AND.scavenging) then
655	radius(igcm_ccn_mass) = radius(igcm_dust_mass)
656	alpha_lift(igcm_ccn_mass) = 1e-30
657	alpha_devil(igcm_ccn_mass) = 1e-30
658	rho_q(igcm_ccn_mass) = rho_dust
659
660	radius(igcm_ccn_number) = radius(igcm_ccn_mass)
661	alpha_lift(igcm_ccn_number) = alpha_lift(igcm_ccn_mass)
662	alpha_devil(igcm_ccn_number) = alpha_devil(igcm_ccn_mass)
663	rho_q(igcm_ccn_number) = rho_q(igcm_ccn_mass)
664	endif ! of if (water.AND.doubleq.AND.scavenging)
665
666	c Submicron dust mode:
667	c --------------------
668
669	if (submicron) then
670	radius(igcm_dust_submicron)=0.1e-6
671	rho_q(igcm_dust_submicron)=rho_dust
672	if (doubleq) then
673	c If doubleq is also active, we use the population ratio:
674	alpha_lift(igcm_dust_submicron) =
675	& alpha_lift(igcm_dust_number)popratio
676	& rho_q(igcm_dust_submicron)4./3.pi*
677	& radius(igcm_dust_submicron)**3.
678	alpha_devil(igcm_dust_submicron)=1.e-30
679	else
680	alpha_lift(igcm_dust_submicron)=1e-6
681	alpha_devil(igcm_dust_submicron)=1.e-30
682	endif ! (doubleq)
683	end if ! (submicron)
684
685	c Initialization for water vapor
686	c ------------------------------
687	if(water) then
688	radius(igcm_h2o_vap)=0.
689	alpha_lift(igcm_h2o_vap) =0.
690	alpha_devil(igcm_h2o_vap)=0.
691	if(water.and.(nq.ge.2)) then
692	radius(igcm_h2o_ice)=3.e-6
693	rho_q(igcm_h2o_ice)=rho_ice
694	alpha_lift(igcm_h2o_ice) =0.
695	alpha_devil(igcm_h2o_ice)=0.
696	elseif(water.and.(nq.lt.2)) then
697	write(,) 'nq is too low : nq=', nq
698	write(,) 'water= ',water
699	endif
700
701	end if ! (water)
702
703	c Initialization for hdo vapor
704	c ------------------------------
705	if (hdo) then
706	radius(igcm_hdo_vap)=0.
707	alpha_lift(igcm_hdo_vap) =0.
708	alpha_devil(igcm_hdo_vap)=0.
709	if(water.and.(nq.ge.2)) then
710	radius(igcm_hdo_ice)=3.e-6
711	rho_q(igcm_hdo_ice)=rho_ice
712	alpha_lift(igcm_hdo_ice) =0.
713	alpha_devil(igcm_hdo_ice)=0.
714	elseif(hdo.and.(nq.lt.6)) then
715	write(,) 'nq is too low : nq=', nq
716	write(,) 'hdo= ',hdo
717	endif
718
719	end if ! (hdo)
720
721
722	! Initialisation for CO2 clouds
723	if (co2clouds ) then
724	radius(igcm_ccnco2_mass) = radius(igcm_dust_mass)
725	alpha_lift(igcm_ccnco2_mass) = 1e-30
726	alpha_devil(igcm_ccnco2_mass) = 1e-30
727	rho_q(igcm_ccnco2_mass) = rho_dust
728	radius(igcm_ccnco2_number) = radius(igcm_ccnco2_mass)
729	alpha_lift(igcm_ccnco2_number) = alpha_lift(igcm_ccnco2_mass)
730	alpha_devil(igcm_ccnco2_number) = alpha_devil(igcm_ccnco2_mass)
731	rho_q(igcm_ccnco2_number) = rho_q(igcm_ccnco2_mass)
732
733	radius(igcm_co2)=0.
734	alpha_lift(igcm_co2) =0.
735	alpha_devil(igcm_co2)=0.
736	radius(igcm_co2_ice)=1.e-8
737	rho_q(igcm_co2_ice)=rho_ice_co2
738	alpha_lift(igcm_co2_ice) =0.
739	alpha_devil(igcm_co2_ice)=0.
740
741	endif
742
743	c Output for records:
744	c ~~~~~~~~~~~~~~~~~~
745	write(,)
746	Write(,) '******** initracer : dust transport parameters :'
747	write(,) 'alpha_lift = ', alpha_lift
748	write(,) 'alpha_devil = ', alpha_devil
749	write(,) 'radius = ', radius
750	if(doubleq) then
751	write(,) 'reff_lift (um) = ', reff_lift
752	write(,) 'size distribution variance = ', varian
753	write(,) 'r3n_q , ref_r0 : ', r3n_q , ref_r0
754	end if
755
756	!
757	! some extra (possibly redundant) sanity checks for tracers:
758	! ---------------------------------------------------------
759
760	if (doubleq) then
761	! verify that we indeed have dust_mass and dust_number tracers
762	if (igcm_dust_mass.eq.0) then
763	write(,) "initracer: error !!"
764	write(,) " cannot use doubleq option without ",
765	& "a dust_mass tracer !"
766	call abort_physic("initracer","doubleq issue",1)
767	endif
768	if (igcm_dust_number.eq.0) then
769	write(,) "initracer: error !!"
770	write(,) " cannot use doubleq option without ",
771	& "a dust_number tracer !"
772	call abort_physic("initracer","doubleq issue",1)
773	endif
774	endif
775
776	if ((.not.doubleq).and.(dustbin.gt.0)) then
777	! verify that we indeed have 'dustbin' dust tracers
778	count=0
779	do iq=1,dustbin
780	if (igcm_dustbin(iq).ne.0) then
781	count=count+1
782	endif
783	enddo
784	if (count.ne.dustbin) then
785	write(,) "initracer: error !!"
786	write(,) " dustbin is set to ",dustbin,
787	& " but we only have the following dust tracers:"
788	do iq=1,count
789	write(,)" ",trim(noms(igcm_dustbin(iq)))
790	enddo
791	call abort_physic("initracer","dustbin issue",1)
792	endif
793	endif
794
795	if (water) then
796	! verify that we indeed have h2o_vap and h2o_ice tracers
797	if (igcm_h2o_vap.eq.0) then
798	write(,) "initracer: error !!"
799	write(,) " cannot use water option without ",
800	& "an h2o_vap tracer !"
801	call abort_physic("initracer","water cycle issue",1)
802	endif
803	if (igcm_h2o_ice.eq.0) then
804	write(,) "initracer: error !!"
805	write(,) " cannot use water option without ",
806	& "an h2o_ice tracer !"
807	call abort_physic("initracer","water cycle issue",1)
808	endif
809	endif
810
811	if (hdo) then
812	! verify that we indeed have hdo_vap and hdo_ice tracers
813	if (igcm_hdo_vap.eq.0) then
814	write(,) "initracer: error !!"
815	write(,) " cannot use hdo option without ",
816	& "an hdo_vap tracer !"
817	call abort_physic("initracer","hdo cycle issue",1)
818	endif
819	if (igcm_hdo_ice.eq.0) then
820	write(,) "initracer: error !!"
821	write(,) " cannot use hdo option without ",
822	& "an hdo_ice tracer !"
823	call abort_physic("initracer","hdo cycle issue",1)
824	endif
825	endif
826
827
828	if (co2clouds) then
829	!verify that we have co2_ice and co2 tracers
830	if (igcm_co2 .eq. 0) then
831	write(,) "initracer: error !!"
832	write(,) " cannot use co2 clouds option without ",
833	& "a co2 tracer !"
834	call abort_physic("initracer","co2 clouds issue",1)
835	end if
836	if (igcm_co2_ice .eq. 0) then
837	write(,) "initracer: error !!"
838	write(,) " cannot use co2 clouds option without ",
839	& "a co2_ice tracer !"
840	call abort_physic("initracer","co2 clouds issue",1)
841	end if
842	! TODO: ajouter les tests avec co2useh2o, et meteo_flux, ainsi que les traceurs ccnco2_number
843	endif
844
845	if (rdstorm) then
846	! verify that we indeed have stormdust_mass and stormdust_number tracers
847	if (igcm_stormdust_mass.eq.0) then
848	write(,) "initracer: error !!"
849	write(,) " cannot use rdstorm option without ",
850	& "a stormdust_mass tracer !"
851	call abort_physic("initracer","rdstorm issue",1)
852	endif
853	if (igcm_stormdust_number.eq.0) then
854	write(,) "initracer: error !!"
855	write(,) " cannot use rdstorm option without ",
856	& "a stormdust_number tracer !"
857	call abort_physic("initracer","rdstorm issue",1)
858	endif
859	endif
860
861	if (slpwind) then
862	! verify that we indeed have topdust_mass and topdust_number tracers
863	if (igcm_topdust_mass.eq.0) then
864	write(,) "initracer: error !!"
865	write(,) " cannot use slpwind option without ",
866	& "a topdust_mass tracer !"
867	call abort_physic("initracer","slpwind issue",1)
868	endif
869	if (igcm_topdust_number.eq.0) then
870	write(,) "initracer: error !!"
871	write(,) " cannot use slpwind option without ",
872	& "a topdust_number tracer !"
873	call abort_physic("initracer","slpwind issue",1)
874	endif
875	endif
876
877	if (callnlte) then ! NLTE requirements
878	if (nltemodel.ge.1) then
879	! check that co2, co, o and n2 tracers are available
880	if (igcm_co2.eq.0) then
881	write(,) "initracer: error !!"
882	write(,) " with nltemodel>0, we need the co2 tracer!"
883	call abort_physic("initracer","missing co2 tracer",1)
884	endif
885	if (igcm_co.eq.0) then
886	write(,) "initracer: error !!"
887	write(,) " with nltemodel>0, we need the co tracer!"
888	call abort_physic("initracer","missing co tracer",1)
889	endif
890	if (igcm_o.eq.0) then
891	write(,) "initracer: error !!"
892	write(,) " with nltemodel>0, we need the o tracer!"
893	call abort_physic("initracer","missing o tracer",1)
894	endif
895	if (igcm_n2.eq.0) then
896	write(,) "initracer: error !!"
897	write(,) " with nltemodel>0, we need the n2 tracer!"
898	call abort_physic("initracer","missing n2 tracer",1)
899	endif
900	endif
901	endif
902
903	if (scavenging) then
904	! verify that we indeed have ccn_mass and ccn_number tracers
905	if (igcm_ccn_mass.eq.0 .and. igcm_ccnco2_mass.eq.0) then
906	write(,) "initracer: error !!"
907	write(,) " cannot use scavenging option without ",
908	& "a ccn_mass or ccnco2_mass tracer !"
909	call abort_physic("initracer","scavenging issue",1)
910	endif
911	if (igcm_ccn_number.eq.0 .and. igcm_ccnco2_number.eq.0 ) then
912	write(,) "initracer: error !!"
913	write(,) " cannot use scavenging option without ",
914	& "a ccn_number or ccnco2_number tracer !"
915	call abort_physic("initracer","scavenging issue",1)
916	endif
917	endif ! of if (scavenging)
918
919	if (photochem .or. callthermos) then
920	! verify that we indeed have the chemistry tracers
921	if (igcm_co2.eq.0) then
922	write(,) "initracer: error !!"
923	write(,) " cannot use chemistry option without ",
924	& "a co2 tracer !"
925	call abort_physic("initracer","missing co2 tracer",1)
926	endif
927	if (igcm_co.eq.0) then
928	write(,) "initracer: error !!"
929	write(,) " cannot use chemistry option without ",
930	& "a co tracer !"
931	call abort_physic("initracer","missing co tracer",1)
932	endif
933	if (igcm_o.eq.0) then
934	write(,) "initracer: error !!"
935	write(,) " cannot use chemistry option without ",
936	& "a o tracer !"
937	call abort_physic("initracer","missing o tracer",1)
938	endif
939	if (igcm_o1d.eq.0) then
940	write(,) "initracer: error !!"
941	write(,) " cannot use chemistry option without ",
942	& "a o1d tracer !"
943	call abort_physic("initracer","missing o1d tracer",1)
944	endif
945	if (igcm_o2.eq.0) then
946	write(,) "initracer: error !!"
947	write(,) " cannot use chemistry option without ",
948	& "an o2 tracer !"
949	call abort_physic("initracer","missing o2 tracer",1)
950	endif
951	if (igcm_o3.eq.0) then
952	write(,) "initracer: error !!"
953	write(,) " cannot use chemistry option without ",
954	& "an o3 tracer !"
955	call abort_physic("initracer","missing o3 tracer",1)
956	endif
957	if (igcm_h.eq.0) then
958	write(,) "initracer: error !!"
959	write(,) " cannot use chemistry option without ",
960	& "a h tracer !"
961	call abort_physic("initracer","missing h tracer",1)
962	endif
963	if (igcm_h2.eq.0) then
964	write(,) "initracer: error !!"
965	write(,) " cannot use chemistry option without ",
966	& "a h2 tracer !"
967	call abort_physic("initracer","missing h2 tracer",1)
968	endif
969	if (igcm_oh.eq.0) then
970	write(,) "initracer: error !!"
971	write(,) " cannot use chemistry option without ",
972	& "an oh tracer !"
973	call abort_physic("initracer","missing oh tracer",1)
974	endif
975	if (igcm_ho2.eq.0) then
976	write(,) "initracer: error !!"
977	write(,) " cannot use chemistry option without ",
978	& "a ho2 tracer !"
979	call abort_physic("initracer","missing ho2 tracer",1)
980	endif
981	if (igcm_h2o2.eq.0) then
982	write(,) "initracer: error !!"
983	write(,) " cannot use chemistry option without ",
984	& "a h2o2 tracer !"
985	call abort_physic("initracer","missing h2o2 tracer",1)
986	endif
987	if (igcm_n2.eq.0) then
988	write(,) "initracer: error !!"
989	write(,) " cannot use chemistry option without ",
990	& "a n2 tracer !"
991	call abort_physic("initracer","missing n2 tracer",1)
992	endif
993	if (igcm_ar.eq.0) then
994	write(,) "initracer: error !!"
995	write(,) " cannot use chemistry option without ",
996	& "an ar tracer !"
997	call abort_physic("initracer","missing ar tracer",1)
998	endif
999	endif ! of if (photochem .or. callthermos)
1000
1001	end

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

Download in other formats:

Original Format