Context Navigation

cvltr_spl.F90 @ 5411

Last change on this file since 5411 was 5160, checked in by abarral, 5 months ago
Put .h into modules
File size: 30.2 KB

Rev	Line
[5099]	1
[2147]	2	! $Id $
[5099]	3
[2147]	4	SUBROUTINE cvltr_spl(pdtime, da, phi,phi2,d1a,dam, mpIN,epIN, &
	5	sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm, &
	6	pmflxrIN,pmflxsIN,ev,te,wdtrainA,wdtrainM, &
	7	paprs,it,tr,upd,dnd,inb,icb, &
	8	kk,henry,zrho, ccntrAA_spla,ccntrENV_spla,coefcoli_spla, &
	9	id_prec,id_fine,id_coss, id_codu, id_scdu, &
	10	dtrcv,trsptd,dtrSscav,dtrsat,dtrUscav,qDi,qPr, &
	11	qPa,qMel,qTrdi,dtrcvMA,Mint, &
	12	zmfd1a,zmfphi2,zmfdam)
	13	USE IOIPSL
	14	USE dimphy
[5101]	15	USE infotrac_phy, ONLY: nbtr
[5140]	16	USE lmdz_conema3
[5142]	17	USE lmdz_YOECUMF
[5144]	18	USE lmdz_yomcst
[5160]	19	USE lmdz_chem, ONLY: idms, iso2, iso4, ih2s, idmso, imsa, ih2o2, &
	20	n_avogadro, masse_s, masse_so4, rho_water, rho_ice
[5140]	21
[2147]	22	IMPLICIT NONE
	23	!=====================================================================
	24	! Objet : convection des traceurs / KE
	25	! Auteurs: M-A Filiberti and J-Y Grandpeix
	26	! modifiee par R Pilon : lessivage des traceurs / KE
	27	!=====================================================================
	28
	29	! Entree
	30	REAL,INTENT(IN) :: pdtime
	31	REAL,DIMENSION(klon,klev),INTENT(IN) :: da
	32	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi
	33	! RomP
	34	REAL,DIMENSION(klon,klev),INTENT(IN) :: d1a,dam ! matrices pour simplifier
	35	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2 ! l'ecriture des tendances
[5099]	36
[2147]	37	REAL,DIMENSION(klon,klev),INTENT(IN) :: mpIN
	38	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs ! pression aux 1/2 couches (bas en haut)
	39	! REAL,DIMENSION(klon,klev),INTENT(IN) :: pplay ! pression aux 1/2 couches (bas en haut)
	40	REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: tr ! q de traceur (bas en haut)
	41	INTEGER,INTENT(IN) :: it
	42	REAL,DIMENSION(klon,klev),INTENT(IN) :: upd ! saturated updraft mass flux
	43	REAL,DIMENSION(klon,klev),INTENT(IN) :: dnd ! saturated downdraft mass flux
[5099]	44
[2147]	45	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainA ! masses precipitantes de l'asc adiab
	46	REAL,DIMENSION(klon,klev),INTENT(IN) :: wdtrainM ! masses precipitantes des melanges
	47	!JE REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxrIN ! vprecip: eau
	48	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxrIN ! vprecip: eau
	49	!JE REAL,DIMENSION(klon,klev),INTENT(IN) :: pmflxsIN ! vprecip: neige
	50	REAL,DIMENSION(klon,klev+1),INTENT(IN) :: pmflxsIN ! vprecip: neige
	51	REAL,DIMENSION(klon,klev),INTENT(IN) :: ev ! evaporation cv30_routine
	52	REAL,DIMENSION(klon,klev),INTENT(IN) :: epIN
	53	REAL,DIMENSION(klon,klev),INTENT(IN) :: te
	54	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij ! fraction dair de lenv
	55	REAL,DIMENSION(klon,klev),INTENT(IN) :: wght_cvfd ! weights of the layers feeding convection
	56	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij ! contenu en eau condensée spécifique/conc deau condensée massique
	57	REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm ! eau condensee precipitee dans mel masse dair sat
	58	REAL,DIMENSION(klon,klev),INTENT(IN) :: eplaMm ! eau condensee precipitee dans aa masse dair sat
	59
	60	REAL,DIMENSION(klon,klev),INTENT(IN) :: clw ! contenu en eau condensée dans lasc adiab
	61	REAL,DIMENSION(klon),INTENT(IN) :: sigd
	62	INTEGER,DIMENSION(klon),INTENT(IN) :: icb,inb
	63	! Sortie
	64	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcv ! tendance totale (bas en haut)
	65	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrcvMA ! M-A Filiberti
	66	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: trsptd ! tendance du transport
	67	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrSscav ! tendance du lessivage courant sat
	68	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrsat ! tendance trsp+sat scav
	69	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: dtrUscav ! tendance du lessivage courant unsat
[5099]	70
[2147]	71	! Variables locales
	72	INTEGER :: i,j,k
	73	REAL,DIMENSION(klon,klev) :: dxpres ! difference de pression entre niveau (j+1) et (j)
	74	REAL :: pdtimeRG ! pas de temps * gravite
	75	REAL,DIMENSION(klon,nbtr) :: qfeed ! tracer concentration feeding convection
	76	! variables pour les courants satures
	77	REAL,DIMENSION(klon,klev,klev) :: zmd
	78	REAL,DIMENSION(klon,klev,klev) :: za
	79	REAL,DIMENSION(klon,klev,nbtr) :: zmfd,zmfa
	80	REAL,DIMENSION(klon,klev,nbtr) :: zmfp,zmfu
	81
	82	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfd1a
	83	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfdam
	84	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: zmfphi2
	85
	86	! RomP ! les variables sont nettoyees des valeurs aberrantes
	87	REAL,DIMENSION(klon,klev) :: Pa, Pm ! pluie AA et mélanges, var temporaire
	88	REAL,DIMENSION(klon,klev) :: pmflxs,pmflxr ! pmflxrIN,pmflxsIN sans valeur aberante
	89	REAL,DIMENSION(klon,klev) :: mp ! flux de masse
	90	REAL,DIMENSION(klon,klev) :: ep ! fraction d'eau convertie en precipitation
	91	REAL,DIMENSION(klon,klev) :: evap ! evaporation : variable temporaire
	92	REAL,DIMENSION(klon,klev) :: rho !environmental density
	93
	94	REAL,DIMENSION(klon,klev) :: kappa ! denominateur du au calcul de la matrice
	95	! pour obtenir qd et qp
	96	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qTrdi ! traceurs descente air insature transport MA
	97	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qDi ! traceurs descente insaturees
	98	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPr ! traceurs colonne precipitante
	99	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qPa ! traceurs dans les precip issues lasc. adiab.
	100	REAL,DIMENSION(klon,klev,nbtr),INTENT(OUT) :: qMel ! traceurs dans les precip issues des melanges
	101	REAL,DIMENSION(klon,klev,nbtr) :: qMeltmp ! variable temporaire
	102	REAL,DIMENSION(klon,klev,nbtr) :: qpmMint
	103	REAL,DIMENSION(klon,klev),INTENT(OUT) :: Mint
	104	! tendances
	105	REAL :: tdcvMA ! terme de transport de traceur (schema Marie Angele)
	106	REAL :: trsptrac ! terme de transport de traceur par l'air
	107	REAL :: scavtrac ! terme de lessivage courant sature
	108	REAL :: uscavtrac ! terme de lessivage courant insature
	109	! impaction
	110	!!! Correction apres discussion Romain P. / Olivier B.
	111	!!! REAL,PARAMETER :: rdrop=2.5e-3 ! rayon des gouttes d'eau
	112	REAL,PARAMETER :: rdrop=1.e-3 ! rayon des gouttes d'eau
	113	!!!
	114	REAL,DIMENSION(klon,klev) :: imp ! coefficient d'impaction
	115	! parametres lessivage
	116	REAL :: ccntrAA_coef ! \alpha_a : fract aerosols de l'AA convertis en CCN
	117	REAL :: ccntrENV_coef ! \beta_m : fract aerosols de l'env convertis en CCN
	118	REAL :: coefcoli ! coefficient de collision des gouttes sur les aerosols
[5099]	119
[2147]	120	LOGICAL,DIMENSION(klon,klev) :: NO_precip
	121	! LOGICAL :: scavON
	122	! var tmp tests
	123	REAL :: conserv
[5158]	124	REAL :: conservMA
[2147]	125	! JE SPLA adaptation
	126
	127	INTEGER :: id_prec,id_fine,id_coss, id_codu, id_scdu
	128	REAL,DIMENSION(nbtr) :: ccntrAA_spla,ccntrENV_spla,coefcoli_spla
	129	REAL,DIMENSION(klon,klev) :: ccntrAA_coef3d
	130	REAL,DIMENSION(klon,klev) :: ccntrENV_coef3d
	131	REAL,DIMENSION(klon,klev) :: coefcoli3d
	132
	133	REAL,DIMENSION(nbtr) :: henry !--cste de Henry mol/l/atm
	134	REAL,DIMENSION(nbtr) :: kk !--coefficient de var avec T (K)
	135	REAL :: henry_t !--constante de Henry a T t (mol/l/atm)
	136	REAL :: f_a !--rapport de la phase aqueuse a la phase gazeuse
	137
	138	REAL, PARAMETER :: ph=5.
	139	REAL :: K1, K2
	140	REAL,DIMENSION(klon,klev) :: zrho
	141	REAL, PARAMETER :: qliq=1.e-3 ! CONVECTIVE ONLY
	142
	143	! Je end 20140105
	144
	145	!JE20140724<<
	146	!JE SPLA adapt:
[5099]	147
	148
[2147]	149	!! coefficient lessivage
	150	! ccntrAA_coef = 0.
	151	! ccntrENV_coef = 0.
	152	! coefcoli = 0.
[5099]	153
[2147]	154	!!$OMP MASTER
[5101]	155	! CALL getin('ccntrAA_coef',ccntrAA_coef)
	156	! CALL getin('ccntrENV_coef',ccntrENV_coef)
	157	! CALL getin('coefcoli',coefcoli)
[2147]	158	!!$OMP END MASTER
	159	!!$OMP BARRIER
	160	!!! JE
	161	! do j=1,klev
	162	! do i=1,klon
	163	! imp(i,j)=0.
	164	! ccntrAA_coef3d(i,j)=ccntrAA_coef
	165	! ccntrENV_coef3d(i,j)=ccntrENV_coef
	166	! coefcoli3d(i,j)=coefcoli
	167	! enddo
	168	! enddo
[5099]	169
	170
[2147]	171	!!! for SPLA
	172	!!!JEend
[5103]	173	! PRINT*,'cvltr coef lessivage convectif', ccntrAA_coef,ccntrENV_coef,coefcoli
[5099]	174
[2147]	175	!JE20140724>>
	176
	177	! scavON=.TRUE.
[5116]	178	! IF(scavON) THEN
[2147]	179	! ccntrAA_coef = 1.
	180	! ccntrENV_coef = 1.
	181	! coefcoli = 1.
	182	! else
	183	! ccntrAA_coef = 0.
	184	! ccntrENV_coef = 0.
	185	! coefcoli = 0.
[5117]	186	! ENDIF
[2147]	187
	188	! ======================================================
	189	! calcul de l'impaction
	190	! ======================================================
	191	!initialisation
[5158]	192	DO j=1,klev
	193	DO i=1,klon
[2147]	194	imp(i,j)=0.
	195	enddo
	196	enddo
	197	!JE init 20140103
	198	!! impaction sur la surface de la colonne de la descente insaturee
	199	!! On prend la moyenne des precip entre le niveau i+1 et i
	200	!! I=3/4* (P(1+1)+P(i))/2 / (sigdrrho_l)
	201	!! 1000kg/m3= densité de l'eau
	202	!! 0.75e-3 = 3/4 /1000
	203	!! Par la suite, I est tout le temps multiplié par sig_d pour avoir l'impaction sur la surface de la maille
	204	!! on le néglige ici pour simplifier le code
	205	! do j=1,klev-1
	206	! do i=1,klon
	207	! imp(i,j) = coefcoli0.75e-3/rdrop &
	208	! 0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
	209	!! rho(i,j)=pplay(i,j)/(rd*te(i,j))
	210	! enddo
	211	! enddo
	212	!JE end 20140103
	213
	214	! initialisation pour flux de traceurs, td et autre
	215	trsptrac = 0.
	216	scavtrac = 0.
	217	uscavtrac = 0.
	218	qfeed(:,it) = 0. !RL
	219	DO j=1,klev
	220	DO i=1,klon
	221	zmfd(i,j,it)=0.
	222	zmfa(i,j,it)=0.
	223	zmfu(i,j,it)=0.
	224	zmfp(i,j,it)=0.
	225	zmfphi2(i,j,it)=0.
	226	zmfd1a(i,j,it)=0.
	227	zmfdam(i,j,it)=0.
	228	qDi(i,j,it)=0.
	229	qPr(i,j,it)=0.
	230	qPa(i,j,it)=0.
	231	qMel(i,j,it)=0.
	232	qMeltmp(i,j,it)=0.
	233	qTrdi(i,j,it)=0.
	234	kappa(i,j)=0.
	235	trsptd(i,j,it)=0.
	236	dtrsat(i,j,it)=0.
	237	dtrSscav(i,j,it)=0.
	238	dtrUscav(i,j,it)=0.
	239	dtrcv(i,j,it)=0.
	240	dtrcvMA(i,j,it)=0.
	241	evap(i,j)=0.
	242	dxpres(i,j)=0.
	243	qpmMint(i,j,it)=0.
	244	Mint(i,j)=0.
	245	END DO
	246	END DO
	247
	248	! suppression des valeurs très faibles (~1e-320)
	249	! multiplication de levaporation pour lavoir par unite de temps
	250	! et par unite de surface de la maille
	251	! -> cv30_unsat : evap : masse evaporee/s/(m2 de la descente)
	252	DO j=1,klev
	253	DO i=1,klon
[5116]	254	IF(ev(i,j)<1.e-16) THEN
[2147]	255	evap(i,j)=0.
	256	else
	257	evap(i,j)=ev(i,j)*sigd(i)
	258	endif
	259	END DO
	260	END DO
	261
	262	DO j=1,klev
	263	DO i=1,klon
[5116]	264	IF(j<klev) THEN
	265	IF(epIN(i,j)<1.e-32) THEN
[2147]	266	ep(i,j)=0.
	267	else
	268	ep(i,j)=epIN(i,j)
	269	endif
	270	else
	271	ep(i,j)=epmax
	272	endif
[5116]	273	IF(mpIN(i,j)<1.e-32) THEN
[2147]	274	mp(i,j)=0.
	275	else
	276	mp(i,j)=mpIN(i,j)
	277	endif
[5116]	278	IF(pmflxsIN(i,j)<1.e-32) THEN
[2147]	279	pmflxs(i,j)=0.
	280	else
	281	pmflxs(i,j)=pmflxsIN(i,j)
	282	endif
[5116]	283	IF(pmflxrIN(i,j)<1.e-32) THEN
[2147]	284	pmflxr(i,j)=0.
	285	else
	286	pmflxr(i,j)=pmflxrIN(i,j)
	287	endif
[5116]	288	IF(wdtrainA(i,j)<1.e-32) THEN
[2147]	289	Pa(i,j)=0.
	290	else
	291	Pa(i,j)=wdtrainA(i,j)
	292	endif
[5116]	293	IF(wdtrainM(i,j)<1.e-32) THEN
[2147]	294	Pm(i,j)=0.
	295	else
	296	Pm(i,j)=wdtrainM(i,j)
	297	endif
	298	END DO
	299	END DO
	300
	301	!==========================================
	302	DO j = klev-1,1,-1
	303	DO i = 1,klon
[5082]	304	NO_precip(i,j) = (pmflxr(i,j+1)+pmflxs(i,j+1))<1.e-10&
[5117]	305	.AND.Pa(i,j)<1.e-10.AND.Pm(i,j)<1.e-10
[2147]	306	END DO
	307	END DO
	308	!==============================================================================
	309	! JE SPLA: Calc of ccntrAA_coef,ccntrENV_coef, coefcoli for SPLA aerosols and
	310	! precursors. From SPLA code inscav_spl.F
[5160]	311	!PRINT *,'Using SPLA values for cvltr_spl ccntr and coefcoli params'
[5099]	312
	313
[2147]	314	!IF (it.EQ.2) THEN !--fine aerosol
	315	! DO j=1,klev
	316	! DO i=1,klon
	317	! ccntrAA_coef3d(i,j)=0.7
	318	! ccntrENV_coef3d(i,j)=0.7
	319	! coefcoli3d(i,j)=0.001
	320	! ENDDO
	321	! ENDDO
	322	!ELSEIF (it.EQ.3) THEN !-- Coarse Sea salt aerosol
	323	! DO j=1,klev
	324	! DO i=1,klon
	325	! ccntrAA_coef3d(i,j)=1.0
	326	! ccntrENV_coef3d(i,j)=1.0
	327	! coefcoli3d(i,j)=0.001
	328	! ENDDO
	329	! ENDDO
[5099]	330
[2147]	331	!ELSEIF (it.EQ.4) THEN !--Coarse Dust aerosol
	332	! DO j=1,klev
	333	! DO i=1,klon
	334	! ccntrAA_coef3d(i,j)=0.7
	335	! ccntrENV_coef3d(i,j)=0.7
	336	! coefcoli3d(i,j)=0.001
[5099]	337
[2147]	338	! ENDDO
	339	! ENDDO
	340	! Gas precursor: Henry's law
	341
[5082]	342	IF (it == id_prec) THEN
[2147]	343	DO k=1, klev
	344	DO i=1, klon
	345	henry_t=henry(it)exp(-kk(it)(1./298.-1./te(i,k))) !--mol/l/atm
	346	K1=1.2e-2exp(-2010(1/298.-1/te(i,k)))
	347	K2=6.6e-8exp(-1510(1/298.-1/te(i,k)))
	348	henry_t=henry_t(1. + K1/10.(-ph) + K1K2/(10.(-ph))2)
	349	f_a=henry_t/101.325Rte(i,k)qliqzrho(i,k)/rho_water
	350	! scav(i,k)=f_a/(1.+f_a)
	351	ccntrAA_coef3d(i,k)= f_a/(1.+f_a)
	352	ccntrENV_coef3d(i,k)= f_a/(1.+f_a)
	353	coefcoli3d(i,k)=0.0
	354	ENDDO
	355	ENDDO
	356	! CALL minmaxqfi2(clw,1.e33,-1.e33,'minmax clw')
	357	ELSE
	358	DO j=1,klev
	359	DO i=1,klon
	360	ccntrAA_coef3d(i,j)=ccntrAA_spla(it)
	361	ccntrENV_coef3d(i,j)=ccntrENV_spla(it)
	362	coefcoli3d(i,j)=coefcoli_spla(it)
	363	ENDDO
	364	ENDDO
	365
	366
	367	ENDIF
	368
	369	! JE end SPLA modifs in ccn parameters
	370	!==============================================================================
	371
	372
	373
	374
	375	!JE init 20140103
	376	! impaction sur la surface de la colonne de la descente insaturee
	377	! On prend la moyenne des precip entre le niveau i+1 et i
	378	! I=3/4* (P(1+1)+P(i))/2 / (sigdrrho_l)
	379	! 1000kg/m3= densité de l'eau
	380	! 0.75e-3 = 3/4 /1000
	381	! Par la suite, I est tout le temps multiplié par sig_d pour avoir l'impaction sur la surface de la maille
	382	! on le néglige ici pour simplifier le code
[5158]	383	DO j=1,klev-1
	384	DO i=1,klon
[2147]	385	!JE imp(i,j) = coefcoli0.75e-3/rdrop &
	386	!JE 0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
	387	imp(i,j) = coefcoli3d(i,j)0.75e-3/rdrop &
	388	0.5*(pmflxr(i,j+1)+pmflxs(i,j+1)+pmflxr(i,j)+pmflxs(i,j))
	389	! rho(i,j)=pplay(i,j)/(rd*te(i,j))
	390	enddo
	391	enddo
	392	!JE end 20140103
	393
	394
	395	! =========================================
	396	! calcul des tendances liees au downdraft
	397	! =========================================
	398	!cdir collapse
	399	DO k=1,klev
	400	DO j=1,klev
	401	DO i=1,klon
	402	zmd(i,j,k)=0.
	403	za (i,j,k)=0.
	404	END DO
	405	END DO
	406	END DO
	407	! calcul de la matrice d echange
	408	! matrice de distribution de la masse entrainee en k
	409	! commmentaire RomP : mp > 0
	410	DO k=1,klev-1
	411	DO i=1,klon
	412	zmd(i,k,k)=max(0.,mp(i,k)-mp(i,k+1)) ! ~ mk(k)
	413	END DO
	414	END DO
	415	DO k=2,klev
	416	DO j=k-1,1,-1
	417	DO i=1,klon
[5116]	418	IF(mp(i,j+1)>1.e-10) THEN
[2147]	419	zmd(i,j,k)=zmd(i,j+1,k)min(1.,mp(i,j)/mp(i,j+1)) !det ~ mk(j)=mk(j+1)mp(i,j)/mp(i,j+1)
	420	ENDif
	421	END DO
	422	END DO
	423	END DO
	424	DO k=1,klev
	425	DO j=1,klev-1
	426	DO i=1,klon
	427	za(i,j,k)=max(0.,zmd(i,j+1,k)-zmd(i,j,k))
	428	END DO
	429	END DO
	430	END DO
	431	!!!!! quantite de traceur dans la descente d'air insaturee : 4 juin 2012
	432	DO k=1,klev
	433	DO j=1,klev-1
	434	DO i=1,klon
[5116]	435	IF(mp(i,j+1)>1.e-10) THEN
[2147]	436	qTrdi(i,j+1,it)=qTrdi(i,j+1,it)+(zmd(i,j+1,k)/mp(i,j+1))*tr(i,k,it)
	437	else
	438	qTrdi(i,j,it)=0.!tr(i,j,it)
	439	endif
	440	ENDDO
	441	ENDDO
	442	ENDDO
	443	!!!!!
[5099]	444
[2147]	445	! rajout du terme lie a l ascendance induite
[5099]	446
[2147]	447	DO j=2,klev
	448	DO i=1,klon
	449	za(i,j,j-1)=za(i,j,j-1)+mp(i,j)
	450	END DO
	451	END DO
[5099]	452
[2147]	453	! tendance courants insatures ! sans lessivage ancien schema
[5099]	454
[2147]	455	DO k=1,klev
	456	DO j=1,klev
	457	DO i=1,klon
	458	zmfd(i,j,it)=zmfd(i,j,it)+za(i,j,k)*(tr(i,k,it)-tr(i,j,it))
	459	END DO
	460	END DO
	461	END DO
[5099]	462
[2147]	463	! =========================================
	464	! calcul des tendances liees aux courants satures j <-> z ; k <-> z'
	465	! =========================================
	466	!RL
	467	! Feeding concentrations
	468	DO j=1,klev
	469	DO i=1,klon
	470	qfeed(i,it)=qfeed(i,it)+wght_cvfd(i,j)*tr(i,j,it)
	471	END DO
	472	END DO
	473	!RL
[5099]	474
[2147]	475	DO j=1,klev
	476	DO i=1,klon
	477	!RL
	478	!! zmfa(i,j,it)=da(i,j)*(tr(i,1,it)-tr(i,j,it)) ! da
	479	zmfa(i,j,it)=da(i,j)*(qfeed(i,it)-tr(i,j,it)) ! da
	480	!RL
	481	END DO
	482	END DO
[5099]	483
[2147]	484	DO k=1,klev
	485	DO j=1,klev
	486	DO i=1,klon
	487	zmfp(i,j,it)=zmfp(i,j,it)+phi(i,j,k)*(tr(i,k,it)-tr(i,j,it)) ! phi
	488	END DO
	489	END DO
	490	END DO
	491	! RomP ajout des matrices liees au lessivage
	492	DO j=1,klev
	493	DO i=1,klon
	494	zmfd1a(i,j,it)=d1a(i,j)*tr(i,1,it) ! da1
	495	zmfdam(i,j,it)=dam(i,j)*tr(i,1,it) ! dam
	496	END DO
	497	END DO
	498	DO k=1,klev
	499	DO j=1,klev
	500	DO i=1,klon
	501	zmfphi2(i,j,it)=zmfphi2(i,j,it)+phi2(i,j,k)*tr(i,k,it) ! psi
	502	END DO
	503	END DO
	504	END DO
	505	DO j=1,klev-1
	506	DO i=1,klon
	507	zmfu(i,j,it)=max(0.,upd(i,j+1)+dnd(i,j+1))*(tr(i,j+1,it)-tr(i,j,it))
	508	END DO
	509	END DO
	510	DO j=2,klev
	511	DO i=1,klon
	512	zmfu(i,j,it)=zmfu(i,j,it)+min(0.,upd(i,j)+dnd(i,j))*(tr(i,j,it)-tr(i,j-1,it))
	513	END DO
	514	END DO
	515	! ===================================================
	516	! calcul des tendances liees aux courants insatures
	517	! ===================================================
	518	! pression
	519	DO k=1, klev
	520	DO i=1, klon
	521	dxpres(i,k)=paprs(i,k)-paprs(i,k+1)
	522	ENDDO
	523	ENDDO
	524	pdtimeRG=pdtime*RG
	525
	526	! q_pa et q_pm traceurs issues des courants satures se retrouvant dans les precipitations
	527	DO j=1,klev
	528	DO i=1,klon
[5117]	529	IF(j>=icb(i).AND.j<=inb(i)) THEN
[5116]	530	IF(clw(i,j)>1.e-16) THEN
[2147]	531	! qPa(i,j,it)=ccntrAA_coef*tr(i,1,it)/clw(i,j)
	532	qPa(i,j,it)=ccntrAA_coef3d(i,j)*tr(i,1,it)/clw(i,j)
	533	else
	534	qPa(i,j,it)=0.
	535	endif
	536	endif
	537	END DO
	538	END DO
	539
	540	! calcul de q_pm en 2 parties :
	541	! 1) calcul de sa valeur pour un niveau z' donne
	542	! 2) integration sur la verticale sur z'
	543	DO j=1,klev
	544	DO k=1,j-1
	545	DO i=1,klon
[5117]	546	IF(k>=icb(i).AND.k<=inb(i).AND.&
[5116]	547	j<=inb(i)) THEN
	548	IF(elij(i,k,j)>1.e-16) THEN
[2147]	549	!JE qMeltmp(i,j,it)=((1-ep(i,k))ccntrAA_coeftr(i,1,it)&
	550	!JE *(1.-sij(i,k,j)) +ccntrENV_coef&
	551	!JE tr(i,k,it)sij(i,k,j)) / elij(i,k,j)
	552	qMeltmp(i,j,it)=((1-ep(i,k))ccntrAA_coef3d(i,k)tr(i,1,it)&
	553	*(1.-sij(i,k,j)) +ccntrENV_coef3d(i,k)&
	554	tr(i,k,it)sij(i,k,j)) / elij(i,k,j)
	555	else
	556	qMeltmp(i,j,it)=0.
	557	endif
	558	qpmMint(i,j,it)=qpmMint(i,j,it) + qMeltmp(i,j,it)*epmlmMm(i,j,k)
	559	Mint(i,j)=Mint(i,j) + epmlmMm(i,j,k)
	560	endif ! end if dans nuage
	561	END DO
	562	END DO
	563	END DO
	564
	565	DO j=1,klev
	566	DO i=1,klon
[5116]	567	IF(Mint(i,j)>1.e-16) THEN
[2147]	568	qMel(i,j,it)=qpmMint(i,j,it)/Mint(i,j)
	569	else
	570	qMel(i,j,it)=0.
	571	endif
	572	END DO
	573	END DO
	574
	575	! calcul de q_d et q_p traceurs de la descente precipitante
	576	DO j=klev-1,1,-1
	577	DO i=1,klon
[5117]	578	IF(mp(i,j+1)>mp(i,j).AND.mp(i,j+1)>1.e-10) then ! detrainement
[2147]	579	kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	580	(-mp(i,j+1)-imp(i,j)/RG*dxpres(i,j))&
	581	+ (imp(i,j)/RGdxpres(i,j))(evap(i,j)/RG*dxpres(i,j)))
	582
[5117]	583	elseif(mp(i,j)>mp(i,j+1).AND.mp(i,j)>1.e-10) then! entrainement
[5116]	584	IF(j==1) THEN
[2147]	585	kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	586	(-mp(i,2)-imp(i,j)/RG*dxpres(i,j))&
	587	+ (imp(i,j)/RGdxpres(i,j))(evap(i,j)/RG*dxpres(i,j)))
	588	else
	589	kappa(i,j)=((pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	590	(-mp(i,j)-imp(i,j)/RG*dxpres(i,j))&
	591	+ (imp(i,j)/RGdxpres(i,j))(evap(i,j)/RG*dxpres(i,j)))
	592	endif
	593	else
	594	kappa(i,j)=1.
	595	endif
	596	ENDDO
	597	ENDDO
	598
	599	DO j=klev-1,1,-1
	600	DO i=1,klon
[5117]	601	IF (abs(kappa(i,j))<1.e-25) then !si denominateur nul (il peut y avoir des mp!=0)
[2147]	602	kappa(i,j)=1.
[5116]	603	IF(j==1) THEN
[2147]	604	qDi(i,j,it)=qDi(i,j+1,it) !orig tr(i,j,it) ! mp(1)=0 donc tout vient de la couche supérieure
[5117]	605	elseif(mp(i,j+1)>mp(i,j).AND.mp(i,j+1)>1.e-10) THEN
[2147]	606	qDi(i,j,it)=qDi(i,j+1,it)
[5117]	607	elseif(mp(i,j)>mp(i,j+1).AND.mp(i,j)>1.e-10) then! entrainement
[2147]	608	qDi(i,j,it)=(-mp(i,j+1)(qDi(i,j+1,it)-tr(i,j,it))-mp(i,j)tr(i,j,it))/(-mp(i,j))
	609	else ! si mp (i)=0 et mp(j+1)=0
	610	qDi(i,j,it)=tr(i,j,it) ! orig 0.
	611	endif
	612
[5116]	613	IF(NO_precip(i,j)) THEN
[2147]	614	qPr(i,j,it)=0.
	615	else
	616	qPr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	617	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)&
	618	+imp(i,j)/RGdxpres(i,j)qDi(i,j,it))/&
	619	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))
	620	endif
	621	else ! denominateur non nul
	622	kappa(i,j)=1./kappa(i,j)
	623	! calcul de qd et qp
	624	!!jyg (20130119) correction pour le sommet du nuage
[5116]	625	!! IF(j.ge.inb(i)) then !au-dessus du nuage, sommet inclu
	626	IF(j>inb(i)) then !au-dessus du nuage
[2147]	627	qDi(i,j,it)=tr(i,j,it) ! pas de descente => environnement = descente insaturee
	628	qPr(i,j,it)=0.
	629
	630	! vvv premiere couche du modele ou mp(1)=0 ! det tout le temps vvv
[5116]	631	elseif(j==1) THEN
	632	IF(mp(i,2)>1.e-10) then !mp(2) non nul -> detrainement (car mp(1) = 0) !ent pas possible
	633	IF(NO_precip(i,j)) then !pas de precip en (i)
[2147]	634	qDi(i,j,it)=qDi(i,j+1,it)
	635	qPr(i,j,it)=0.
	636	else
	637	qDi(i,j,it)=kappa(i,j)*(&
	638	(-evap(i,j)/RGdxpres(i,j))((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	639	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)) +&
	640	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	641	(-mp(i,j+1)*qDi(i,j+1,it)))
	642
	643	qPr(i,j,it)=kappa(i,j)*(&
	644	(-mp(i,j+1)-imp(i,j)/RGdxpres(i,j))&
	645	((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	646	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it))&
	647	+(-mp(i,j+1)qDi(i,j+1,it)) (imp(i,j)/RG*dxpres(i,j)))
	648	endif
	649
	650	else !mp(2) nul -> plus de descente insaturee -> pluie agit sur environnement
	651	qDi(i,j,it)=tr(i,j,it) ! orig 0.
[5116]	652	IF(NO_precip(i,j)) THEN
[2147]	653	qPr(i,j,it)=0.
	654	else
	655	qPr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	656	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)&
	657	+imp(i,j)/RGdxpres(i,j)tr(i,j,it))/&
	658	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))
	659	endif
	660
	661	endif !mp(2) nul ou non
	662
[5117]	663	! vvv (j!=1.AND.j.lt.inb(i)) en-dessous du sommet nuage vvv
[2147]	664	else
	665	!------------------------------------------------------------- detrainement
[5117]	666	IF(mp(i,j+1)>mp(i,j).AND.mp(i,j+1)>1.e-10) then !mp(i,j).gt.1.e-10) THEN
[5116]	667	IF(NO_precip(i,j)) THEN
[2147]	668	qDi(i,j,it)=qDi(i,j+1,it)
	669	qPr(i,j,it)=0.
	670	else
	671	qDi(i,j,it)=kappa(i,j)*(&
	672	(-evap(i,j)/RGdxpres(i,j))((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	673	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)) +&
	674	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	675	(-mp(i,j+1)*qDi(i,j+1,it)))
[5099]	676
[2147]	677	qPr(i,j,it)=kappa(i,j)*(&
	678	(-mp(i,j+1)-imp(i,j)/RGdxpres(i,j))&
	679	((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	680	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it))&
	681	+(-mp(i,j+1)qDi(i,j+1,it)) (imp(i,j)/RG*dxpres(i,j)))
	682	endif !precip
	683	!------------------------------------------------------------- entrainement
[5117]	684	elseif(mp(i,j)>mp(i,j+1).AND.mp(i,j)>1.e-10) THEN
[5116]	685	IF(NO_precip(i,j)) THEN
[2147]	686	qDi(i,j,it)=(-mp(i,j+1)(qDi(i,j+1,it)-tr(i,j,it))-mp(i,j)tr(i,j,it))/(-mp(i,j))
	687	qPr(i,j,it)=0.
	688	else
	689	qDi(i,j,it)=kappa(i,j)*(&
	690	(-evap(i,j)/RGdxpres(i,j))((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	691	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)) +&
	692	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))*&
	693	(-mp(i,j+1)(qDi(i,j+1,it)-tr(i,j,it))-mp(i,j)tr(i,j,it)))
[5099]	694
[2147]	695	qPr(i,j,it)=kappa(i,j)*(&
	696	(-mp(i,j)-imp(i,j)/RGdxpres(i,j))&
	697	((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	698	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it))&
	699	+(-mp(i,j+1)(qDi(i,j+1,it)-tr(i,j,it))-mp(i,j)tr(i,j,it))*&
	700	(imp(i,j)/RG*dxpres(i,j)))
	701	endif !precip
[5117]	702	!------------------------------------------------------------- END IF ! ent/det
[2147]	703	else !mp nul
	704	qDi(i,j,it)=tr(i,j,it) ! orig 0.
[5116]	705	IF(NO_precip(i,j)) THEN
[2147]	706	qPr(i,j,it)=0.
	707	else
	708	qPr(i,j,it)=((pmflxr(i,j+1)+pmflxs(i,j+1))*qPr(i,j+1,it)+&
	709	Pa(i,j)qPa(i,j,it)+Pm(i,j)qMel(i,j,it)&
	710	+imp(i,j)/RGdxpres(i,j)tr(i,j,it))/&
	711	(pmflxr(i,j+1)+pmflxs(i,j+1)+Pa(i,j)+Pm(i,j))
	712	endif
	713	endif ! mp nul ou non
	714	endif ! condition sur j
	715	endif ! kappa
	716	ENDDO
	717	ENDDO
	718
	719	!! print test descente insaturee
	720	! DO j=klev,1,-1
	721	! DO i=1,klon
[5117]	722	! IF(it.EQ.3) THEN
[5116]	723	! WRITE(*,'(I2,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12)') j,&
[2147]	724	!! 'zmfdam',zmfdam(i,j,it),'zmfpsi',zmfphi2(i,j,it),&
	725	! 'zmfdam+zmfpsi',zmfdam(i,j,it)+zmfphi2(i,j,it),'qpmMint',qpmMint(i,j,it),&
	726	! 'Pm',Pm(i,j),'Mint',Mint(i,j),&
	727	!! 'zmfa',zmfa(i,j,it),'zmfp',zmfp(i,j,it),&
	728	! 'zmfdam',zmfdam(i,j,it),'zmfpsi',zmfphi2(i,j,it),'zmfd1a',zmfd1a(i,j,it)
	729	!! 'Pa',Pa(i,j),'eplaMm',eplaMm(i,j)
	730	!! 'zmfd1a=da1qa',zmfd1a(i,j,it),'PaqPa',wdtrainA(i,j)*qPa(i,j,it),'da1',d1a(i,j)
	731	! endif
	732	! ENDDO
	733	! ENDDO
	734
	735
	736	! ===================================================
	737	! calcul final des tendances
	738	! ===================================================
	739
	740	DO k=klev-1,1,-1
	741	DO i=1, klon
	742	! transport
	743	tdcvMA=zmfd(i,k,it)+zmfu(i,k,it)+zmfa(i,k,it)+zmfp(i,k,it) ! double comptage des downdraft insatures
	744	trsptrac=zmfu(i,k,it)+zmfa(i,k,it)+zmfp(i,k,it)
	745	! lessivage courants satures
	746	!JE scavtrac=-ccntrAA_coef*zmfd1a(i,k,it)&
	747	!JE -zmfphi2(i,k,it)*ccntrENV_coef&
	748	!JE -zmfdam(i,k,it)*ccntrAA_coef
	749	scavtrac=-ccntrAA_coef3d(i,k)*zmfd1a(i,k,it)&
	750	-zmfphi2(i,k,it)*ccntrENV_coef3d(i,k)&
	751	-zmfdam(i,k,it)*ccntrAA_coef3d(i,k)
	752	! lessivage courants insatures
[5117]	753	IF(k<=inb(i).AND.k>1) then ! tendances dans le nuage
[2147]	754	!------------------------------------------------------------- detrainement
[5117]	755	IF(mp(i,k+1)>mp(i,k).AND.mp(i,k+1)>1.e-10) THEN
[2147]	756	uscavtrac= (-mp(i,k)+mp(i,k+1))*(qDi(i,k,it)-tr(i,k,it))&
	757	+ mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
[5099]	758
[5117]	759	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,5X,e20.12,82X,a,e20.12)')k,' det incloud',&
[2147]	760	! (-mp(i,k)+mp(i,k+1))(qDi(i,k,it)-tr(i,k,it))pdtimeRG/dxpres(i,k)+&
	761	! mp(i,k)(tr(i,k-1,it)-tr(i,k,it))pdtimeRG/dxpres(i,k),&
	762	! 'mp',mp(i,k)
	763	!------------------------------------------------------------- entrainement
[5117]	764	elseif(mp(i,k)>mp(i,k+1).AND.mp(i,k)>1.e-10) THEN
[2147]	765	uscavtrac= mp(i,k)*(tr(i,k-1,it)-tr(i,k,it))
[5099]	766
[5117]	767	! IF(it.EQ.3) WRITE(,'(I2,1X,a,5X,e20.12,82X,a,e20.12)')k,' ent incloud',uscavtracpdtimeRG/dxpres(i,k), 'mp',mp(i,k)
[2147]	768	!=!------------------------------------------------------------- end ent/det
	769	else ! mp(i,k+1)=0. et mp(i,k)=0. pluie directement sur l environnement
	770
[5116]	771	IF(NO_precip(i,k)) THEN
[2147]	772	uscavtrac=0.
[5117]	773	! IF(it.EQ.3) WRITE(,'(I2,1X,a,e20.12,82X,a,e20.12)')k,' no P ent incloud',uscavtracpdtimeRG/dxpres(i,k), 'mp',mp(i,k)
[2147]	774	else
	775	uscavtrac=-imp(i,k)tr(i,k,it)dxpres(i,k)/RG+evap(i,k)qPr(i,k,it)dxpres(i,k)/RG
[5117]	776	! IF(it.EQ.3) WRITE(,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtracpdtimeRG/dxpres(i,k), 'mp',mp(i,k)
[2147]	777	!!JE adds
[5117]	778	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'imp',imp(i,k)
	779	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'tr',tr(i,k,it)
	780	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'evap',evap(i,k)
	781	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'qPr',qPr(i,k,it)
	782	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,3X,e20.12,82X,a,e20.12)')k,' P env incloud',uscavtrac, 'dxpres',dxpres(i,k)
[2147]	783	!!Je end
	784
	785	endif
	786	endif ! mp/det/ent
	787	!------------------------------------------------------------- premiere couche
[5082]	788	elseif(k==1) then ! mp(1)=0.
[5116]	789	IF(mp(i,2)>1.e-10) then !detrainement
[2147]	790	uscavtrac= (-0.+mp(i,2))*(qDi(i,k,it)-tr(i,k,it)) !&
	791	! + mp(i,2)*(0.-tr(i,k,it))
[5099]	792
[5117]	793	! IF(it.EQ.3) WRITE(*,'(I2,1X,a,e20.12,84X,a,e20.12)')k,' 1 det',&
[2147]	794	! (-0.+mp(i,2))(qDi(i,k,it)-tr(i,k,it))pdtimeRG/dxpres(i,k)+&
	795	! mp(i,2)(0.-tr(i,k,it))pdtimeRG/dxpres(i,k),&
	796	! 'mp',mp(i,k)
	797	else ! mp(2) = 0 = mp(1) pas de descente insaturee, rien ne se passe s'il ne pleut pas, sinon pluie->env
[5116]	798	IF(NO_precip(i,1)) THEN
[2147]	799	uscavtrac=0.
	800	else
	801	uscavtrac=-imp(i,k)tr(i,k,it)dxpres(i,k)/RG+evap(i,k)qPr(i,k,it)dxpres(i,k)/RG
	802	endif
[5117]	803	! IF(it.EQ.3) WRITE(,'(I2,1X,a,2X,e20.12,82X,a,e20.12)')k,'1 P env incloud',uscavtracpdtimeRG/dxpres(i,k), 'mp',mp(i,k)
[2147]	804	endif
	805
	806	else ! k > INB au-dessus du nuage
	807	uscavtrac=0.
	808	endif
	809
	810	! ===== tendances finales ======
	811	trsptd(i,k,it)=trsptrac*pdtimeRG/dxpres(i,k) ! td transport sans eau dans courants satures
	812	dtrSscav(i,k,it)=scavtrac*pdtimeRG/dxpres(i,k) ! td du lessivage dans courants satures
	813	dtrUscav(i,k,it)=uscavtrac*pdtimeRG/dxpres(i,k) ! td courant insat
	814	dtrsat(i,k,it)=(trsptrac+scavtrac)*pdtimeRG/dxpres(i,k) ! td courant sat
	815	dtrcv(i,k,it)=(trsptrac+scavtrac+uscavtrac)*pdtimeRG/dxpres(i,k)!dtrsat(i,k,it)+dtrUscav(i,k,it) td conv
	816	!!!!!!
	817	dtrcvMA(i,k,it)=tdcvMA*pdtimeRG/dxpres(i,k) ! MA tendance convection
	818	ENDDO
	819	ENDDO
	820
	821	! test de conservation du traceur
[5103]	822	!PRINT*,"_____________________________________________________________"
	823	!PRINT*," "
[2147]	824	! conserv=0.
	825	! conservMA=0.
	826	! DO k= klev-1,1,-1
	827	! DO i=1, klon
	828	! conserv=conserv+dtrcv(i,k,it)* &
	829	! (paprs(i,k)-paprs(i,k+1))/RG
	830	! conservMA=conservMA+dtrcvMA(i,k,it)* &
	831	! (paprs(i,k)-paprs(i,k+1))/RG
[5099]	832
[5117]	833	! IF(it.EQ.3) WRITE(*,'(I2,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12,2X,a,e20.12)') k,&
[2147]	834	! 'MA td ',dtrcvMA(i,k,it)*dxpres(i,k)/RG,&
	835	! ' td',dtrcv(i,k,it)*dxpres(i,k)/RG,' conservMA ',conservMA,'conserv ',conserv
	836	!!
	837	! ENDDO
	838	! ENDDO
[5160]	839	! IF(it.EQ.3) PRINT *,'it',it,'conserv ',conserv,'conservMA ',conservMA
[2147]	840
	841	END SUBROUTINE cvltr_spl

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

Download in other formats:

Original Format