Context Navigation

calltherm.F90 @ 5450

Last change on this file since 5450 was 5390, checked in by yann meurdesoif, 3 weeks ago

Remove UTF8 character that inihibit fortran parsing with GPU morphosis
Add missing END SUBROUTINE instead of simple END, that inhibit correct parsing with regulat expression parser (quick and dirty parsing)

Property copyright set to
Name of program: LMDZ
Creation date: 1984
Version: LMDZ5
License: CeCILL version 2
Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
See the license file in the root directory
Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 18.6 KB

Rev	Line
[878]	1	!
[1403]	2	! $Id: calltherm.F90 5390 2024-12-05 16:09:25Z aborella $
[878]	3	!
	4	subroutine calltherm(dtime &
	5	& ,pplay,paprs,pphi,weak_inversion &
[4690]	6	& ,u_seri_,v_seri_,t_seri_,q_seri_,t_env,q_env,zqsat,debut &
[878]	7	& ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
[4837]	8	& ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth,&
[4089]	9	& ratqsdiff,zqsatth,ale_bl,alp_bl,lalim_conv,wght_th, &
[1638]	10	& zmax0,f0,zw2,fraca,ztv,zpspsk,ztla,zthl &
	11	!!! nrlmd le 10/04/2012
	12	& ,pbl_tke,pctsrf,omega,airephy &
	13	& ,zlcl_th,fraca0,w0,w_conv,therm_tke_max0,env_tke_max0 &
[4843]	14	& ,n2,s2,strig,zcong,ale_bl_stat &
[1638]	15	& ,therm_tke_max,env_tke_max &
	16	& ,alp_bl_det,alp_bl_fluct_m,alp_bl_fluct_tke &
	17	& ,alp_bl_conv,alp_bl_stat &
	18	!!! fin nrlmd le 10/04/2012
[4089]	19	& ,zqla,ztva &
	20	#ifdef ISO
	21	& ,xt_seri,d_xt_ajs &
	22	#ifdef DIAGISO
	23	& ,q_the,xt_the &
	24	#endif
	25	#endif
	26	& )
[878]	27
[940]	28	USE dimphy
[1785]	29	USE indice_sol_mod
[2311]	30	USE print_control_mod, ONLY: prt_level,lunout
[4590]	31	USE lmdz_thermcell_alp, ONLY: thermcell_alp
	32	USE lmdz_thermcell_main, ONLY: thermcell_main
	33	USE lmdz_thermcell_old, ONLY: thermcell, thermcell_2002, thermcell_eau, calcul_sec, thermcell_sec
[4089]	34	#ifdef ISO
[4143]	35	use infotrac_phy, ONLY: ntiso
[4089]	36	#ifdef ISOVERIF
	37	USE isotopes_mod, ONLY: iso_eau,iso_HDO
	38	USE isotopes_verif_mod, ONLY: iso_verif_aberrant_enc_vect2D, &
	39	iso_verif_egalite_vect2D
	40	#endif
	41	#endif
[5282]	42	USE clesphys_mod_h
[5291]	43	USE thermcell_old_mod_h, ONLY: r_aspect_thermals, l_mix_thermals, w2di_thermals
[878]	44	implicit none
	45
[1638]	46
[973]	47	!IM 140508
[1776]	48	INTEGER, SAVE :: itap
	49	!$OMP THREADPRIVATE(itap)
[878]	50	REAL dtime
	51	LOGICAL debut
[973]	52	LOGICAL logexpr0, logexpr2(klon,klev), logexpr1(klon)
	53	REAL fact(klon)
	54	INTEGER nbptspb
	55
[4678]	56	REAL, DIMENSION(klon,klev), INTENT(IN) :: u_seri_,v_seri_
[4837]	57	REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri_,q_seri_
	58	REAL, DIMENSION(klon,klev), INTENT(IN) :: t_env,q_env
[4678]	59	REAL, DIMENSION(klon,klev) :: u_seri,v_seri
	60	REAL, DIMENSION(klon,klev) :: t_seri,q_seri
	61	REAL, DIMENSION(klon,klev) :: qmemoire
[878]	62	REAL weak_inversion(klon)
	63	REAL paprs(klon,klev+1)
	64	REAL pplay(klon,klev)
	65	REAL pphi(klon,klev)
	66	real zlev(klon,klev+1)
[879]	67	!test: on sort lentr et a* pour alimenter KE
	68	REAL wght_th(klon,klev)
	69	INTEGER lalim_conv(klon)
[1026]	70	REAL zw2(klon,klev+1),fraca(klon,klev+1)
[878]	71
	72	!FH Update Thermiques
	73	REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
	74	REAL d_u_ajs(klon,klev),d_v_ajs(klon,klev)
[973]	75	real fm_therm(klon,klev+1)
	76	real entr_therm(klon,klev),detr_therm(klon,klev)
[878]	77
	78	!********************************************************
	79	! declarations
[1403]	80	LOGICAL flag_bidouille_stratocu
[973]	81	real fmc_therm(klon,klev+1),zqasc(klon,klev)
[878]	82	real zqla(klon,klev)
[1790]	83	real ztv(klon,klev),ztva(klon,klev)
[1403]	84	real zpspsk(klon,klev)
	85	real ztla(klon,klev)
	86	real zthl(klon,klev)
[878]	87	real wmax_sec(klon)
[4843]	88	real zcong(klon)
[878]	89	real zmax_sec(klon)
	90	real f_sec(klon)
[973]	91	real detrc_therm(klon,klev)
	92	! FH WARNING : il semble que ces save ne servent a rien
	93	! save fmc_therm, detrc_therm
[878]	94	real clwcon0(klon,klev)
	95	real zqsat(klon,klev)
	96	real zw_sec(klon,klev+1)
	97	integer lmix_sec(klon)
	98	integer lmax(klon)
	99	real ratqscth(klon,klev)
	100	real ratqsdiff(klon,klev)
	101	real zqsatth(klon,klev)
[879]	102	!nouvelles variables pour la convection
[4089]	103	real ale_bl(klon)
	104	real alp_bl(klon)
	105	real ale(klon)
	106	real alp(klon)
[879]	107	!RC
[927]	108	!on garde le zmax du pas de temps precedent
	109	real zmax0(klon), f0(klon)
[1638]	110
	111	!!! nrlmd le 10/04/2012
	112	real pbl_tke(klon,klev+1,nbsrf)
	113	real pctsrf(klon,nbsrf)
	114	real omega(klon,klev)
	115	real airephy(klon)
	116	real zlcl_th(klon),fraca0(klon),w0(klon),w_conv(klon)
	117	real therm_tke_max0(klon),env_tke_max0(klon)
[4843]	118	real n2(klon),s2(klon),strig(klon)
[1638]	119	real ale_bl_stat(klon)
	120	real therm_tke_max(klon,klev),env_tke_max(klon,klev)
	121	real alp_bl_det(klon),alp_bl_fluct_m(klon),alp_bl_fluct_tke(klon),alp_bl_conv(klon),alp_bl_stat(klon)
	122	!!! fin nrlmd le 10/04/2012
	123
[878]	124	!********************************************************
	125
[4089]	126	real, dimension(klon) :: pcon
	127	real, dimension(klon,klev) :: rhobarz,wth3
	128	integer,dimension(klon) :: lalim
	129	real, dimension(klon,klev+1) :: fm
	130	real, dimension(klon,klev) :: alim_star
	131	real, dimension(klon) :: zmax
[878]	132
[4089]	133
	134
	135
[878]	136	! variables locales
	137	REAL d_t_the(klon,klev), d_q_the(klon,klev)
	138	REAL d_u_the(klon,klev),d_v_the(klon,klev)
	139	!
[973]	140	real zfm_therm(klon,klev+1),zdt
	141	real zentr_therm(klon,klev),zdetr_therm(klon,klev)
	142	! FH A VERIFIER : SAVE INUTILES
[940]	143	! save zentr_therm,zfm_therm
[973]	144
[1403]	145	character (len=20) :: modname='calltherm'
	146	character (len=80) :: abort_message
	147
[4089]	148	integer i,k,isplit
[973]	149	logical, save :: first=.true.
[4089]	150	logical :: new_thermcell
	151
	152	#ifdef ISO
[4143]	153	REAL xt_seri(ntiso,klon,klev),xtmemoire(ntiso,klon,klev)
	154	REAL d_xt_ajs(ntiso,klon,klev)
	155	real d_xt_the(ntiso,klon,klev)
[4089]	156	#ifdef DIAGISO
	157	real q_the(klon,klev)
[4143]	158	real xt_the(ntiso,klon,klev)
[4089]	159	#endif
	160	real qprec(klon,klev)
	161	integer ixt
	162	#endif
	163
	164
[987]	165	!$OMP THREADPRIVATE(first)
[878]	166	!********************************************************
[973]	167	if (first) then
	168	itap=0
	169	first=.false.
	170	endif
[878]	171
[4678]	172	u_seri(:,:)=u_seri_(:,:)
	173	v_seri(:,:)=v_seri_(:,:)
	174	t_seri(:,:)=t_seri_(:,:)
	175	q_seri(:,:)=q_seri_(:,:)
	176
[973]	177	! Incrementer le compteur de la physique
	178	itap = itap + 1
	179
[878]	180	! Modele du thermique
	181	! ===================
	182	! print*,'thermiques: WARNING on passe t au lieu de t_seri'
	183
	184
[973]	185	! On prend comme valeur initiale des thermiques la valeur du pas
	186	! de temps precedent
	187	zfm_therm(:,:)=fm_therm(:,:)
	188	zdetr_therm(:,:)=detr_therm(:,:)
	189	zentr_therm(:,:)=entr_therm(:,:)
	190
	191	! On reinitialise les flux de masse a zero pour le cumul en
	192	! cas de splitting
[878]	193	fm_therm(:,:)=0.
	194	entr_therm(:,:)=0.
[973]	195	detr_therm(:,:)=0.
	196
[4089]	197	ale_bl(:)=0.
	198	alp_bl(:)=0.
[938]	199	if (prt_level.ge.10) then
	200	print*,'thermV4 nsplit: ',nsplit_thermals,' weak_inversion'
	201	endif
[878]	202
	203	! tests sur les valeurs negatives de l'eau
[973]	204	logexpr0=prt_level.ge.10
	205	nbptspb=0
[878]	206	do k=1,klev
	207	do i=1,klon
[1146]	208	! Attention teste abderr 19-03-09
	209	! logexpr2(i,k)=.not.q_seri(i,k).ge.0.
	210	logexpr2(i,k)=.not.q_seri(i,k).ge.1.e-15
[973]	211	if (logexpr2(i,k)) then
[4089]	212	#ifdef ISO
	213	qprec(i,k)=q_seri(i,k)
	214	#endif
[938]	215	q_seri(i,k)=1.e-15
[973]	216	nbptspb=nbptspb+1
[4089]	217	#ifdef ISO
[4143]	218	do ixt=1,ntiso
[4089]	219	xt_seri(ixt,i,k)=1.e-15*(xt_seri(ixt,i,k)/qprec(i,k))
	220	! xt_seri(ixt,i,k)=1.e-15*(Rdefault(index_iso(ixt)))
	221	enddo
	222	#endif
[878]	223	endif
[973]	224	! if (logexpr0) &
	225	! & print*,'WARN eau<0 avant therm i=',i,' k=',k &
	226	! & ,' dq,q',d_q_the(i,k),q_seri(i,k)
[878]	227	enddo
	228	enddo
[973]	229	if(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	230
[4089]	231
	232	new_thermcell=iflag_thermals>=15.and.iflag_thermals<=18
	233	#ifdef ISO
	234	if (.not.new_thermcell) then
[5217]	235	CALL abort_physic('calltherm 234','isos pas prevus ici',1)
[4089]	236	endif
	237	#ifdef ISOVERIF
	238	if (iso_eau.gt.0) then
	239	call iso_verif_egalite_vect2D( &
	240	& xt_seri,q_seri, &
[4143]	241	& 'calltherm 174',ntiso,klon,klev)
[4089]	242	endif !if (iso_eau.gt.0) then
	243	#endif
	244	#endif
[1403]	245	zdt=dtime/REAL(nsplit_thermals)
[4089]	246
	247
[878]	248	do isplit=1,nsplit_thermals
	249
[1943]	250	if (iflag_thermals>=1000) then
	251	CALL thermcell_2002(klon,klev,zdt,iflag_thermals &
[878]	252	& ,pplay,paprs,pphi &
	253	& ,u_seri,v_seri,t_seri,q_seri &
	254	& ,d_u_the,d_v_the,d_t_the,d_q_the &
[1943]	255	& ,zfm_therm,zentr_therm,fraca,zw2 &
[878]	256	& ,r_aspect_thermals,30.,w2di_thermals &
[1517]	257	& ,tau_thermals)
[878]	258	else if (iflag_thermals.eq.2) then
	259	CALL thermcell_sec(klon,klev,zdt &
	260	& ,pplay,paprs,pphi,zlev &
	261	& ,u_seri,v_seri,t_seri,q_seri &
	262	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	263	& ,zfm_therm,zentr_therm &
	264	& ,r_aspect_thermals,30.,w2di_thermals &
[1517]	265	& ,tau_thermals)
[878]	266	else if (iflag_thermals.eq.3) then
	267	CALL thermcell(klon,klev,zdt &
	268	& ,pplay,paprs,pphi &
	269	& ,u_seri,v_seri,t_seri,q_seri &
	270	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	271	& ,zfm_therm,zentr_therm &
	272	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[1517]	273	& ,tau_thermals)
[878]	274	else if (iflag_thermals.eq.10) then
	275	CALL thermcell_eau(klon,klev,zdt &
	276	& ,pplay,paprs,pphi &
	277	& ,u_seri,v_seri,t_seri,q_seri &
	278	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	279	& ,zfm_therm,zentr_therm &
	280	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[1517]	281	& ,tau_thermals)
[878]	282	else if (iflag_thermals.eq.11) then
[1403]	283	abort_message = 'cas non prevu dans calltherm'
[2311]	284	CALL abort_physic (modname,abort_message,1)
[878]	285	else if (iflag_thermals.eq.12) then
	286	CALL calcul_sec(klon,klev,zdt &
	287	& ,pplay,paprs,pphi,zlev &
	288	& ,u_seri,v_seri,t_seri,q_seri &
	289	& ,zmax_sec,wmax_sec,zw_sec,lmix_sec &
	290	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	291	& ,tau_thermals)
[1403]	292	else if (iflag_thermals==13.or.iflag_thermals==14) then
[4089]	293	abort_message = 'thermcellV0_main enleve svn>2084'
	294	CALL abort_physic (modname,abort_message,1)
	295	else if (new_thermcell) then
[1403]	296	CALL thermcell_main(itap,klon,klev,zdt &
	297	& ,pplay,paprs,pphi,debut &
[4690]	298	& ,u_seri,v_seri,t_seri,q_seri,t_env,q_env &
[1403]	299	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	300	& ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax &
	301	& ,ratqscth,ratqsdiff,zqsatth &
	302	& ,zmax0,f0,zw2,fraca,ztv,zpspsk &
[4089]	303	& ,ztla,zthl,ztva &
[4843]	304	& ,pcon,rhobarz,wth3,wmax_sec,lalim,fm,alim_star,zmax,zcong &
[4089]	305	#ifdef ISO
	306	& ,xt_seri,d_xt_the &
	307	#endif
	308	& )
	309
	310	CALL thermcell_alp(klon,klev,zdt & ! in
	311	& ,pplay,paprs & ! in
	312	& ,zfm_therm,zentr_therm,lmax & ! in
	313	& ,pbl_tke,pctsrf,omega,airephy & ! in
	314	& ,zw2,fraca & ! in
	315	& ,pcon,rhobarz,wth3,wmax_sec,lalim,fm,alim_star,zmax & ! in
[4843]	316	& ,zcong,ale,alp,lalim_conv,wght_th & ! out
[4089]	317	& ,zlcl_th,fraca0,w0,w_conv,therm_tke_max0,env_tke_max0 &! out
[4843]	318	& ,n2,s2,strig,ale_bl_stat & ! out
[4089]	319	& ,therm_tke_max,env_tke_max & ! out
	320	& ,alp_bl_det,alp_bl_fluct_m,alp_bl_fluct_tke & ! out
	321	& ,alp_bl_conv,alp_bl_stat & ! out
	322	& )
	323
[1403]	324	if (prt_level.gt.10) write(lunout,*)'Apres thermcell_main OK'
	325	else
	326	abort_message = 'Cas des thermiques non prevu'
[2311]	327	CALL abort_physic (modname,abort_message,1)
[878]	328	endif
	329
[1638]	330	! Attention : les noms sont contre intuitif.
	331	! flag_bidouille_stratocu est .true. si on ne fait pas de bidouille.
	332	! Il aurait mieux valu avoir un nobidouille_stratocu
	333	! Et pour simplifier :
	334	! nobidouille_stratocu=.not.(iflag_thermals==13.or.iflag_thermals=15)
	335	! Ce serait bien de changer, mai en prenant le temps de vérifier que ca
	336	! fait bien ce qu'on croit.
[878]	337
[1738]	338	flag_bidouille_stratocu=iflag_thermals<=12.or.iflag_thermals==14.or.iflag_thermals==16.or.iflag_thermals==18
[1638]	339
[1943]	340	! Calcul a posteriori du niveau max des thermiques pour les schémas qui
	341	! ne la sortent pas.
	342	if (iflag_thermals<=12.or.iflag_thermals>=1000) then
	343	lmax(:)=1
[1638]	344	do k=1,klev-1
	345	zdetr_therm(:,k)=zentr_therm(:,k)+zfm_therm(:,k)-zfm_therm(:,k+1)
	346	enddo
[1943]	347	do k=1,klev-1
	348	do i=1,klon
	349	if (zfm_therm(i,k+1)>0.) lmax(i)=k
	350	enddo
	351	enddo
[1638]	352	endif
	353
[973]	354	fact(:)=0.
[878]	355	DO i=1,klon
[1403]	356	logexpr1(i)=flag_bidouille_stratocu.or.weak_inversion(i).gt.0.5
	357	IF(logexpr1(i)) fact(i)=1./REAL(nsplit_thermals)
[973]	358	ENDDO
[878]	359
[973]	360	DO k=1,klev
[878]	361	! transformation de la derivee en tendance
[973]	362	d_t_the(:,k)=d_t_the(:,k)dtimefact(:)
	363	d_u_the(:,k)=d_u_the(:,k)dtimefact(:)
	364	d_v_the(:,k)=d_v_the(:,k)dtimefact(:)
	365	d_q_the(:,k)=d_q_the(:,k)dtimefact(:)
	366	fm_therm(:,k)=fm_therm(:,k) &
	367	& +zfm_therm(:,k)*fact(:)
	368	entr_therm(:,k)=entr_therm(:,k) &
	369	& +zentr_therm(:,k)*fact(:)
[1026]	370	detr_therm(:,k)=detr_therm(:,k) &
	371	& +zdetr_therm(:,k)*fact(:)
[4089]	372	#ifdef ISO
[4143]	373	do ixt=1,ntiso
[4089]	374	d_xt_the(ixt,:,k)=d_xt_the(ixt,:,k)dtimefact(:)
	375	enddo
	376	#endif
[973]	377	ENDDO
	378	fm_therm(:,klev+1)=0.
[878]	379
	380
	381
	382	! accumulation de la tendance
[973]	383	d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
	384	d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
	385	d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
	386	d_q_ajs(:,:)=d_q_ajs(:,:)+d_q_the(:,:)
[4089]	387	#ifdef ISO
	388	d_xt_ajs(:,:,:)=d_xt_ajs(:,:,:)+d_xt_the(:,:,:)
	389	#endif
[878]	390
	391	! incrementation des variables meteo
[973]	392	t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
	393	u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
	394	v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
	395	qmemoire(:,:)=q_seri(:,:)
	396	q_seri(:,:) = q_seri(:,:) + d_q_the(:,:)
[4089]	397	#ifdef ISO
	398	xtmemoire(:,:,:)=xt_seri(:,:,:)
	399	xt_seri(:,:,:) = xt_seri(:,:,:) + d_xt_the(:,:,:)
	400	#ifdef ISOVERIF
	401	! write(,) 'calltherm 350 tmp: ajout d_xt_the'
	402	if (iso_HDO.gt.0) then
	403	! i=479
	404	! k=4
	405	! write(,) 'xt_seri(iso_hdo,i,k),q_seri(i,k)=', &
	406	! & xt_seri(iso_hdo,i,k),q_seri(i,k)
	407	! write(,) 'd_xt_the(iso_hdo,i,k),d_q_the(i,k)=', &
	408	! & d_xt_the(iso_hdo,i,k),d_q_the(i,k)
	409	call iso_verif_aberrant_enc_vect2D( &
	410	& xt_seri,q_seri, &
[4143]	411	& 'calltherm 353, apres ajout d_xt_the',ntiso,klon,klev)
[4089]	412	endif
	413	#endif
	414	#endif
[1403]	415	if (prt_level.gt.10) write(lunout,*)'Apres apres thermcell_main OK'
[878]	416
[879]	417	DO i=1,klon
	418	fm_therm(i,klev+1)=0.
[4089]	419	ale_bl(i)=ale_bl(i)+ale(i)/REAL(nsplit_thermals)
	420	! write(22,*)'ALE CALLTHERM',ale_bl(i),ale(i)
	421	alp_bl(i)=alp_bl(i)+alp(i)/REAL(nsplit_thermals)
	422	! write(23,*)'ALP CALLTHERM',alp_bl(i),alp(i)
	423	if(prt_level.GE.10) print*,'calltherm i alp_bl alp ale_bl ale',i,alp_bl(i),alp(i),ale_bl(i),ale(i)
[879]	424	ENDDO
	425
[973]	426	!IM 060508 marche pas comme cela !!! enddo ! isplit
	427
[878]	428	! tests sur les valeurs negatives de l'eau
[973]	429	nbptspb=0
[878]	430	DO k = 1, klev
	431	DO i = 1, klon
[973]	432	logexpr2(i,k)=.not.q_seri(i,k).ge.0.
	433	if (logexpr2(i,k)) then
	434	q_seri(i,k)=1.e-15
	435	nbptspb=nbptspb+1
[4089]	436	#ifdef ISO
[4143]	437	do ixt=1,ntiso
[4089]	438	xt_seri(ixt,i,k)=1.e-15*(xtmemoire(ixt,i,k)/qmemoire(i,k))
	439	enddo
	440	#endif
[973]	441	! if (prt_level.ge.10) then
	442	! print*,'WARN eau<0 apres therm i=',i,' k=',k &
	443	! & ,' dq,q',d_q_the(i,k),q_seri(i,k), &
	444	! & 'fm=',zfm_therm(i,k),'entr=',entr_therm(i,k)
[938]	445	endif
[878]	446	ENDDO
	447	ENDDO
[4089]	448	#ifdef ISO
	449	#ifdef ISOVERIF
	450	if (iso_HDO.gt.0) then
	451	call iso_verif_aberrant_enc_vect2D( &
	452	& xt_seri,q_seri, &
[4143]	453	& 'calltherm 393, apres bidouille q<0',ntiso,klon,klev)
[4089]	454	endif
	455	#endif
	456	#endif
	457
[973]	458	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	459	! tests sur les valeurs de la temperature
[973]	460	nbptspb=0
[878]	461	DO k = 1, klev
	462	DO i = 1, klon
[973]	463	logexpr2(i,k)=t_seri(i,k).lt.50..or.t_seri(i,k).gt.370.
	464	if (logexpr2(i,k)) nbptspb=nbptspb+1
	465	! if ((t_seri(i,k).lt.50.) .or. &
	466	! & (t_seri(i,k).gt.370.)) then
	467	! print*,'WARN temp apres therm i=',i,' k=',k &
	468	! & ,' t_seri',t_seri(i,k)
[878]	469	! CALL abort
[973]	470	! endif
[878]	471	ENDDO
	472	ENDDO
[973]	473	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	474	enddo ! isplit
	475
	476	!
	477	!***************************************************************
	478	! calcul du flux ascencant conservatif
	479	! print*,'<<<<calcul flux ascendant conservatif'
	480
	481	fmc_therm=0.
	482	do k=1,klev
	483	do i=1,klon
	484	if (entr_therm(i,k).gt.0.) then
	485	fmc_therm(i,k+1)=fmc_therm(i,k)+entr_therm(i,k)
	486	else
	487	fmc_therm(i,k+1)=fmc_therm(i,k)
	488	endif
	489	detrc_therm(i,k)=(fmc_therm(i,k+1)-fm_therm(i,k+1)) &
	490	& -(fmc_therm(i,k)-fm_therm(i,k))
	491	enddo
	492	enddo
	493
	494
	495	!****************************************************************
	496	! calcul de l'humidite dans l'ascendance
	497	! print*,'<<<<calcul de lhumidite dans thermique'
	498	!CR:on ne le calcule que pour le cas sec
	499	if (iflag_thermals.le.11) then
	500	do i=1,klon
	501	zqasc(i,1)=q_seri(i,1)
	502	do k=2,klev
	503	if (fmc_therm(i,k+1).gt.1.e-6) then
	504	zqasc(i,k)=(fmc_therm(i,k)*zqasc(i,k-1) &
	505	& +entr_therm(i,k)*q_seri(i,k))/fmc_therm(i,k+1)
	506	!CR:test on asseche le thermique
	507	! zqasc(i,k)=zqasc(i,k)/2.
	508	! else
	509	! zqasc(i,k)=q_seri(i,k)
	510	endif
	511	enddo
	512	enddo
	513
	514
	515	! calcul de l'eau condensee dans l'ascendance
	516	! print*,'<<<<calcul de leau condensee dans thermique'
	517	do i=1,klon
	518	do k=1,klev
	519	clwcon0(i,k)=zqasc(i,k)-zqsat(i,k)
	520	if (clwcon0(i,k).lt.0. .or. &
	521	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
	522	clwcon0(i,k)=0.
	523	endif
	524	enddo
	525	enddo
	526	else
	527	do i=1,klon
	528	do k=1,klev
	529	clwcon0(i,k)=zqla(i,k)
	530	if (clwcon0(i,k).lt.0. .or. &
	531	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
	532	clwcon0(i,k)=0.
	533	endif
	534	enddo
	535	enddo
	536	endif
	537	!*******************************************************************
	538
	539
[1428]	540	!jyg Protection contre les temperatures nulles
	541	do i=1,klon
	542	do k=1,klev
	543	if (ztla(i,k) .lt. 1.e-10) fraca(i,k) =0.
	544	enddo
	545	enddo
	546
	547
[878]	548	return
	549
[5390]	550	end subroutine calltherm

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

Download in other formats:

Original Format