Context Navigation

initracer.F @ 3741

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

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

Download in other formats:

Original Format