Context Navigation

initracer.F @ 648

Last change on this file since 648 was 648, checked in by emillour, 13 years ago

Mars GCM:

some syntax corrections in thermcall_main_mars, vdif_cd, pbl_parameters which cause problems when compiling with some strict compilers (g95, gfortran)
added an initialisation of 'varian' in initracer for cases when using conrath dust; because that value can be is used elsewhere (e.g. surfacearea)

JYC+EM

File size: 21.5 KB

Line
1	SUBROUTINE initracer(qsurf,co2ice)
2
3	IMPLICIT NONE
4	c=======================================================================
5	c subject:
6	c --------
7	c Initialization related to tracer
8	c (transported dust, water, chemical species, ice...)
9	c
10	c Name of the tracer
11	c
12	c Test of dimension :
13	c Initialize COMMON tracer in tracer.h, using tracer names provided
14	c by the dynamics in "advtrac.h"
15	c
16	c Old conventions: (not used any more)
17	c
18	c If water=T : q(iq=nqmx) is the water mass mixing ratio
19	c and q(iq=nqmx-1) is the ice mass mixing ratio
20
21	c If there is transported dust, it uses iq=1 to iq=dustbin
22	c If there is no transported dust : dustbin=0
23	c If doubleq=T : q(iq=1) is the dust mass mixing ratio
24	c q(iq=2) is the dust number mixing ratio
25
26	c
27	c author: F.Forget
28	c ------
29	c Modifs: Franck Montmessin, Sebastien Lebonnois (june 2003)
30	c Ehouarn Millour (oct. 2008) identify tracers by their names
31	c=======================================================================
32
33
34	#include "dimensions.h"
35	#include "dimphys.h"
36	#include "comcstfi.h"
37	#include "callkeys.h"
38	#include "tracer.h"
39	#include "advtrac.h"
40	#include "comgeomfi.h"
41
42	#include "surfdat.h"
43
44	real qsurf(ngridmx,nqmx) ! tracer on surface (e.g. kg.m-2)
45	real co2ice(ngridmx) ! co2 ice mass on surface (e.g. kg.m-2)
46	integer iq,ig,count
47	real r0_lift , reff_lift, nueff_lift
48	c Ratio of small over large dust particles (used when both
49	c doubleq and the submicron mode are active); In Montmessin
50	c et al. (2002), a value of 25 has been deduced;
51	real, parameter :: popratio = 25.
52	character(len=20) :: txt ! to store some text
53
54	c-----------------------------------------------------------------------
55	c radius(nqmx) ! aerosol particle radius (m)
56	c rho_q(nqmx) ! tracer densities (kg.m-3)
57	c alpha_lift(nqmx) ! saltation vertical flux/horiz flux ratio (m-1)
58	c alpha_devil(nqmx) ! lifting coeeficient by dust devil
59	c rho_dust ! Mars dust density
60	c rho_ice ! Water ice density
61	c nuice_ref ! Effective variance nueff of the
62	c ! water-ice size distributions
63	c doubleq ! if method with mass (iq=1) and number(iq=2) mixing ratio
64	c varian ! Characteristic variance of log-normal distribution
65	c-----------------------------------------------------------------------
66
67	! Initialization: get tracer names from the dynamics and check if we are
68	! using 'old' tracer convention ('q01',q02',...)
69	! or new convention (full tracer names)
70	! check if tracers have 'old' names
71	count=0
72	do iq=1,nqmx
73	txt=" "
74	write(txt,'(a1,i2.2)') 'q',iq
75	if (txt.eq.tnom(iq)) then
76	count=count+1
77	endif
78	enddo ! of do iq=1,nqmx
79
80	if (count.eq.nqmx) then
81	write(,) "initracer: tracers seem to follow old naming ",
82	& "convention (q01,q02,...)"
83	write(,) "you should run newstart to rename them"
84	stop
85	endif
86
87	! copy tracer names from dynamics
88	do iq=1,nqmx
89	noms(iq)=tnom(iq)
90	enddo
91
92	c------------------------------------------------------------
93	c NAME and molar mass of the tracer
94	c------------------------------------------------------------
95
96	! Identify tracers by their names: (and set corresponding values of mmol)
97	! 0. initialize tracer indexes to zero:
98	do iq=1,nqmx
99	igcm_dustbin(iq)=0
100	enddo
101	igcm_dust_mass=0
102	igcm_dust_number=0
103	igcm_ccn_mass=0
104	igcm_ccn_number=0
105	igcm_dust_submicron=0
106	igcm_h2o_vap=0
107	igcm_h2o_ice=0
108	igcm_co2=0
109	igcm_co=0
110	igcm_o=0
111	igcm_o1d=0
112	igcm_o2=0
113	igcm_o3=0
114	igcm_h=0
115	igcm_h2=0
116	igcm_oh=0
117	igcm_ho2=0
118	igcm_h2o2=0
119	igcm_ch4=0
120	igcm_n2=0
121	igcm_ar=0
122	igcm_ar_n2=0
123	igcm_n=0
124	igcm_no=0
125	igcm_no2=0
126	igcm_n2d=0
127	igcm_co2plus=0
128	igcm_oplus=0
129	igcm_o2plus=0
130	igcm_coplus=0
131	igcm_cplus=0
132	igcm_nplus=0
133	igcm_noplus=0
134	igcm_n2plus=0
135	igcm_hplus=0
136	igcm_hco2plus=0
137	igcm_elec=0
138
139	! 1. find dust tracers
140	count=0
141	if (dustbin.gt.0) then
142	do iq=1,nqmx
143	txt=" "
144	write(txt,'(a4,i2.2)')'dust',count+1
145	if (noms(iq).eq.txt) then
146	count=count+1
147	igcm_dustbin(count)=iq
148	mmol(iq)=100.
149	endif
150	enddo !do iq=1,nqmx
151	endif ! of if (dustbin.gt.0)
152	if (doubleq) then
153	do iq=1,nqmx
154	if (noms(iq).eq."dust_mass") then
155	igcm_dust_mass=iq
156	count=count+1
157	endif
158	if (noms(iq).eq."dust_number") then
159	igcm_dust_number=iq
160	count=count+1
161	endif
162	enddo
163	endif ! of if (doubleq)
164	if (scavenging) then
165	do iq=1,nqmx
166	if (noms(iq).eq."ccn_mass") then
167	igcm_ccn_mass=iq
168	count=count+1
169	endif
170	if (noms(iq).eq."ccn_number") then
171	igcm_ccn_number=iq
172	count=count+1
173	endif
174	enddo
175	endif ! of if (scavenging)
176	if (submicron) then
177	do iq=1,nqmx
178	if (noms(iq).eq."dust_submicron") then
179	igcm_dust_submicron=iq
180	mmol(iq)=100.
181	count=count+1
182	endif
183	enddo
184	endif ! of if (submicron)
185	! 2. find chemistry and water tracers
186	do iq=1,nqmx
187	if (noms(iq).eq."co2") then
188	igcm_co2=iq
189	mmol(igcm_co2)=44.
190	count=count+1
191	endif
192	if (noms(iq).eq."co") then
193	igcm_co=iq
194	mmol(igcm_co)=28.
195	count=count+1
196	endif
197	if (noms(iq).eq."o") then
198	igcm_o=iq
199	mmol(igcm_o)=16.
200	count=count+1
201	endif
202	if (noms(iq).eq."o1d") then
203	igcm_o1d=iq
204	mmol(igcm_o1d)=16.
205	count=count+1
206	endif
207	if (noms(iq).eq."o2") then
208	igcm_o2=iq
209	mmol(igcm_o2)=32.
210	count=count+1
211	endif
212	if (noms(iq).eq."o3") then
213	igcm_o3=iq
214	mmol(igcm_o3)=48.
215	count=count+1
216	endif
217	if (noms(iq).eq."h") then
218	igcm_h=iq
219	mmol(igcm_h)=1.
220	count=count+1
221	endif
222	if (noms(iq).eq."h2") then
223	igcm_h2=iq
224	mmol(igcm_h2)=2.
225	count=count+1
226	endif
227	if (noms(iq).eq."oh") then
228	igcm_oh=iq
229	mmol(igcm_oh)=17.
230	count=count+1
231	endif
232	if (noms(iq).eq."ho2") then
233	igcm_ho2=iq
234	mmol(igcm_ho2)=33.
235	count=count+1
236	endif
237	if (noms(iq).eq."h2o2") then
238	igcm_h2o2=iq
239	mmol(igcm_h2o2)=34.
240	count=count+1
241	endif
242	if (noms(iq).eq."n2") then
243	igcm_n2=iq
244	mmol(igcm_n2)=28.
245	count=count+1
246	endif
247	if (noms(iq).eq."ch4") then
248	igcm_ch4=iq
249	mmol(igcm_ch4)=16.
250	count=count+1
251	endif
252	if (noms(iq).eq."ar") then
253	igcm_ar=iq
254	mmol(igcm_ar)=40.
255	count=count+1
256	endif
257	if (noms(iq).eq."n") then
258	igcm_n=iq
259	mmol(igcm_n)=14.
260	count=count+1
261	endif
262	if (noms(iq).eq."no") then
263	igcm_no=iq
264	mmol(igcm_no)=30.
265	count=count+1
266	endif
267	if (noms(iq).eq."no2") then
268	igcm_no2=iq
269	mmol(igcm_no2)=46.
270	count=count+1
271	endif
272	if (noms(iq).eq."n2d") then
273	igcm_n2d=iq
274	mmol(igcm_n2d)=28.
275	count=count+1
276	endif
277	if (noms(iq).eq."co2plus") then
278	igcm_co2plus=iq
279	mmol(igcm_co2plus)=44.
280	count=count+1
281	endif
282	if (noms(iq).eq."oplus") then
283	igcm_oplus=iq
284	mmol(igcm_oplus)=16.
285	count=count+1
286	endif
287	if (noms(iq).eq."o2plus") then
288	igcm_o2plus=iq
289	mmol(igcm_o2plus)=32.
290	count=count+1
291	endif
292	if (noms(iq).eq."coplus") then
293	igcm_coplus=iq
294	mmol(igcm_coplus)=28.
295	count=count+1
296	endif
297	if (noms(iq).eq."cplus") then
298	igcm_cplus=iq
299	mmol(igcm_cplus)=12.
300	count=count+1
301	endif
302	if (noms(iq).eq."nplus") then
303	igcm_nplus=iq
304	mmol(igcm_nplus)=14.
305	count=count+1
306	endif
307	if (noms(iq).eq."noplus") then
308	igcm_noplus=iq
309	mmol(igcm_noplus)=30.
310	count=count+1
311	endif
312	if (noms(iq).eq."n2plus") then
313	igcm_n2plus=iq
314	mmol(igcm_n2plus)=28.
315	count=count+1
316	endif
317	if (noms(iq).eq."hplus") then
318	igcm_hplus=iq
319	mmol(igcm_hplus)=1.
320	count=count+1
321	endif
322	if (noms(iq).eq."hco2plus") then
323	igcm_hco2plus=iq
324	mmol(igcm_hco2plus)=45.
325	count=count+1
326	endif
327	if (noms(iq).eq."elec") then
328	igcm_elec=iq
329	mmol(igcm_elec)=1./1822.89
330	count=count+1
331	endif
332	if (noms(iq).eq."h2o_vap") then
333	igcm_h2o_vap=iq
334	mmol(igcm_h2o_vap)=18.
335	count=count+1
336	endif
337	if (noms(iq).eq."h2o_ice") then
338	igcm_h2o_ice=iq
339	mmol(igcm_h2o_ice)=18.
340	count=count+1
341	endif
342	! Other stuff: e.g. for simulations using co2 + neutral gaz
343	if (noms(iq).eq."Ar_N2") then
344	igcm_ar_n2=iq
345	mmol(igcm_ar_n2)=30.
346	count=count+1
347	endif
348
349	enddo ! of do iq=1,nqmx
350
351	! check that we identified all tracers:
352	if (count.ne.nqmx) then
353	write(,) "initracer: found only ",count," tracers"
354	write(,) " expected ",nqmx
355	do iq=1,count
356	write(,)' ',iq,' ',trim(noms(iq))
357	enddo
358	stop
359	else
360	write(,) "initracer: found all expected tracers, namely:"
361	do iq=1,nqmx
362	write(,)' ',iq,' ',trim(noms(iq))
363	enddo
364	endif
365
366	! if water cycle but iceparty=.false., there will nevertheless be
367	! water ice at the surface (iceparty is not used anymore, but this
368	! part is still relevant, as we want to stay compatible with the
369	! older versions).
370	if (water.and.(igcm_h2o_ice.eq.0)) then
371	igcm_h2o_ice=igcm_h2o_vap ! so that qsurf(i_h2o_ice) is identified
372	! even though there is no q(i_h2o_ice)
373	else
374	! surface ice qsurf(i_h2o_ice) was loaded twice by phyetat0,
375	! as qsurf(i_h2o_vap) & as qsurf(i_h2o_ice), so to be clean:
376	if (igcm_h2o_vap.ne.0) then
377	qsurf(1:ngridmx,igcm_h2o_vap)=0
378	endif
379	endif
380
381	c------------------------------------------------------------
382	c Initialisation tracers ....
383	c------------------------------------------------------------
384	call zerophys(nqmx,rho_q)
385
386	rho_dust=2500. ! Mars dust density (kg.m-3)
387	rho_ice=920. ! Water ice density (kg.m-3)
388	nuice_ref=0.1 ! Effective variance nueff of the
389	! water-ice size distribution
390	!!!nuice_sed=0.45 ! Sedimentation effective variance
391	! of the water-ice size distribution
392
393	if (doubleq) then
394	c "doubleq" technique
395	c -------------------
396	c (transport of mass and number mixing ratio)
397	c iq=1: Q mass mixing ratio, iq=2: N number mixing ratio
398
399	if( (nqmx.lt.2).or.(water.and.(nqmx.lt.4)) ) then
400	write(,)'initracer: nqmx is too low : nqmx=', nqmx
401	write(,)'water= ',water,' doubleq= ',doubleq
402	end if
403
404	nueff_lift = 0.5
405	varian=sqrt(log(1.+nueff_lift))
406
407	rho_q(igcm_dust_mass)=rho_dust
408	rho_q(igcm_dust_number)=rho_dust
409
410	c Intermediate calcul for computing geometric mean radius r0
411	c as a function of mass and number mixing ratio Q and N
412	c (r0 = (r3n_q * Q/ N)^(1/3))
413	r3n_q = exp(-4.5varian2)(3./4.)/(pi*rho_dust)
414
415	c Intermediate calcul for computing effective radius reff
416	c from geometric mean radius r0
417	c (reff = ref_r0 * r0)
418	ref_r0 = exp(2.5varian*2)
419
420	c lifted dust :
421	c '''''''''''
422	reff_lift = 3.0e-6 !3.e-6 !Effective radius of lifted dust (m)
423	alpha_devil(igcm_dust_mass)=9.e-9 ! dust devil lift mass coeff
424	c alpha_lift(igcm_dust_mass)=3.0e-15 ! Lifted mass coeff
425	alpha_lift(igcm_dust_mass)=1.e-6 !1.e-6 !Lifted mass coeff
426
427	r0_lift = reff_lift/ref_r0
428	alpha_devil(igcm_dust_number)=r3n_q*
429	& alpha_devil(igcm_dust_mass)/r0_lift**3
430	alpha_lift(igcm_dust_number)=r3n_q*
431	& alpha_lift(igcm_dust_mass)/r0_lift**3
432
433	radius(igcm_dust_mass) = reff_lift
434	radius(igcm_dust_number) = reff_lift
435
436	write(,) "initracer: doubleq_param reff_lift:", reff_lift
437	write(,) "initracer: doubleq_param nueff_lift:", nueff_lift
438	write(,) "initracer: doubleq_param alpha_lift:",
439	& alpha_lift(igcm_dust_mass)
440	else
441
442	! initialize varian, which may be used (e.g. by surfacearea)
443	! even with conrath dust
444	nueff_lift = 0.5
445	varian=sqrt(log(1.+nueff_lift))
446
447	if (dustbin.gt.1) then
448	print*,'initracer: STOP!',
449	$ ' properties of dust need to be set in initracer !!!'
450	stop
451
452	else if (dustbin.eq.1) then
453
454	c This will be used for 1 dust particle size:
455	c ------------------------------------------
456	radius(igcm_dustbin(1))=3.e-6
457	alpha_lift(igcm_dustbin(1))=0.0e-6
458	alpha_devil(igcm_dustbin(1))=7.65e-9
459	rho_q(igcm_dustbin(1))=rho_dust
460
461	endif
462	end if ! (doubleq)
463
464
465	c Scavenging of dust particles by H2O clouds:
466	c ------------------------------------------
467	c Initialize the two tracers used for the CCNs
468	if (water.AND.doubleq.AND.scavenging) then
469	radius(igcm_ccn_mass) = radius(igcm_dust_mass)
470	alpha_lift(igcm_ccn_mass) = 1e-30
471	alpha_devil(igcm_ccn_mass) = 1e-30
472	rho_q(igcm_ccn_mass) = rho_dust
473
474	radius(igcm_ccn_number) = radius(igcm_ccn_mass)
475	alpha_lift(igcm_ccn_number) = alpha_lift(igcm_ccn_mass)
476	alpha_devil(igcm_ccn_number) = alpha_devil(igcm_ccn_mass)
477	rho_q(igcm_ccn_number) = rho_q(igcm_ccn_mass)
478	endif ! of if (water.AND.doubleq.AND.scavenging)
479
480	c Submicron dust mode:
481	c --------------------
482
483	if (submicron) then
484	radius(igcm_dust_submicron)=0.1e-6
485	rho_q(igcm_dust_submicron)=rho_dust
486	if (doubleq) then
487	c If doubleq is also active, we use the population ratio:
488	alpha_lift(igcm_dust_submicron) =
489	& alpha_lift(igcm_dust_number)popratio
490	& rho_q(igcm_dust_submicron)4./3.pi*
491	& radius(igcm_dust_submicron)**3.
492	alpha_devil(igcm_dust_submicron)=1.e-30
493	else
494	alpha_lift(igcm_dust_submicron)=1e-6
495	alpha_devil(igcm_dust_submicron)=1.e-30
496	endif ! (doubleq)
497	end if ! (submicron)
498
499	c Initialization for photochemistry:
500	c ---------------------------------
501	if (photochem) then
502	! initialize chemistry+water (water will be correctly initialized below)
503	! by initializing everything which is not dust ...
504	do iq=1,nqmx
505	txt=noms(iq)
506	if (txt(1:4).ne."dust") then
507	radius(iq)=0.
508	alpha_lift(iq) =0.
509	alpha_devil(iq)=0.
510	endif
511	enddo ! do iq=1,nqmx
512	endif
513
514	c Initialization for water vapor
515	c ------------------------------
516	if(water) then
517	radius(igcm_h2o_vap)=0.
518	alpha_lift(igcm_h2o_vap) =0.
519	alpha_devil(igcm_h2o_vap)=0.
520	if(water.and.(nqmx.ge.2)) then
521	radius(igcm_h2o_ice)=3.e-6
522	rho_q(igcm_h2o_ice)=rho_ice
523	alpha_lift(igcm_h2o_ice) =0.
524	alpha_devil(igcm_h2o_ice)=0.
525	elseif(water.and.(nqmx.lt.2)) then
526	write(,) 'nqmx is too low : nqmx=', nqmx
527	write(,) 'water= ',water
528	endif
529
530	end if ! (water)
531
532	c Output for records:
533	c ~~~~~~~~~~~~~~~~~~
534	write(,)
535	Write(,) '******** initracer : dust transport parameters :'
536	write(,) 'alpha_lift = ', alpha_lift
537	write(,) 'alpha_devil = ', alpha_devil
538	write(,) 'radius = ', radius
539	if(doubleq) then
540	write(,) 'reff_lift (um) = ', reff_lift
541	write(,) 'size distribution variance = ', varian
542	write(,) 'r3n_q , ref_r0 : ', r3n_q , ref_r0
543	end if
544
545	!
546	! some extra (possibly redundant) sanity checks for tracers:
547	! ---------------------------------------------------------
548
549	if (doubleq) then
550	! verify that we indeed have dust_mass and dust_number tracers
551	if (igcm_dust_mass.eq.0) then
552	write(,) "initracer: error !!"
553	write(,) " cannot use doubleq option without ",
554	& "a dust_mass tracer !"
555	stop
556	endif
557	if (igcm_dust_number.eq.0) then
558	write(,) "initracer: error !!"
559	write(,) " cannot use doubleq option without ",
560	& "a dust_number tracer !"
561	stop
562	endif
563	endif
564
565	if ((.not.doubleq).and.(dustbin.gt.0)) then
566	! verify that we indeed have 'dustbin' dust tracers
567	count=0
568	do iq=1,dustbin
569	if (igcm_dustbin(iq).ne.0) then
570	count=count+1
571	endif
572	enddo
573	if (count.ne.dustbin) then
574	write(,) "initracer: error !!"
575	write(,) " dusbin is set to ",dustbin,
576	& " but we only have the following dust tracers:"
577	do iq=1,count
578	write(,)" ",trim(noms(igcm_dustbin(iq)))
579	enddo
580	stop
581	endif
582	endif
583
584	if (water) then
585	! verify that we indeed have h2o_vap and h2o_ice tracers
586	if (igcm_h2o_vap.eq.0) then
587	write(,) "initracer: error !!"
588	write(,) " cannot use water option without ",
589	& "an h2o_vap tracer !"
590	stop
591	endif
592	if (igcm_h2o_ice.eq.0) then
593	write(,) "initracer: error !!"
594	write(,) " cannot use water option without ",
595	& "an h2o_ice tracer !"
596	stop
597	endif
598	endif
599
600	if (scavenging) then
601	! verify that we indeed have ccn_mass and ccn_number tracers
602	if (igcm_ccn_mass.eq.0) then
603	write(,) "initracer: error !!"
604	write(,) " cannot use scavenging option without ",
605	& "a ccn_mass tracer !"
606	stop
607	endif
608	if (igcm_ccn_number.eq.0) then
609	write(,) "initracer: error !!"
610	write(,) " cannot use scavenging option without ",
611	& "a ccn_number tracer !"
612	stop
613	endif
614	endif ! of if (scavenging)
615
616	if (photochem .or. callthermos) then
617	! verify that we indeed have the chemistry tracers
618	if (igcm_co2.eq.0) then
619	write(,) "initracer: error !!"
620	write(,) " cannot use chemistry option without ",
621	& "a co2 tracer !"
622	stop
623	endif
624	if (igcm_co.eq.0) then
625	write(,) "initracer: error !!"
626	write(,) " cannot use chemistry option without ",
627	& "a co tracer !"
628	stop
629	endif
630	if (igcm_o.eq.0) then
631	write(,) "initracer: error !!"
632	write(,) " cannot use chemistry option without ",
633	& "a o tracer !"
634	stop
635	endif
636	if (igcm_o1d.eq.0) then
637	write(,) "initracer: error !!"
638	write(,) " cannot use chemistry option without ",
639	& "a o1d tracer !"
640	stop
641	endif
642	if (igcm_o2.eq.0) then
643	write(,) "initracer: error !!"
644	write(,) " cannot use chemistry option without ",
645	& "an o2 tracer !"
646	stop
647	endif
648	if (igcm_o3.eq.0) then
649	write(,) "initracer: error !!"
650	write(,) " cannot use chemistry option without ",
651	& "an o3 tracer !"
652	stop
653	endif
654	if (igcm_h.eq.0) then
655	write(,) "initracer: error !!"
656	write(,) " cannot use chemistry option without ",
657	& "a h tracer !"
658	stop
659	endif
660	if (igcm_h2.eq.0) then
661	write(,) "initracer: error !!"
662	write(,) " cannot use chemistry option without ",
663	& "a h2 tracer !"
664	stop
665	endif
666	if (igcm_oh.eq.0) then
667	write(,) "initracer: error !!"
668	write(,) " cannot use chemistry option without ",
669	& "an oh tracer !"
670	stop
671	endif
672	if (igcm_ho2.eq.0) then
673	write(,) "initracer: error !!"
674	write(,) " cannot use chemistry option without ",
675	& "a ho2 tracer !"
676	stop
677	endif
678	if (igcm_h2o2.eq.0) then
679	write(,) "initracer: error !!"
680	write(,) " cannot use chemistry option without ",
681	& "a h2o2 tracer !"
682	stop
683	endif
684	if (igcm_n2.eq.0) then
685	write(,) "initracer: error !!"
686	write(,) " cannot use chemistry option without ",
687	& "a n2 tracer !"
688	stop
689	endif
690	if (igcm_ar.eq.0) then
691	write(,) "initracer: error !!"
692	write(,) " cannot use chemistry option without ",
693	& "an ar tracer !"
694	stop
695	endif
696	endif ! of if (photochem .or. callthermos)
697
698	end

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

Download in other formats:

Original Format