Context Navigation

cnuages3D.F @ 3094

Last change on this file since 3094 was 808, checked in by slebonnois, 12 years ago
SL: Many changes for VENUS (related to newstart) and TITAN (related to clouds). Please read DOC/chantiers/commit_importants.log (cf v808).
File size: 13.6 KB

Line
1	subroutine cnuages(
2	& tq,tqc1,tqc2,tqc3,tqcn,gaz1,gaz2,gaz3, ! aerosol/glace/gas
3	& ddt)
4
5	c
6	c
7	c °
8	c °
9	c SERT A APPELE LA ROUTINE MICROPHYSIQUE DES NUAGES
10	c
11	c °
12	c ICI ON NE FAIT QUE LA NUCLEATION/CONDENSATION
13	c ET GESTION°DES NOYAUX. LA SEDIMENTATION EST DANS
14	c SNUAGES.F
15	c ° °
16	c
17	c °
18	c °
19	c ° °
20	c
21	c ° °
22	c °
23	c ° °
24	c ° \\|/ °
25	c (@ @)°
26	c-----------------oOo--O--oOo--------------------------
27	c
28	c
29	c
30	c Interface entre physiq.F et les routines n_<nom_compose>.F
31	c
32	c Date: 3 Nov 2003
33	c
34	c EN ENTREE/SORTIE DE LA ROUTINE
35	c ------------------------------------
36	c
37	c Les aerosols, noyaux (tq,tqcn) sont en nbre/m^2 dans la colonne
38	c Les condensats (tqc1,tqc2) sont en volume/m^2 dans la colonne
39	c Le gaz (gaz1, gaz2) est en fraction molaire
40	c
41	c EN APPEL DES ROUTINES NUAGES
42	c ------------------------------------
43	c
44	c Les aerosols et noyaux doivent etre en nombre /kg d'air
45	c Les condensats doivent etre en volume / kg d'air
46	c Le gaz en kg/kg d'air
47	c
48	c LES TENDANCES ET DIFFERENCES SONT HOMOGENES AUX QUANTITES
49	c ------------------------------------------------------------
50	c
51	c
52	c------------------------------------------------------
53
54	use dimphy
55	IMPLICIT NONE
56	#include "dimensions.h"
57	#include "microtab.h"
58	#include "varmuphy.h"
59
60	integer NG1,NG,NL
61	parameter (NG1=1,NG=NG1,NL=llm)
62
63	c*************************************
64	c declaration des variables internes *
65	c*************************************
66
67	c INTERNE! *
68	c-----------------------*
69
70	real tqc1(NG,NL,nrad)
71	real tqc2(NG,NL,nrad)
72	real tqc3(NG,NL,nrad)
73	real tqcn(NG,NL,nrad)
74	real tq(NG,NL,nrad)
75	*
76	real tdqc1(NG,NL,nrad)
77	real tdqc2(NG,NL,nrad)
78	real tdqc3(NG,NL,nrad)
79	real tdq(NG,NL,nrad,ntype-2+1)
80	real tdqcn(NG,NL,nrad,ntype-2+1)
81	*
82	real gaz1(NG,NL)
83	real gaz2(NG,NL)
84	real gaz3(NG,NL)
85	real dgaz1(NG,NL)
86	real dgaz2(NG,NL)
87	real dgaz3(NG,NL)
88	*
89	real ppch4(NG,NL)
90	real ppc2h6(NG,NL)
91	real ppn2(NG,NL)
92	*
93	real pmixch4(NL)
94	real pmixc2h6(NL)
95	real pmixn2(NL)
96	c composition initiale estimée (interne)
97	real xprime1(NG,NL)
98	real xprime2(NG,NL)
99	real xprime3(NG,NL)
100	c composition calculée (output)
101	real x1(NL)
102	real x2(NL)
103	real x3(NL)
104	c moyenne "glissante" pondéréee (output + mémoire)
105	real x1o(NL)
106	real x2o(NL)
107	real x3o(NL)
108	real icefrac(NL)
109	real dmn2(NL+1)
110
111	real ppch4t,ppc2h6t,ppn2t
112	real psatch4,psatc2h6,psatn2
113	real xprime(3),x(3),frac
114	real melange
115	real sum,sum0
116
117	* RAPPEL: NG=1
118
119	real ddt
120	real masspaer
121
122	common/mixing/x1,x2,x3,icefrac,
123	& pmixch4,pmixc2h6,pmixn2,
124	& x1o,x2o,x3o
125
126
127
128	c FORMAT MICRO DES NUAGES
129	c------------------------*
130
131	real especes(NG,NL,3*nrad+1)
132	real condens(NG,NL,nrad)
133	real gg(NL),xmair
134	real effg ! effg est une fonction(z), z en m.
135
136	integer jsup,jinf,h,i,j,k,ndim
137	integer ival1,ival2,ival3
138
139	integer iprem
140
141	save iprem,xprime
142	data iprem/0/
143
144	ndim=3*nrad+1
145
146	do j=1,NL
147	gg(j)=effg(z(j))
148	enddo
149
150	*********************************************
151	*********************************************
152	* Appel de la condensation du methane
153	*********************************************
154	*********************************************
155
156
157
158	*=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
159	* Bilan avant sur le methane *
160	*=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
161
162	do i=1,ng1 ! ng1=1 !!
163	do j=1,NL
164
165	* RAZ des delta d'especes communes.
166	*----------------------------------
167	do k=1,nrad
168	tdqcn(i,j,k,1)= 0.
169	tdq( i,j,k,1) = 0.
170	enddo
171
172	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
173
174	do k=1,nrad
175	especes(i,j,k)=tq(i,j,k) /xmair ! aerosols, noyaux,
176	especes(i,j,k+nrad)=tqc1(i,j,k) /xmair ! methane condense,
177	condens(i,j,k)=(tqc2(i,j,k)+tqc3(i,j,k))/xmair! autre(s) condensat(s)
178	especes(i,j,k+2*nrad)=tqcn(i,j,k)/xmair ! nombre/kg & volume/kg
179	enddo
180	especes(i,j,3nrad+1)=gaz1(i,j)mch4/mair ! methane gazeux kg/kg
181	enddo
182	enddo
183
184
185	1001 format(7(1x,e12.6),' avN2CH4C2H6')
186	1003 format(7(1x,e12.6),' miN2CH4C2H6')
187	1002 format(7(1x,e12.6),' apN2CH4C2H6')
188	*
189	call n_methane(ng1,ndim,nrad,ddt,
190	& p,t,r_e,especes,condens)
191
192
193	*=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
194	* Bilan apres sur le methane *
195	*=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
196
197	do i=1,ng1
198	do j=1,NL
199
200	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
201
202	* ici ce sont les tendances a sortir de nuages.F pour le methane....
203	*-------------------------------------------------------------------
204	sum=0.
205	do k=1,nrad
206	tdqc1(i,j,k)=(especes(i,j,k+nrad)*xmair-tqc1(i,j,k) )
207	sum =sum+tdqc1(i,j,k)/xmair*rhoi_ch4
208	tqc1(i,j,k) = especes(i,j,k+nrad)*xmair
209	enddo
210
211	dgaz1(i,j)= especes(i,j,3nrad+1)mair/mch4-gaz1(i,j)
212	gaz1(i,j)= especes(i,j,3nrad+1)mair/mch4
213
214	c dgaz1(i,j)=-sum*xmuair/16.
215	c gaz1(i,j)=gaz1(i,j)+dgaz1(i,j)
216
217	* Premiere tendance sur les variables communes (aerosols et noyaux)
218	*------------------------------------------------------------------
219	do k=1,nrad
220	tdqcn(i,j,k,1)=(especes(i,j,k+2nrad)xmair-tqcn(i,j,k))
221	tdq( i,j,k,1) =(especes(i,j,k)*xmair -tq(i,j,k))
222	enddo
223	enddo
224	enddo
225
226
227
228	* attention, si il y a de l'ethane sur les noyaux... il est impossible
229	* de les restituer - en revanche on peut en creer de nouveaux ! !!
230	* Le corrolaire de la condition ci dessus est que si il est impossible de
231	* restituer des noyaux, le nombre d'aeorsols ne peux pas augmenter, il
232	* peut en revanche diminuer
233
234
235	*********************************************
236	*********************************************
237	* Appel de la condensation de l'ethane
238	*********************************************
239	*********************************************
240
241	do i=1,ng1
242	do j=1,NL
243	do k=1,nrad
244	tdqcn(i,j,k,2)= 0.
245	tdq( i,j,k,2) = 0.
246	enddo
247
248	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
249
250	do k=1,nrad
251	especes(i,j,k)=tq(i,j,k)/xmair
252	especes(i,j,k+nrad)=tqc2(i,j,k)/xmair ! ethane condense
253	condens(i,j,k)=(tqc1(i,j,k)+tqc3(i,j,k))/xmair ! autres condensats
254	especes(i,j,k+2*nrad)=tqcn(i,j,k)/xmair
255	enddo
256
257	especes(i,j,3nrad+1)=gaz2(i,j)mc2h6/mair ! ethane gazeux
258
259	enddo
260	enddo
261
262
263
264	call n_ethane(ng1,ndim,nrad,ddt,
265	& p,t,r_e,especes,condens)
266
267	do i=1,ng1
268	do j=1,NL
269	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
270
271	* ici ce sont les tendances a sortir de nuages.F pour l'ethane....
272	*-----------------------------------------------------------------
273
274	sum=0.
275	do k=1,nrad
276	tdqc2(i,j,k)=(especes(i,j,k+nrad)*xmair-tqc2(i,j,k) )
277	sum =sum+tdqc2(i,j,k)/xmair*rhoi_c2h6
278	tqc2(i,j,k) = especes(i,j,k+nrad)*xmair
279	enddo
280
281	dgaz2(i,j)=especes(i,j,3nrad+1)mair/mc2h6 - gaz2(i,j)
282	gaz2(i,j)=especes(i,j,3nrad+1)mair/mc2h6
283
284	c dgaz2(i,j)=-sum*xmuair/30.
285	c gaz2(i,j)=gaz2(i,j)+dgaz2(i,j)
286
287
288	* Deuxieme tendance sur les variables communes (aerosols et noyaux)
289	*------------------------------------------------------------------
290
291	do k=1,nrad
292	tdqcn(i,j,k,2)=(especes(i,j,k+2nrad)xmair-tqcn(i,j,k))
293	tdq(i,j,k,2) =(especes(i,j,k)*xmair -tq(i,j,k))
294	enddo
295
296	enddo
297	enddo
298
299
300	*********************************************
301	*********************************************
302	* Appel de la condensation de l'acethylene
303	*********************************************
304	*********************************************
305
306	do i=1,ng1
307	do j=1,NL
308	do k=1,nrad
309	tdqcn(i,j,k,3)= 0.
310	tdq( i,j,k,3) = 0.
311	enddo
312	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
313
314	do k=1,nrad
315	especes(i,j,k)=tq(i,j,k)/xmair
316	especes(i,j,k+nrad)=tqc3(i,j,k)/xmair ! acethylene condense
317	condens(i,j,k)=(tqc1(i,j,k)+tqc2(i,j,k))/xmair ! autres condensats
318	especes(i,j,k+2*nrad)=tqcn(i,j,k)/xmair
319	enddo
320
321	especes(i,j,3nrad+1)=gaz3(i,j)mc2h2/mair ! acethylene gazeux
322
323	enddo
324	enddo
325
326	call n_acethylene(ng1,ndim,nrad,ddt,
327	& p,t,r_e,especes,condens)
328
329
330	do i=1,ng1
331	do j=1,NL
332
333	xmair=(pb(j+1)-pb(j))/gg(j)/dzb(j)
334
335	* ici ce sont les tendances a sortir de nuages.F pour l'ethane....
336	*-----------------------------------------------------------------
337
338	sum=0.
339	do k=1,nrad
340	tdqc3(i,j,k)=(especes(i,j,k+nrad)*xmair-tqc3(i,j,k) )
341	sum =sum+tdqc3(i,j,k)/xmair*rhoi_c2h2
342	tqc3(i,j,k) = especes(i,j,k+nrad)*xmair
343	enddo
344
345	dgaz3(i,j)=especes(i,j,3nrad+1)mair/mc2h2 - gaz3(i,j)
346	gaz3(i,j)=especes(i,j,3nrad+1)mair/mc2h2
347
348	c dgaz3(i,j)=-sum*xmuair/26.
349	c gaz3(i,j)=gaz3(i,j)+dgaz3(i,j)
350
351	* Troisieme tendance sur les variables communes (aerosols et noyaux)
352	*------------------------------------------------------------------
353
354
355	do k=1,nrad
356	tdqcn(i,j,k,3)=(especes(i,j,k+2nrad)xmair-tqcn(i,j,k))
357	tdq(i,j,k,3) =(especes(i,j,k)*xmair -tq(i,j,k))
358	enddo
359
360	enddo
361	enddo
362
363
364	* FIN DES APPELS DE NUAGES ET BILAN DES TENDANCES...
365	*------------------------------------------------------------------
366
367
368	do i=1,ng1
369	do j=1,NL
370	do k=1,nrad
371
372	* Ici on test l'activité nuageuse : si on a l'association
373	* tdqcX(i,j,k) = 0 et tqcX(i,j,k) = 0 alors ivalX reste à 0 (pas d'ativité)
374	* sinon ivalX passe à 1 (activité)
375	*------------------------------------------------------------------------------------
376
377	ival1=0
378	ival2=0
379	ival3=0
380
381	if(tdqc1(i,j,k).ne.0. .or. tqc1(i,j,k).gt.0.) ival1=1
382	if(tdqc2(i,j,k).ne.0. .or. tqc2(i,j,k).gt.0.) ival2=1
383	if(tdqc3(i,j,k).ne.0. .or. tqc3(i,j,k).gt.0.) ival3=1
384
385	* Ici on definit la tendances des noyaux en faisant deux choses:
386	* -1 On ecarte les cas tdqcn(i,j,k,X)=0 si ils sont associés à une
387	* absence d'activité nuageuse de l'espèce (tdqcX(i,j,k)=0.)
388	* -2 Sélectionne la tendance la plus élevée. Si aucune activité nuageuse
389	* n'exits dans cette case (ivalX=0 pour les 3 especes), alors on
390	* retrouve la valeur -1.e40 que l'on mets alors à 0.
391	*----------------------------------------------------------------------
392	c23456789012345678901234567890123456789012345678901234567890123456789012
393
394	tdqcn(i,j,k,ntype-1)=-1.e40 ! plus petite valeur possible
395
396	if(ival1.eq.1.and.
397	& tdqcn(i,j,k,1).ge.tdqcn(i,j,k,ntype-1)) ! Si activité de l'espece 1
398	& tdqcn(i,j,k,ntype-1)=tdqcn(i,j,k,1)
399	if(ival2.eq.1.and.
400	& tdqcn(i,j,k,2).ge.tdqcn(i,j,k,ntype-1)) ! Si activité de l'espece 2
401	& tdqcn(i,j,k,ntype-1)=tdqcn(i,j,k,2)
402	if(ival3.eq.1.and.
403	& tdqcn(i,j,k,3).ge.tdqcn(i,j,k,ntype-1)) ! Si activité de l'espece 3
404	& tdqcn(i,j,k,ntype-1)=tdqcn(i,j,k,3)
405
406	if(tdqcn(i,j,k,ntype-1).le.-0.99e39)
407	& tdqcn(i,j,k,ntype-1)=0.
408
409	tdq(i,j,k,ntype-1)=1.e40 ! plus grande valeur possible
410
411	if(ival1.eq.1 .and. tdq(i,j,k,1).le.tdq(i,j,k,ntype-1)) ! Si activité de l'espece 1
412	& tdq(i,j,k,ntype-1)=tdq(i,j,k,1)
413	if(ival2.eq.1 .and. tdq(i,j,k,2).le.tdq(i,j,k,ntype-1)) ! Si activité de l'espece 2
414	& tdq(i,j,k,ntype-1)=tdq(i,j,k,2)
415	if(ival3.eq.1 .and. tdq(i,j,k,3).le.tdq(i,j,k,ntype-1)) ! Si activité de l'espece 3
416	& tdq(i,j,k,ntype-1)=tdq(i,j,k,3)
417
418	if(tdq(i,j,k,ntype-1).ge.0.99e39) tdq(i,j,k,ntype-1)=0.
419
420	tqcn(i,j,k)=tqcn(i,j,k)+tdqcn(i,j,k,ntype-1) ! Alors on ajoute les tendances (positive pour qcn ?)
421	tq(i,j,k) =tq(i,j,k) +tdq(i,j,k,ntype-1) !
422
423	if(tqcn(i,j,k).le.0.) tqcn(i,j,k)=0. ! et on régularise les tableaux noyaux et aerosols.
424	if(tq(i,j,k) .le.0.) tq(i,j,k)=0. !
425
426	enddo
427	enddo
428	enddo
429
430	continue
431
432	1202 format(i2,1x,i2,6(1x,e12.4) )
433	return
434	end

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

Download in other formats:

Original Format