Context Navigation

fisrtilp.F @ 1239

Last change on this file since 1239 was 1091, checked in by yann meurdesoif, 16 years ago
optimisation Othman Bouizi : fisrtilp
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 16.8 KB

Line
1	!
2	! $Header$
3	!
4	c
5	SUBROUTINE fisrtilp(dtime,paprs,pplay,t,q,ptconv,ratqs,
6	s d_t, d_q, d_ql, rneb, radliq, rain, snow,
7	s pfrac_impa, pfrac_nucl, pfrac_1nucl,
8	s frac_impa, frac_nucl,
9	s prfl, psfl, rhcl)
10
11	c
12	USE dimphy
13	IMPLICIT none
14	c======================================================================
15	c Auteur(s): Z.X. Li (LMD/CNRS)
16	c Date: le 20 mars 1995
17	c Objet: condensation et precipitation stratiforme.
18	c schema de nuage
19	c======================================================================
20	c======================================================================
21	cym#include "dimensions.h"
22	cym#include "dimphy.h"
23	#include "YOMCST.h"
24	#include "tracstoke.h"
25	#include "fisrtilp.h"
26	c
27	c Arguments:
28	c
29	REAL dtime ! intervalle du temps (s)
30	REAL paprs(klon,klev+1) ! pression a inter-couche
31	REAL pplay(klon,klev) ! pression au milieu de couche
32	REAL t(klon,klev) ! temperature (K)
33	REAL q(klon,klev) ! humidite specifique (kg/kg)
34	REAL d_t(klon,klev) ! incrementation de la temperature (K)
35	REAL d_q(klon,klev) ! incrementation de la vapeur d'eau
36	REAL d_ql(klon,klev) ! incrementation de l'eau liquide
37	REAL rneb(klon,klev) ! fraction nuageuse
38	REAL radliq(klon,klev) ! eau liquide utilisee dans rayonnements
39	REAL rhcl(klon,klev) ! humidite relative en ciel clair
40	REAL rain(klon) ! pluies (mm/s)
41	REAL snow(klon) ! neige (mm/s)
42	REAL prfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
43	REAL psfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s)
44	cAA
45	c Coeffients de fraction lessivee : pour OFF-LINE
46	c
47	REAL pfrac_nucl(klon,klev)
48	REAL pfrac_1nucl(klon,klev)
49	REAL pfrac_impa(klon,klev)
50	c
51	c Fraction d'aerosols lessivee par impaction et par nucleation
52	c POur ON-LINE
53	c
54	REAL frac_impa(klon,klev)
55	REAL frac_nucl(klon,klev)
56	real zct(klon),zcl(klon)
57	cAA
58	c
59	c Options du programme:
60	c
61	REAL seuil_neb ! un nuage existe vraiment au-dela
62	PARAMETER (seuil_neb=0.001)
63
64	INTEGER ninter ! sous-intervals pour la precipitation
65	INTEGER ncoreczq
66	PARAMETER (ninter=5)
67	LOGICAL evap_prec ! evaporation de la pluie
68	PARAMETER (evap_prec=.TRUE.)
69	REAL ratqs(klon,klev) ! determine la largeur de distribution de vapeur
70	logical ptconv(klon,klev) ! determine la largeur de distribution de vapeur
71
72	real zpdf_sig(klon),zpdf_k(klon),zpdf_delta(klon)
73	real Zpdf_a(klon),zpdf_b(klon),zpdf_e1(klon),zpdf_e2(klon)
74	real erf
75	c
76	LOGICAL cpartiel ! condensation partielle
77	PARAMETER (cpartiel=.TRUE.)
78	REAL t_coup
79	PARAMETER (t_coup=234.0)
80	c
81	c Variables locales:
82	c
83	INTEGER i, k, n, kk
84	REAL zqs(klon), zdqs(klon), zdelta, zcor, zcvm5
85	REAL zrfl(klon), zrfln(klon), zqev, zqevt
86	REAL zoliq(klon), zcond(klon), zq(klon), zqn(klon), zdelq
87	REAL ztglace, zt(klon)
88	INTEGER nexpo ! exponentiel pour glace/eau
89	REAL zdz(klon),zrho(klon),ztot(klon), zrhol(klon)
90	REAL zchau(klon),zfroi(klon),zfice(klon),zneb(klon)
91	c
92	LOGICAL appel1er
93	SAVE appel1er
94	c$OMP THREADPRIVATE(appel1er)
95	c
96	c---------------------------------------------------------------
97	c
98	cAA Variables traceurs:
99	cAA Provisoire !!! Parametres alpha du lessivage
100	cAA A priori on a 4 scavenging # possibles
101	c
102	REAL a_tr_sca(4)
103	save a_tr_sca
104	c$OMP THREADPRIVATE(a_tr_sca)
105	c
106	c Variables intermediaires
107	c
108	REAL zalpha_tr
109	REAL zfrac_lessi
110	REAL zprec_cond(klon)
111	cAA
112	REAL zmair, zcpair, zcpeau
113	C Pour la conversion eau-neige
114	REAL zlh_solid(klon), zm_solid
115	cIM
116	cym INTEGER klevm1
117	c---------------------------------------------------------------
118	c
119	c Fonctions en ligne:
120	c
121	REAL fallvs,fallvc ! vitesse de chute pour crystaux de glace
122	REAL zzz
123	#include "YOETHF.h"
124	#include "FCTTRE.h"
125	fallvc (zzz) = 3.29/2.0 * ((zzz)*0.16) ffallv_con
126	fallvs (zzz) = 3.29/2.0 * ((zzz)*0.16) ffallv_lsc
127	c
128	DATA appel1er /.TRUE./
129	cym
130	zdelq=0.0
131
132	IF (appel1er) THEN
133	c
134	PRINT*, 'fisrtilp, ninter:', ninter
135	PRINT*, 'fisrtilp, evap_prec:', evap_prec
136	PRINT*, 'fisrtilp, cpartiel:', cpartiel
137	IF (ABS(dtime/FLOAT(ninter)-360.0).GT.0.001) THEN
138	PRINT*, 'fisrtilp: Ce n est pas prevu, voir Z.X.Li', dtime
139	PRINT*, 'Je prefere un sous-intervalle de 6 minutes'
140	c CALL abort
141	ENDIF
142	appel1er = .FALSE.
143	c
144	cAA initialiation provisoire
145	a_tr_sca(1) = -0.5
146	a_tr_sca(2) = -0.5
147	a_tr_sca(3) = -0.5
148	a_tr_sca(4) = -0.5
149	c
150	cAA Initialisation a 1 des coefs des fractions lessivees
151	c
152	DO k = 1, klev
153	DO i = 1, klon
154	pfrac_nucl(i,k)=1.
155	pfrac_1nucl(i,k)=1.
156	pfrac_impa(i,k)=1.
157	ENDDO
158	ENDDO
159
160	ENDIF ! test sur appel1er
161	c
162	cMAf Initialisation a 0 de zoliq
163	DO i = 1, klon
164	zoliq(i)=0.
165	ENDDO
166	c Determiner les nuages froids par leur temperature
167	c nexpo regle la raideur de la transition eau liquide / eau glace.
168	c
169	ztglace = RTT - 15.0
170	nexpo = 6
171	ccc nexpo = 1
172	c
173	c Initialiser les sorties:
174	c
175	DO k = 1, klev+1
176	DO i = 1, klon
177	prfl(i,k) = 0.0
178	psfl(i,k) = 0.0
179	ENDDO
180	ENDDO
181
182	DO k = 1, klev
183	DO i = 1, klon
184	d_t(i,k) = 0.0
185	d_q(i,k) = 0.0
186	d_ql(i,k) = 0.0
187	rneb(i,k) = 0.0
188	radliq(i,k) = 0.0
189	frac_nucl(i,k) = 1.
190	frac_impa(i,k) = 1.
191	ENDDO
192	ENDDO
193	DO i = 1, klon
194	rain(i) = 0.0
195	snow(i) = 0.0
196	ENDDO
197	c
198	c Initialiser le flux de precipitation a zero
199	c
200	DO i = 1, klon
201	zrfl(i) = 0.0
202	zneb(i) = seuil_neb
203	ENDDO
204	c
205	c
206	cAA Pour plus de securite
207
208	zalpha_tr = 0.
209	zfrac_lessi = 0.
210
211	cAA----------------------------------------------------------
212	c
213	ncoreczq=0
214	c Boucle verticale (du haut vers le bas)
215	c
216	cIM : klevm1
217	cym klevm1=klev-1
218	DO 9999 k = klev, 1, -1
219	c
220	cAA----------------------------------------------------------
221	c
222	DO i = 1, klon
223	zt(i)=t(i,k)
224	zq(i)=q(i,k)
225	ENDDO
226	c
227	c Calculer la varition de temp. de l'air du a la chaleur sensible
228	C transporter par la pluie.
229	C Il resterait a rajouter cet effet de la chaleur sensible sur les
230	C flux de surface, du a la diff. de temp. entre le 1er niveau et la
231	C surface.
232	C
233	IF(k.LE.klevm1) THEN
234	DO i = 1, klon
235	cIM
236	zmair=(paprs(i,k)-paprs(i,k+1))/RG
237	zcpair=RCPD(1.0+RVTMP2zq(i))
238	zcpeau=RCPD*RVTMP2
239	zt(i) = ( (t(i,k+1)+d_t(i,k+1))zrfl(i)dtime*zcpeau
240	$ + zmairzcpairzt(i) )
241	$ / (zmairzcpair + zrfl(i)dtime*zcpeau)
242	C C WRITE (6,*) 'cppluie ', zt(i)-(t(i,k+1)+d_t(i,k+1))
243	ENDDO
244	ENDIF
245	c
246	c
247	c Calculer l'evaporation de la precipitation
248	c
249
250
251	IF (evap_prec) THEN
252	DO i = 1, klon
253	IF (zrfl(i) .GT.0.) THEN
254	IF (thermcep) THEN
255	zdelta=MAX(0.,SIGN(1.,RTT-zt(i)))
256	zqs(i)= R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
257	zqs(i)=MIN(0.5,zqs(i))
258	zcor=1./(1.-RETV*zqs(i))
259	zqs(i)=zqs(i)*zcor
260	ELSE
261	IF (zt(i) .LT. t_coup) THEN
262	zqs(i) = qsats(zt(i)) / pplay(i,k)
263	ELSE
264	zqs(i) = qsatl(zt(i)) / pplay(i,k)
265	ENDIF
266	ENDIF
267	zqev = MAX (0.0, (zqs(i)-zq(i))*zneb(i) )
268	zqevt = coef_eva * (1.0-zq(i)/zqs(i)) * SQRT(zrfl(i))
269	. * (paprs(i,k)-paprs(i,k+1))/pplay(i,k)zt(i)RD/RG
270	zqevt = MAX(0.0,MIN(zqevt,zrfl(i)))
271	. * RG*dtime/(paprs(i,k)-paprs(i,k+1))
272	zqev = MIN (zqev, zqevt)
273	zrfln(i) = zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1))
274	. /RG/dtime
275
276	c pour la glace, on r�vapore toute la pr�ip dans la couche du dessous
277	c la glace venant de la couche du dessus est simplement dans la couche
278	c du dessous.
279
280	IF (zt(i) .LT. t_coup.and.reevap_ice) zrfln(i)=0.
281
282	zq(i) = zq(i) - (zrfln(i)-zrfl(i))
283	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
284	zt(i) = zt(i) + (zrfln(i)-zrfl(i))
285	. * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime
286	. * RLVTT/RCPD/(1.0+RVTMP2*zq(i))
287	zrfl(i) = zrfln(i)
288	ENDIF
289	ENDDO
290	ENDIF
291	c
292	c Calculer Qs et L/Cp*dQs/dT:
293	c
294	IF (thermcep) THEN
295	DO i = 1, klon
296	zdelta = MAX(0.,SIGN(1.,RTT-zt(i)))
297	zcvm5 = R5LESRLVTT(1.-zdelta) + R5IESRLSTTzdelta
298	zcvm5 = zcvm5 /RCPD/(1.0+RVTMP2*zq(i))
299	zqs(i) = R2ES*FOEEW(zt(i),zdelta)/pplay(i,k)
300	zqs(i) = MIN(0.5,zqs(i))
301	zcor = 1./(1.-RETV*zqs(i))
302	zqs(i) = zqs(i)*zcor
303	zdqs(i) = FOEDE(zt(i),zdelta,zcvm5,zqs(i),zcor)
304	ENDDO
305	ELSE
306	DO i = 1, klon
307	IF (zt(i).LT.t_coup) THEN
308	zqs(i) = qsats(zt(i))/pplay(i,k)
309	zdqs(i) = dqsats(zt(i),zqs(i))
310	ELSE
311	zqs(i) = qsatl(zt(i))/pplay(i,k)
312	zdqs(i) = dqsatl(zt(i),zqs(i))
313	ENDIF
314	ENDDO
315	ENDIF
316	c
317	c Determiner la condensation partielle et calculer la quantite
318	c de l'eau condensee:
319	c
320	IF (cpartiel) THEN
321
322	c print*,'Dans partiel k=',k
323	c
324	c Calcul de l'eau condensee et de la fraction nuageuse et de l'eau
325	c nuageuse a partir des PDF de Sandrine Bony.
326	c rneb : fraction nuageuse
327	c zqn : eau totale dans le nuage
328	c zcond : eau condensee moyenne dans la maille.
329	c on prend en compte le r�hauffement qui diminue la partie condensee
330	c
331	c Version avec les raqts
332
333	if (iflag_pdf.eq.0) then
334
335	do i=1,klon
336	zdelq = min(ratqs(i,k),0.99) * zq(i)
337	rneb(i,k) = (zq(i)+zdelq-zqs(i)) / (2.0*zdelq)
338	zqn(i) = (zq(i)+zdelq+zqs(i))/2.0
339	enddo
340
341	else
342	c
343	c Version avec les nouvelles PDFs.
344	do i=1,klon
345	if(zq(i).lt.1.e-15) then
346	ncoreczq=ncoreczq+1
347	zq(i)=1.e-15
348	endif
349	enddo
350	do i=1,klon
351	zpdf_sig(i)=ratqs(i,k)*zq(i)
352	zpdf_k(i)=-sqrt(log(1.+(zpdf_sig(i)/zq(i))**2))
353	zpdf_delta(i)=log(zq(i)/zqs(i))
354	zpdf_a(i)=zpdf_delta(i)/(zpdf_k(i)*sqrt(2.))
355	zpdf_b(i)=zpdf_k(i)/(2.*sqrt(2.))
356	zpdf_e1(i)=zpdf_a(i)-zpdf_b(i)
357	zpdf_e1(i)=sign(min(abs(zpdf_e1(i)),5.),zpdf_e1(i))
358	zpdf_e1(i)=1.-erf(zpdf_e1(i))
359	zpdf_e2(i)=zpdf_a(i)+zpdf_b(i)
360	zpdf_e2(i)=sign(min(abs(zpdf_e2(i)),5.),zpdf_e2(i))
361	zpdf_e2(i)=1.-erf(zpdf_e2(i))
362	if (zpdf_e1(i).lt.1.e-10) then
363	rneb(i,k)=0.
364	zqn(i)=zqs(i)
365	else
366	rneb(i,k)=0.5*zpdf_e1(i)
367	zqn(i)=zq(i)*zpdf_e2(i)/zpdf_e1(i)
368	endif
369
370	enddo
371
372	endif ! iflag_pdf
373
374	DO i=1,klon
375	IF (rneb(i,k) .LE. 0.0) THEN
376	zqn(i) = 0.0
377	rneb(i,k) = 0.0
378	zcond(i) = 0.0
379	rhcl(i,k)=zq(i)/zqs(i)
380	ELSE IF (rneb(i,k) .GE. 1.0) THEN
381	zqn(i) = zq(i)
382	rneb(i,k) = 1.0
383	zcond(i) = MAX(0.0,zqn(i)-zqs(i))
384	rhcl(i,k)=1.0
385	ELSE
386	zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)
387	rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
388	ENDIF
389	ENDDO
390	! do i=1,klon
391	! IF (rneb(i,k) .LE. 0.0) zqn(i) = 0.0
392	! IF (rneb(i,k) .GE. 1.0) zqn(i) = zq(i)
393	! rneb(i,k) = MAX(0.0,MIN(1.0,rneb(i,k)))
394	!c zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1.+zdqs(i))
395	!c On ne divise pas par 1+zdqs pour forcer a avoir l'eau predite par
396	!c la convection.
397	!c ATTENTION !!! Il va falloir verifier tout ca.
398	! zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)
399	!c print*,'ZDQS ',zdqs(i)
400	!c--Olivier
401	! rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i)
402	! IF (rneb(i,k) .LE. 0.0) rhcl(i,k)=zq(i)/zqs(i)
403	! IF (rneb(i,k) .GE. 1.0) rhcl(i,k)=1.0
404	!c--fin
405	! ENDDO
406	ELSE
407	DO i = 1, klon
408	IF (zq(i).GT.zqs(i)) THEN
409	rneb(i,k) = 1.0
410	ELSE
411	rneb(i,k) = 0.0
412	ENDIF
413	zcond(i) = MAX(0.0,zq(i)-zqs(i))/(1.+zdqs(i))
414	ENDDO
415	ENDIF
416	c
417	DO i = 1, klon
418	zq(i) = zq(i) - zcond(i)
419	c zt(i) = zt(i) + zcond(i) * RLVTT/RCPD
420	zt(i) = zt(i) + zcond(i) * RLVTT/RCPD/(1.0+RVTMP2*zq(i))
421	ENDDO
422	c
423	c Partager l'eau condensee en precipitation et eau liquide nuageuse
424	c
425	DO i = 1, klon
426	IF (rneb(i,k).GT.0.0) THEN
427	zoliq(i) = zcond(i)
428	zrho(i) = pplay(i,k) / zt(i) / RD
429	zdz(i) = (paprs(i,k)-paprs(i,k+1)) / (zrho(i)*RG)
430	zfice(i) = 1.0 - (zt(i)-ztglace) / (273.13-ztglace)
431	zfice(i) = MIN(MAX(zfice(i),0.0),1.0)
432	zfice(i) = zfice(i)**nexpo
433	zneb(i) = MAX(rneb(i,k), seuil_neb)
434	radliq(i,k) = zoliq(i)/FLOAT(ninter+1)
435	ENDIF
436	ENDDO
437	c
438	DO n = 1, ninter
439	DO i = 1, klon
440	IF (rneb(i,k).GT.0.0) THEN
441	zrhol(i) = zrho(i) * zoliq(i) / zneb(i)
442
443	if (ptconv(i,k)) then
444	zcl(i)=cld_lc_con
445	zct(i)=1./cld_tau_con
446	else
447	zcl(i)=cld_lc_lsc
448	zct(i)=1./cld_tau_lsc
449	endif
450	c quantit�d'eau ��minier.
451	zchau(i) = zct(i)dtime/FLOAT(ninter) zoliq(i)
452	. (1.0-EXP(-(zoliq(i)/zneb(i)/zcl(i))2)) (1.-zfice(i))
453	c meme chose pour la glace.
454	if (ptconv(i,k)) then
455	zfroi(i) = dtime/FLOAT(ninter)/zdz(i)*zoliq(i)
456	. fallvc(zrhol(i)) zfice(i)
457	else
458	zfroi(i) = dtime/FLOAT(ninter)/zdz(i)*zoliq(i)
459	. fallvs(zrhol(i)) zfice(i)
460	endif
461	ztot(i) = zchau(i) + zfroi(i)
462	IF (zneb(i).EQ.seuil_neb) ztot(i) = 0.0
463	ztot(i) = MIN(MAX(ztot(i),0.0),zoliq(i))
464	zoliq(i) = MAX(zoliq(i)-ztot(i), 0.0)
465	radliq(i,k) = radliq(i,k) + zoliq(i)/FLOAT(ninter+1)
466	ENDIF
467	ENDDO
468	ENDDO
469	c
470	DO i = 1, klon
471	IF (rneb(i,k).GT.0.0) THEN
472	d_ql(i,k) = zoliq(i)
473	zrfl(i) = zrfl(i)+ MAX(zcond(i)-zoliq(i),0.0)
474	. * (paprs(i,k)-paprs(i,k+1))/(RG*dtime)
475	ENDIF
476	IF (zt(i).LT.RTT) THEN
477	psfl(i,k)=zrfl(i)
478	ELSE
479	prfl(i,k)=zrfl(i)
480	ENDIF
481	ENDDO
482	c
483	c Calculer les tendances de q et de t:
484	c
485	DO i = 1, klon
486	d_q(i,k) = zq(i) - q(i,k)
487	d_t(i,k) = zt(i) - t(i,k)
488	ENDDO
489	c
490	cAA--------------- Calcul du lessivage stratiforme -------------
491
492	DO i = 1,klon
493	c
494	zprec_cond(i) = MAX(zcond(i)-zoliq(i),0.0)
495	. * (paprs(i,k)-paprs(i,k+1))/RG
496	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
497	cAA lessivage nucleation LMD5 dans la couche elle-meme
498	if (t(i,k) .GE. ztglace) THEN
499	zalpha_tr = a_tr_sca(3)
500	else
501	zalpha_tr = a_tr_sca(4)
502	endif
503	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
504	pfrac_nucl(i,k)=pfrac_nucl(i,k)(1.-zneb(i)zfrac_lessi)
505	frac_nucl(i,k)= 1.-zneb(i)*zfrac_lessi
506	c
507	c nucleation avec un facteur -1 au lieu de -0.5
508	zfrac_lessi = 1. - EXP(-zprec_cond(i)/zneb(i))
509	pfrac_1nucl(i,k)=pfrac_1nucl(i,k)(1.-zneb(i)zfrac_lessi)
510	ENDIF
511	c
512	ENDDO ! boucle sur i
513	c
514	cAA Lessivage par impaction dans les couches en-dessous
515	DO kk = k-1, 1, -1
516	DO i = 1, klon
517	IF (rneb(i,k).GT.0.0.and.zprec_cond(i).gt.0.) THEN
518	if (t(i,kk) .GE. ztglace) THEN
519	zalpha_tr = a_tr_sca(1)
520	else
521	zalpha_tr = a_tr_sca(2)
522	endif
523	zfrac_lessi = 1. - EXP(zalpha_tr*zprec_cond(i)/zneb(i))
524	pfrac_impa(i,kk)=pfrac_impa(i,kk)(1.-zneb(i)zfrac_lessi)
525	frac_impa(i,kk)= 1.-zneb(i)*zfrac_lessi
526	ENDIF
527	ENDDO
528	ENDDO
529	c
530	cAA----------------------------------------------------------
531	c FIN DE BOUCLE SUR K
532	9999 CONTINUE
533	c
534	cAA-----------------------------------------------------------
535	c
536	c Pluie ou neige au sol selon la temperature de la 1ere couche
537	c
538	DO i = 1, klon
539	IF ((t(i,1)+d_t(i,1)) .LT. RTT) THEN
540	snow(i) = zrfl(i)
541	zlh_solid(i) = RLSTT-RLVTT
542	ELSE
543	rain(i) = zrfl(i)
544	zlh_solid(i) = 0.
545	ENDIF
546	ENDDO
547	C
548	C For energy conservation : when snow is present, the solification
549	c latent heat is considered.
550	DO k = 1, klev
551	DO i = 1, klon
552	zcpair=RCPD(1.0+RVTMP2(q(i,k)+d_q(i,k)))
553	zmair=(paprs(i,k)-paprs(i,k+1))/RG
554	zm_solid = (prfl(i,k)-prfl(i,k+1)+psfl(i,k)-psfl(i,k+1))*dtime
555	d_t(i,k) = d_t(i,k) + zlh_solid(i) zm_solid / (zcpairzmair)
556	END DO
557	END DO
558	c
559
560	if (ncoreczq>0) then
561	print*,'WARNING : ZQ dans fisrtilp ',ncoreczq,' val < 1.e-15.'
562	endif
563	RETURN
564	END

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

Download in other formats:

Original Format