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