Context Navigation

calltherm.F90 @ 996

Last change on this file since 996 was 987, checked in by Laurent Fairhead, 16 years ago
Du nettoyage sur le parallelisme, inclusion de nouvelles interfaces pour OPA9 et ORCHIDEE YM LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 12.3 KB

Rev	Line
[878]	1	!
	2	! $Header$
	3	!
	4	subroutine calltherm(dtime &
	5	& ,pplay,paprs,pphi,weak_inversion &
	6	& ,u_seri,v_seri,t_seri,q_seri,zqsat,debut &
	7	& ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
[973]	8	& ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth, &
[927]	9	& ratqsdiff,zqsatth,Ale_bl,Alp_bl,lalim_conv,wght_th, &
	10	& zmax0,f0)
[878]	11
[940]	12	USE dimphy
[878]	13	implicit none
	14	#include "dimensions.h"
[940]	15	!#include "dimphy.h"
[878]	16	#include "thermcell.h"
[938]	17	#include "iniprint.h"
[878]	18
	19	! A inclure eventuellement dans les fichiers de configuration
[973]	20	data r_aspect_thermals,l_mix_thermals/2.,30./
	21	data w2di_thermals/1/
[878]	22
[973]	23	!IM 140508
	24	INTEGER itap
[878]	25	REAL dtime
	26	LOGICAL debut
[973]	27	LOGICAL logexpr0, logexpr2(klon,klev), logexpr1(klon)
	28	REAL fact(klon)
	29	INTEGER nbptspb
	30
[878]	31	REAL u_seri(klon,klev),v_seri(klon,klev)
	32	REAL t_seri(klon,klev),q_seri(klon,klev),qmemoire(klon,klev)
	33	REAL weak_inversion(klon)
	34	REAL paprs(klon,klev+1)
	35	REAL pplay(klon,klev)
	36	REAL pphi(klon,klev)
	37	real zlev(klon,klev+1)
[879]	38	!test: on sort lentr et a* pour alimenter KE
	39	REAL wght_th(klon,klev)
	40	INTEGER lalim_conv(klon)
[878]	41
	42	!FH Update Thermiques
	43	REAL d_t_ajs(klon,klev), d_q_ajs(klon,klev)
	44	REAL d_u_ajs(klon,klev),d_v_ajs(klon,klev)
[973]	45	real fm_therm(klon,klev+1)
	46	real entr_therm(klon,klev),detr_therm(klon,klev)
[878]	47
	48	!********************************************************
	49	! declarations
[973]	50	real fmc_therm(klon,klev+1),zqasc(klon,klev)
[878]	51	real zqla(klon,klev)
	52	real wmax_sec(klon)
	53	real zmax_sec(klon)
	54	real f_sec(klon)
[973]	55	real detrc_therm(klon,klev)
	56	! FH WARNING : il semble que ces save ne servent a rien
	57	! save fmc_therm, detrc_therm
[878]	58	real clwcon0(klon,klev)
	59	real zqsat(klon,klev)
	60	real zw_sec(klon,klev+1)
	61	integer lmix_sec(klon)
	62	integer lmax(klon)
	63	real ratqscth(klon,klev)
	64	real ratqsdiff(klon,klev)
	65	real zqsatth(klon,klev)
[879]	66	!nouvelles variables pour la convection
	67	real Ale_bl(klon)
	68	real Alp_bl(klon)
	69	real Ale(klon)
	70	real Alp(klon)
	71	!RC
[927]	72	!on garde le zmax du pas de temps precedent
	73	real zmax0(klon), f0(klon)
[878]	74	!********************************************************
	75
	76
	77	! variables locales
	78	REAL d_t_the(klon,klev), d_q_the(klon,klev)
	79	REAL d_u_the(klon,klev),d_v_the(klon,klev)
	80	!
[973]	81	real zfm_therm(klon,klev+1),zdt
	82	real zentr_therm(klon,klev),zdetr_therm(klon,klev)
	83	! FH A VERIFIER : SAVE INUTILES
[940]	84	! save zentr_therm,zfm_therm
[973]	85
[878]	86	integer i,k
[973]	87	logical, save :: first=.true.
[987]	88	!$OMP THREADPRIVATE(first)
[878]	89	!********************************************************
[973]	90	if (first) then
	91	itap=0
	92	first=.false.
	93	endif
[878]	94
[973]	95	! Incrementer le compteur de la physique
	96	itap = itap + 1
	97
[878]	98	! Modele du thermique
	99	! ===================
	100	! print*,'thermiques: WARNING on passe t au lieu de t_seri'
	101
	102
[973]	103	! On prend comme valeur initiale des thermiques la valeur du pas
	104	! de temps precedent
	105	zfm_therm(:,:)=fm_therm(:,:)
	106	zdetr_therm(:,:)=detr_therm(:,:)
	107	zentr_therm(:,:)=entr_therm(:,:)
	108
	109	! On reinitialise les flux de masse a zero pour le cumul en
	110	! cas de splitting
[878]	111	fm_therm(:,:)=0.
	112	entr_therm(:,:)=0.
[973]	113	detr_therm(:,:)=0.
	114
[879]	115	Ale_bl(:)=0.
	116	Alp_bl(:)=0.
[938]	117	if (prt_level.ge.10) then
	118	print*,'thermV4 nsplit: ',nsplit_thermals,' weak_inversion'
	119	endif
[878]	120
	121	! tests sur les valeurs negatives de l'eau
[973]	122	logexpr0=prt_level.ge.10
	123	nbptspb=0
[878]	124	do k=1,klev
	125	do i=1,klon
[973]	126	logexpr2(i,k)=.not.q_seri(i,k).ge.0.
	127	if (logexpr2(i,k)) then
[938]	128	q_seri(i,k)=1.e-15
[973]	129	nbptspb=nbptspb+1
[878]	130	endif
[973]	131	! if (logexpr0) &
	132	! & print*,'WARN eau<0 avant therm i=',i,' k=',k &
	133	! & ,' dq,q',d_q_the(i,k),q_seri(i,k)
[878]	134	enddo
	135	enddo
[973]	136	if(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	137
	138	zdt=dtime/float(nsplit_thermals)
	139	do isplit=1,nsplit_thermals
	140
	141	if (iflag_thermals.eq.1) then
	142	CALL thermcell_2002(klon,klev,zdt &
	143	& ,pplay,paprs,pphi &
	144	& ,u_seri,v_seri,t_seri,q_seri &
	145	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	146	& ,zfm_therm,zentr_therm &
	147	& ,r_aspect_thermals,30.,w2di_thermals &
[973]	148	& ,tau_thermals,3)
[878]	149	else if (iflag_thermals.eq.2) then
	150	CALL thermcell_sec(klon,klev,zdt &
	151	& ,pplay,paprs,pphi,zlev &
	152	& ,u_seri,v_seri,t_seri,q_seri &
	153	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	154	& ,zfm_therm,zentr_therm &
	155	& ,r_aspect_thermals,30.,w2di_thermals &
[973]	156	& ,tau_thermals,3)
[878]	157	else if (iflag_thermals.eq.3) then
	158	CALL thermcell(klon,klev,zdt &
	159	& ,pplay,paprs,pphi &
	160	& ,u_seri,v_seri,t_seri,q_seri &
	161	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	162	& ,zfm_therm,zentr_therm &
	163	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	164	& ,tau_thermals,3)
[878]	165	else if (iflag_thermals.eq.10) then
	166	CALL thermcell_eau(klon,klev,zdt &
	167	& ,pplay,paprs,pphi &
	168	& ,u_seri,v_seri,t_seri,q_seri &
	169	& ,d_u_the,d_v_the,d_t_the,d_q_the &
	170	& ,zfm_therm,zentr_therm &
	171	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	172	& ,tau_thermals,3)
[878]	173	else if (iflag_thermals.eq.11) then
	174	stop'cas non prevu dans calltherm'
	175	! CALL thermcell_pluie(klon,klev,zdt &
	176	! & ,pplay,paprs,pphi,zlev &
	177	! & ,u_seri,v_seri,t_seri,q_seri &
	178	! & ,d_u_the,d_v_the,d_t_the,d_q_the &
	179	! & ,zfm_therm,zentr_therm,zqla &
	180	! & ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	181	! & ,tau_thermals,3)
[878]	182	else if (iflag_thermals.eq.12) then
	183	CALL calcul_sec(klon,klev,zdt &
	184	& ,pplay,paprs,pphi,zlev &
	185	& ,u_seri,v_seri,t_seri,q_seri &
	186	& ,zmax_sec,wmax_sec,zw_sec,lmix_sec &
	187	& ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	188	& ,tau_thermals)
[878]	189	! CALL calcul_sec_entr(klon,klev,zdt
	190	! s ,pplay,paprs,pphi,zlev,debut
	191	! s ,u_seri,v_seri,t_seri,q_seri
	192	! s ,zmax_sec,wmax_sec,zw_sec,lmix_sec
	193	! s ,r_aspect_thermals,l_mix_thermals,w2di_thermals
[973]	194	! s ,tau_thermals,3)
[878]	195	! CALL thermcell_pluie_detr(klon,klev,zdt &
	196	! & ,pplay,paprs,pphi,zlev,debut &
	197	! & ,u_seri,v_seri,t_seri,q_seri &
	198	! & ,d_u_the,d_v_the,d_t_the,d_q_the &
	199	! & ,zfm_therm,zentr_therm,zqla,lmax &
	200	! & ,zmax_sec,wmax_sec,zw_sec,lmix_sec &
	201	! & ,ratqscth,ratqsdiff,zqsatth &
	202	! & ,r_aspect_thermals,l_mix_thermals,w2di_thermals &
[973]	203	! & ,tau_thermals)
[878]	204	else if (iflag_thermals.ge.13) then
[973]	205	CALL thermcell_main(itap,klon,klev,zdt &
[878]	206	& ,pplay,paprs,pphi,debut &
	207	& ,u_seri,v_seri,t_seri,q_seri &
	208	& ,d_u_the,d_v_the,d_t_the,d_q_the &
[973]	209	& ,zfm_therm,zentr_therm,zdetr_therm,zqla,lmax &
[878]	210	& ,ratqscth,ratqsdiff,zqsatth &
[973]	211	& ,r_aspect_thermals,l_mix_thermals &
	212	& ,tau_thermals,Ale,Alp,lalim_conv,wght_th &
[927]	213	& ,zmax0,f0)
[878]	214	endif
	215
	216
[973]	217	fact(:)=0.
[878]	218	DO i=1,klon
[973]	219	logexpr1(i)=iflag_thermals.lt.14.or.weak_inversion(i).gt.0.5
	220	IF(logexpr1(i)) fact(i)=1./float(nsplit_thermals)
	221	ENDDO
[878]	222
[973]	223	DO k=1,klev
[878]	224	! transformation de la derivee en tendance
[973]	225	d_t_the(:,k)=d_t_the(:,k)dtimefact(:)
	226	d_u_the(:,k)=d_u_the(:,k)dtimefact(:)
	227	d_v_the(:,k)=d_v_the(:,k)dtimefact(:)
	228	d_q_the(:,k)=d_q_the(:,k)dtimefact(:)
	229	fm_therm(:,k)=fm_therm(:,k) &
	230	& +zfm_therm(:,k)*fact(:)
	231	entr_therm(:,k)=entr_therm(:,k) &
	232	& +zentr_therm(:,k)*fact(:)
	233	ENDDO
	234	fm_therm(:,klev+1)=0.
[878]	235
	236
	237
	238	! accumulation de la tendance
[973]	239	d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
	240	d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
	241	d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
	242	d_q_ajs(:,:)=d_q_ajs(:,:)+d_q_the(:,:)
[878]	243
	244	! incrementation des variables meteo
[973]	245	t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
	246	u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
	247	v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
	248	qmemoire(:,:)=q_seri(:,:)
	249	q_seri(:,:) = q_seri(:,:) + d_q_the(:,:)
[878]	250
[879]	251	DO i=1,klon
[973]	252	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]	253	fm_therm(i,klev+1)=0.
	254	Ale_bl(i)=Ale_bl(i)+Ale(i)/float(nsplit_thermals)
	255	! write(22,*)'ALE CALLTHERM',Ale_bl(i),Ale(i)
	256	Alp_bl(i)=Alp_bl(i)+Alp(i)/float(nsplit_thermals)
	257	! write(23,*)'ALP CALLTHERM',Alp_bl(i),Alp(i)
	258	ENDDO
	259
[973]	260	!IM 060508 marche pas comme cela !!! enddo ! isplit
	261
[878]	262	! tests sur les valeurs negatives de l'eau
[973]	263	nbptspb=0
[878]	264	DO k = 1, klev
	265	DO i = 1, klon
[973]	266	logexpr2(i,k)=.not.q_seri(i,k).ge.0.
	267	if (logexpr2(i,k)) then
	268	q_seri(i,k)=1.e-15
	269	nbptspb=nbptspb+1
	270	! if (prt_level.ge.10) then
	271	! print*,'WARN eau<0 apres therm i=',i,' k=',k &
	272	! & ,' dq,q',d_q_the(i,k),q_seri(i,k), &
	273	! & 'fm=',zfm_therm(i,k),'entr=',entr_therm(i,k)
[938]	274	endif
[878]	275	! stop
	276	ENDDO
	277	ENDDO
[973]	278	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	279	! tests sur les valeurs de la temperature
[973]	280	nbptspb=0
[878]	281	DO k = 1, klev
	282	DO i = 1, klon
[973]	283	logexpr2(i,k)=t_seri(i,k).lt.50..or.t_seri(i,k).gt.370.
	284	if (logexpr2(i,k)) nbptspb=nbptspb+1
	285	! if ((t_seri(i,k).lt.50.) .or. &
	286	! & (t_seri(i,k).gt.370.)) then
	287	! print*,'WARN temp apres therm i=',i,' k=',k &
	288	! & ,' t_seri',t_seri(i,k)
[878]	289	! CALL abort
[973]	290	! endif
[878]	291	ENDDO
	292	ENDDO
[973]	293	IF(nbptspb.GT.0) print*,'Number of points with q_seri(i,k)<=0 ',nbptspb
[878]	294	enddo ! isplit
	295
	296	!
	297	!***************************************************************
	298	! calcul du flux ascencant conservatif
	299	! print*,'<<<<calcul flux ascendant conservatif'
	300
	301	fmc_therm=0.
	302	do k=1,klev
	303	do i=1,klon
	304	if (entr_therm(i,k).gt.0.) then
	305	fmc_therm(i,k+1)=fmc_therm(i,k)+entr_therm(i,k)
	306	else
	307	fmc_therm(i,k+1)=fmc_therm(i,k)
	308	endif
	309	detrc_therm(i,k)=(fmc_therm(i,k+1)-fm_therm(i,k+1)) &
	310	& -(fmc_therm(i,k)-fm_therm(i,k))
	311	enddo
	312	enddo
	313
	314
	315	!****************************************************************
	316	! calcul de l'humidite dans l'ascendance
	317	! print*,'<<<<calcul de lhumidite dans thermique'
	318	!CR:on ne le calcule que pour le cas sec
	319	if (iflag_thermals.le.11) then
	320	do i=1,klon
	321	zqasc(i,1)=q_seri(i,1)
	322	do k=2,klev
	323	if (fmc_therm(i,k+1).gt.1.e-6) then
	324	zqasc(i,k)=(fmc_therm(i,k)*zqasc(i,k-1) &
	325	& +entr_therm(i,k)*q_seri(i,k))/fmc_therm(i,k+1)
	326	!CR:test on asseche le thermique
	327	! zqasc(i,k)=zqasc(i,k)/2.
	328	! else
	329	! zqasc(i,k)=q_seri(i,k)
	330	endif
	331	enddo
	332	enddo
	333
	334
	335	! calcul de l'eau condensee dans l'ascendance
	336	! print*,'<<<<calcul de leau condensee dans thermique'
	337	do i=1,klon
	338	do k=1,klev
	339	clwcon0(i,k)=zqasc(i,k)-zqsat(i,k)
	340	if (clwcon0(i,k).lt.0. .or. &
	341	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
	342	clwcon0(i,k)=0.
	343	endif
	344	enddo
	345	enddo
	346	else
	347	do i=1,klon
	348	do k=1,klev
	349	clwcon0(i,k)=zqla(i,k)
	350	if (clwcon0(i,k).lt.0. .or. &
	351	& (fm_therm(i,k+1)+detrc_therm(i,k)).lt.1.e-6) then
	352	clwcon0(i,k)=0.
	353	endif
	354	enddo
	355	enddo
	356	endif
	357	!*******************************************************************
	358
	359
	360	return
	361
	362	end

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

Download in other formats:

Original Format