Context Navigation

calltherm.F90 @ 1574

Last change on this file since 1574 was 1517, checked in by jghattas, 14 years ago

Errors in argument list not corresponding between CALL and SUBROUTINE in some routines. These errors were found with compile option -qextchk at vargas.

Removed option -qextchk from debug compile options at vargas. This option is not usable with mpi codes.

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 13.2 KB

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

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

Download in other formats:

Original Format