Context Navigation

concentrations.F @ 1747

Last change on this file since 1747 was 1660, checked in by emillour, 8 years ago

Mars GCM:

Added possibility to run with an Helium "he" tracer (to be initialized at constant value of 3.6e-7 kg/kg_air, i.e. the 4ppm of Krasnopolsky 1996 EUVE satellite, using newstart).
corrected case for CH4 in aeronomars/photochemistry.F (was using undefined indexes when there is no CH4).
updated aki/cpi coefficients for Argon used to compute mean atmospheric Cp and thermal conductivity in aeronomars/concentrations.F

JYC+EM

File size: 9.1 KB

Line
1	SUBROUTINE concentrations(ngrid,nlayer,nq,
2	& pplay,pt,pdt,pq,pdq,ptimestep)
3
4	use tracer_mod, only: igcm_co2, igcm_co, igcm_o, igcm_o1d,
5	& igcm_o2, igcm_o3, igcm_h, igcm_h2,
6	& igcm_oh, igcm_ho2, igcm_n2, igcm_ar,
7	& igcm_h2o_vap, igcm_n, igcm_no, igcm_no2,
8	& igcm_n2d, igcm_co2plus, igcm_oplus,
9	& igcm_o2plus, igcm_coplus, igcm_cplus,
10	& igcm_nplus, igcm_noplus, igcm_n2plus,
11	& igcm_hplus, igcm_hco2plus, mmol,
12	& igcm_he
13	use conc_mod, only: mmean, Akknew, rnew, cpnew
14	implicit none
15
16	!=======================================================================
17	! CALCULATION OF MEAN MOLECULAR MASS, Cp, Akk and R
18	!
19	! mmean(ngrid,nlayer) amu
20	! cpnew(ngrid,nlayer) J/kg/K
21	! rnew(ngrid,nlayer) J/kg/K
22	! akknew(ngrid,nlayer) coefficient of thermal concduction
23	!
24	! version: April 2012 - Franck Lefevre
25	!=======================================================================
26
27	! declarations
28
29	#include "callkeys.h"
30	#include "chimiedata.h"
31
32	! input/output
33
34	integer,intent(in) :: ngrid ! number of atmospheric columns
35	integer,intent(in) :: nlayer ! number of atmospheric layers
36	integer,intent(in) :: nq ! number of tracers
37	real,intent(in) :: pplay(ngrid,nlayer)
38	real,intent(in) :: pt(ngrid,nlayer)
39	real,intent(in) :: pdt(ngrid,nlayer)
40	real,intent(in) :: pq(ngrid,nlayer,nq)
41	real,intent(in) :: pdq(ngrid,nlayer,nq)
42	real,intent(in) :: ptimestep
43
44	! local variables
45
46	integer :: i, l, ig, iq
47	integer, save :: nbq
48	integer,allocatable,save :: niq(:)
49	real :: ni(nq), ntot
50	real :: zq(ngrid, nlayer, nq)
51	real :: zt(ngrid, nlayer)
52	real,allocatable,save :: aki(:)
53	real,allocatable,save :: cpi(:)
54
55	logical, save :: firstcall = .true.
56
57	if (firstcall) then
58
59	! allocate local saved arrays:
60	allocate(aki(nq))
61	allocate(cpi(nq))
62	allocate(niq(nq))
63	! find index of chemical tracers to use
64	! initialize thermal conductivity and specific heat coefficients
65	! !? values are estimated
66
67	nbq = 0 ! to count number of tracers used in this subroutine
68
69	if (igcm_co2 /= 0) then
70	nbq = nbq + 1
71	niq(nbq) = igcm_co2
72	aki(nbq) = 3.072e-4
73	cpi(nbq) = 0.834e3
74	end if
75	if (igcm_co /= 0) then
76	nbq = nbq + 1
77	niq(nbq) = igcm_co
78	aki(nbq) = 4.87e-4
79	cpi(nbq) = 1.034e3
80	end if
81	if (igcm_o /= 0) then
82	nbq = nbq + 1
83	niq(nbq) = igcm_o
84	aki(nbq) = 7.59e-4
85	cpi(nbq) = 1.3e3
86	end if
87	if (igcm_o1d /= 0) then
88	nbq = nbq + 1
89	niq(nbq) = igcm_o1d
90	aki(nbq) = 7.59e-4 !?
91	cpi(nbq) = 1.3e3 !?
92	end if
93	if (igcm_o2 /= 0) then
94	nbq = nbq + 1
95	niq(nbq) = igcm_o2
96	aki(nbq) = 5.68e-4
97	cpi(nbq) = 0.9194e3
98	end if
99	if (igcm_o3 /= 0) then
100	nbq = nbq + 1
101	niq(nbq) = igcm_o3
102	aki(nbq) = 3.00e-4 !?
103	cpi(nbq) = 0.800e3 !?
104	end if
105	if (igcm_h /= 0) then
106	nbq = nbq + 1
107	niq(nbq) = igcm_h
108	aki(nbq) = 0.0
109	cpi(nbq) = 20.780e3
110	end if
111	if (igcm_h2 /= 0) then
112	nbq = nbq + 1
113	niq(nbq) = igcm_h2
114	aki(nbq) = 36.314e-4
115	cpi(nbq) = 14.266e3
116	end if
117	if (igcm_oh /= 0) then
118	nbq = nbq + 1
119	niq(nbq) = igcm_oh
120	aki(nbq) = 7.00e-4 !?
121	cpi(nbq) = 1.045e3
122	end if
123	if (igcm_ho2 /= 0) then
124	nbq = nbq + 1
125	niq(nbq) = igcm_ho2
126	aki(nbq) = 0.0
127	cpi(nbq) = 1.065e3 !?
128	end if
129	if (igcm_n2 /= 0) then
130	nbq = nbq + 1
131	niq(nbq) = igcm_n2
132	aki(nbq) = 5.6e-4
133	cpi(nbq) = 1.034e3
134	end if
135	if (igcm_ar /= 0) then
136	nbq = nbq + 1
137	niq(nbq) = igcm_ar
138	aki(nbq) = 3.4e-4
139	cpi(nbq) = 5.2e2
140	! aki/cpi from Leslie A. Guildner,
141	! JOURNAL OF RESEARCH of the National Bureau of Standards-
142	! A. Physics and Chemistry
143	! Vol. 79A, No.2, March-April 1975
144	end if
145	if (igcm_h2o_vap /= 0) then
146	nbq = nbq + 1
147	niq(nbq) = igcm_h2o_vap
148	aki(nbq) = 0.0
149	cpi(nbq) = 1.870e3
150	end if
151	if (igcm_n /= 0) then
152	nbq = nbq + 1
153	niq(nbq) = igcm_n
154	aki(nbq) = 0.0
155	cpi(nbq) = 0.0
156	endif
157	if(igcm_no /= 0) then
158	nbq = nbq + 1
159	niq(nbq) = igcm_no
160	aki(nbq) = 0.0
161	cpi(nbq) = 0.0
162	endif
163	if(igcm_no2 /= 0) then
164	nbq = nbq + 1
165	niq(nbq) = igcm_no2
166	aki(nbq) = 0.0
167	cpi(nbq) = 0.0
168	endif
169	if(igcm_n2d /= 0) then
170	nbq = nbq + 1
171	niq(nbq) = igcm_n2d
172	aki(nbq) = 0.0
173	cpi(nbq) = 0.0
174	endif
175	if(igcm_he /= 0) then
176	nbq = nbq + 1
177	niq(nbq) = igcm_he
178	aki(nbq) = 30.e-4
179	cpi(nbq) = 5.2e3
180	! aki/cpi from Leslie A. Guildner,
181	! JOURNAL OF RESEARCH of the National Bureau of Standards-
182	! A. Physics and Chemistry
183	! Vol. 79A, No.2, March-April 1975
184	endif
185	if(igcm_co2plus /= 0) then
186	nbq = nbq + 1
187	niq(nbq) = igcm_co2plus
188	aki(nbq) = 0.0
189	cpi(nbq) = 0.0
190	endif
191	if(igcm_oplus /= 0) then
192	nbq = nbq + 1
193	niq(nbq) = igcm_oplus
194	aki(nbq) = 0.0
195	cpi(nbq) = 0.0
196	endif
197	if(igcm_o2plus /= 0) then
198	nbq = nbq + 1
199	niq(nbq) = igcm_o2plus
200	aki(nbq) = 0.0
201	cpi(nbq) = 0.0
202	endif
203	if(igcm_coplus /= 0) then
204	nbq = nbq + 1
205	niq(nbq) = igcm_coplus
206	aki(nbq) = 0.0
207	cpi(nbq) = 0.0
208	endif
209	if(igcm_cplus /= 0) then
210	nbq = nbq + 1
211	niq(nbq) = igcm_cplus
212	aki(nbq) = 0.0
213	cpi(nbq) = 0.0
214	endif
215	if(igcm_nplus /= 0) then
216	nbq = nbq + 1
217	niq(nbq) = igcm_nplus
218	aki(nbq) = 0.0
219	cpi(nbq) = 0.0
220	endif
221	if(igcm_noplus /= 0) then
222	nbq = nbq + 1
223	niq(nbq) = igcm_noplus
224	aki(nbq) = 0.0
225	cpi(nbq) = 0.0
226	endif
227	if(igcm_n2plus /= 0) then
228	nbq = nbq + 1
229	niq(nbq) = igcm_n2plus
230	aki(nbq) = 0.0
231	cpi(nbq) = 0.0
232	endif
233	if(igcm_hplus /= 0) then
234	nbq = nbq + 1
235	niq(nbq) = igcm_hplus
236	aki(nbq) = 0.0
237	cpi(nbq) = 0.0
238	endif
239	if(igcm_hco2plus /= 0) then
240	nbq = nbq + 1
241	niq(nbq) = igcm_hco2plus
242	aki(nbq) = 0.0
243	cpi(nbq) = 0.0
244	endif
245
246	! tell the world about it:
247	write(,) "concentrations: firstcall, nbq=",nbq
248	write(,) " niq(1:nbq)=",niq(1:nbq)
249	write(,) " aki(1:nbq)=",aki(1:nbq)
250	write(,) " cpi(1:nbq)=",cpi(1:nbq)
251
252	firstcall = .false.
253
254	end if ! if (firstcall)
255
256	! update temperature
257
258	do l = 1,nlayer
259	do ig = 1,ngrid
260	zt(ig,l) = pt(ig,l) + pdt(ig,l)*ptimestep
261	end do
262	end do
263
264	! update tracers
265
266	do l = 1,nlayer
267	do ig = 1,ngrid
268	do i = 1,nbq
269	iq = niq(i)
270	zq(ig,l,iq) = max(1.e-30, pq(ig,l,iq)
271	$ + pdq(ig,l,iq)*ptimestep)
272	end do
273	end do
274	end do
275
276	! mmean : mean molecular mass
277	! rnew : specific gas constant
278
279	mmean(:,:) = 0.
280
281	do l = 1,nlayer
282	do ig = 1,ngrid
283	do i = 1,nbq
284	iq = niq(i)
285	mmean(ig,l) = mmean(ig,l) + zq(ig,l,iq)/mmol(iq)
286	end do
287	mmean(ig,l) = 1./mmean(ig,l)
288	rnew(ig,l) = 8.314/mmean(ig,l)*1.e3 ! J/kg/K
289	end do
290	end do
291
292	! cpnew : specicic heat
293	! akknew : thermal conductivity cofficient
294
295	cpnew(:,:) = 0.
296	akknew(:,:) = 0.
297
298	do l = 1,nlayer
299	do ig = 1,ngrid
300	ntot = pplay(ig,l)/(1.381e-23zt(ig,l))1.e-6 ! in #/cm3
301	do i = 1,nbq
302	iq = niq(i)
303	ni(iq) = ntotzq(ig,l,iq)mmean(ig,l)/mmol(iq)
304	cpnew(ig,l) = cpnew(ig,l) + ni(iq)*cpi(i)
305	akknew(ig,l) = akknew(ig,l) + ni(iq)*aki(i)
306	end do
307	cpnew(ig,l) = cpnew(ig,l)/ntot
308	akknew(ig,l) = akknew(ig,l)/ntot
309	end do
310	! print*, l, mmean(1,l), cpnew(1,l), rnew(1,l)
311	end do
312
313	return
314	end

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

Download in other formats:

Original Format