Context Navigation

alpale_th.F90 @ 3324

Last change on this file since 3324 was 3209, checked in by jyg, 7 years ago
error in the prvious commit
File size: 11.9 KB

Rev	Line
[2558]	1	SUBROUTINE alpale_th ( dtime, lmax_th, t_seri, cell_area, &
[2513]	2	cin, s2, n2, &
	3	ale_bl_trig, ale_bl_stat, ale_bl, &
[2556]	4	alp_bl, alp_bl_stat, &
[3209]	5	proba_notrig, random_notrig, birth_rate)
[2513]	6
	7	! **************************************************************
	8	! *
	9	! ALPALE_TH *
	10	! *
	11	! *
	12	! written by : Jean-Yves Grandpeix, 11/05/2016 *
	13	! modified by : *
	14	! **************************************************************
	15
	16	USE dimphy
	17	USE ioipsl_getin_p_mod, ONLY : getin_p
	18	USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
	19	!
	20	IMPLICIT NONE
	21
	22	!================================================================
	23	! Auteur(s) : Jean-Yves Grandpeix, 11/05/2016
	24	! Objet : Contribution of the thermal scheme to Ale and Alp
	25	!================================================================
	26
	27	! Input arguments
	28	!----------------
	29	REAL, INTENT(IN) :: dtime
[2558]	30	REAL, DIMENSION(klon), INTENT(IN) :: cell_area
[2513]	31	INTEGER, DIMENSION(klon), INTENT(IN) :: lmax_th
	32	REAL, DIMENSION(klon,klev), INTENT(IN) :: t_seri
	33	REAL, DIMENSION(klon), INTENT(IN) :: ale_bl_stat
	34	REAL, DIMENSION(klon), INTENT(IN) :: cin
	35	REAL, DIMENSION(klon), INTENT(IN) :: s2, n2
	36
	37	REAL, DIMENSION(klon), INTENT(INOUT) :: ale_bl_trig, ale_bl
	38	REAL, DIMENSION(klon), INTENT(INOUT) :: alp_bl
	39	REAL, DIMENSION(klon), INTENT(INOUT) :: alp_bl_stat
[2730]	40	REAL, DIMENSION(klon), INTENT(INOUT) :: proba_notrig
[2513]	41
[2730]	42	REAL, DIMENSION(klon), INTENT(OUT) :: random_notrig
[2556]	43
[3208]	44	REAL, DIMENSION(klon), INTENT(OUT) :: birth_rate
	45
[2513]	46	include "thermcell.h"
	47
	48	! Local variables
	49	!----------------
	50	INTEGER :: i
	51	LOGICAL, SAVE :: first = .TRUE.
	52	REAL, SAVE :: random_notrig_max=1.
[2558]	53	REAL, SAVE :: cv_feed_area
	54	REAL :: birth_number
	55	REAL, DIMENSION(klon) :: ale_bl_ref
[2556]	56	REAL, DIMENSION(klon) :: tau_trig
[2513]	57	!
	58	!$OMP THREADPRIVATE(random_notrig_max)
[2558]	59	!$OMP THREADPRIVATE(cv_feed_area)
[2513]	60	!$OMP THREADPRIVATE(first)
	61	!
[2558]	62	REAL umexp ! expression of (1.-exp(-x))/x valid for all x, especially when x->0
	63	REAL x
	64	umexp(x) = max(sign(1.,x-1.e-3),0.)*(1.-exp(-x))/max(x,1.e-3) + &
[3073]	65	(1.-max(sign(1.,x-1.e-3),0.))(1.-0.5x(1.-x/3.(1.-0.25*x))) !!! correct formula (jyg)
	66	!!! (1.-max(sign(1.,x-1.e-3),0.))(-0.5x(1.-x/3.(1.-0.25*x))) !!! bug introduced by mistake (jyg)
	67	!!! (1.-max(sign(1.,x-1.e-3),0.))(1.-0.5x(1.-x/3.(1.-0.25*x))) !!! initial correct formula (jyg)
[2558]	68	!
	69	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
	70	! JYG, 20160513 : Introduction of the Effective Lifting Power (ELP), which
	71	! takes into account the area (cv_feed_area) covered by thermals contributing
	72	! to each cumulonimbus.
	73	! The use of ELP prevents singularities when the trigger probability tends to
	74	! zero. It is activated by iflag_clos_bl = 3.
	75	! The ELP values are stored in the ALP_bl variable.
	76	!
	77	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
	78	!
[2730]	79	!!
	80	!! The following 3 lines should be commented if one wants to activate the
	81	!! multiplication of no-trigger probabilities between calls to the convection
	82	!! scheme.
	83	!!
	84	do i=1,klon
	85	proba_notrig(i)=1.
	86	enddo
	87	!!
	88	!!
[2558]	89	!---------------------------------------
	90	IF (iflag_clos_bl .LT. 3) THEN
	91	!---------------------------------------
	92	!
	93	! Original code (Nicolas Rochetin)
	94	! --------------------------------
	95
[2513]	96	IF (first) THEN
[2558]	97	random_notrig_max=1.
[2513]	98	CALL getin_p('random_notrig_max',random_notrig_max)
[2515]	99	first=.FALSE.
[2513]	100	ENDIF
	101	!cc nrlmd le 10/04/2012
	102	!-----------Stochastic triggering-----------
	103	if (iflag_trig_bl.ge.1) then
	104	!
	105	IF (prt_level .GE. 10) THEN
[3208]	106	print *,'cin, ale_bl_stat, alp_bl, alp_bl_stat ', &
	107	cin, ale_bl_stat, alp_bl, alp_bl_stat
[2513]	108	ENDIF
	109
	110
	111	!----Initialisations
	112	do i=1,klon
[2730]	113	!!jyg proba_notrig(i)=1.
[2513]	114	random_notrig(i)=1e6ale_bl_stat(i)-int(1e6ale_bl_stat(i))
	115	if ( random_notrig(i) > random_notrig_max ) random_notrig(i)=0.
	116	if ( ale_bl_trig(i) .lt. abs(cin(i))+1.e-10 ) then
	117	tau_trig(i)=tau_trig_shallow
	118	else
	119	tau_trig(i)=tau_trig_deep
	120	endif
	121	enddo
	122	!
	123	IF (prt_level .GE. 10) THEN
	124	print *,'random_notrig, tau_trig ', &
	125	random_notrig, tau_trig
	126	print *,'s_trig,s2,n2 ', &
	127	s_trig,s2,n2
	128	ENDIF
	129
	130	!Option pour re-activer l'ancien calcul de Ale_bl (iflag_trig_bl=2)
	131	IF (iflag_trig_bl.eq.1) then
	132
	133	!----Tirage al\'eatoire et calcul de ale_bl_trig
	134	do i=1,klon
	135	if ( (ale_bl_stat(i) .gt. abs(cin(i))+1.e-10) ) then
[2730]	136	proba_notrig(i)=proba_notrig(i)* &
	137	(1.-exp(-s_trig/s2(i)))*(n2(i)dtime/tau_trig(i))
[2513]	138	! print *, 'proba_notrig(i) ',proba_notrig(i)
	139	if (random_notrig(i) .ge. proba_notrig(i)) then
	140	ale_bl_trig(i)=ale_bl_stat(i)
	141	else
	142	ale_bl_trig(i)=0.
	143	endif
[3208]	144	birth_rate(i) = n2(i)exp(-s_trig/s2(i))/(tau_trig(i)cell_area(i))
	145	!!! birth_rate(i) = max(birth_rate(i),1.e-18)
[2513]	146	else
[2730]	147	!!jyg proba_notrig(i)=1.
[3208]	148	birth_rate(i) = 0.
[2513]	149	random_notrig(i)=0.
	150	ale_bl_trig(i)=0.
	151	endif
	152	enddo
	153
	154	ELSE IF (iflag_trig_bl.ge.2) then
	155
	156	do i=1,klon
	157	if ( (Ale_bl(i) .gt. abs(cin(i))+1.e-10) ) then
[2730]	158	proba_notrig(i)=proba_notrig(i)* &
	159	(1.-exp(-s_trig/s2(i)))*(n2(i)dtime/tau_trig(i))
[2513]	160	! print *, 'proba_notrig(i) ',proba_notrig(i)
	161	if (random_notrig(i) .ge. proba_notrig(i)) then
	162	ale_bl_trig(i)=Ale_bl(i)
	163	else
	164	ale_bl_trig(i)=0.
	165	endif
[3208]	166	birth_rate(i) = n2(i)exp(-s_trig/s2(i))/(tau_trig(i)cell_area(i))
	167	!!! birth_rate(i) = max(birth_rate(i),1.e-18)
[2513]	168	else
[2730]	169	!!jyg proba_notrig(i)=1.
[3208]	170	birth_rate(i) = 0.
[2513]	171	random_notrig(i)=0.
	172	ale_bl_trig(i)=0.
	173	endif
	174	enddo
	175
	176	ENDIF
	177
	178	!
	179	IF (prt_level .GE. 10) THEN
	180	print *,'proba_notrig, ale_bl_trig ', &
	181	proba_notrig, ale_bl_trig
	182	ENDIF
	183
	184	endif !(iflag_trig_bl)
	185
	186	!-----------Statistical closure-----------
	187	if (iflag_clos_bl.eq.1) then
	188
	189	do i=1,klon
	190	!CR: alp probabiliste
	191	if (ale_bl_trig(i).gt.0.) then
	192	alp_bl(i)=alp_bl(i)/(1.-min(proba_notrig(i),0.999))
	193	endif
	194	enddo
	195
	196	else if (iflag_clos_bl.eq.2) then
	197
	198	!CR: alp calculee dans thermcell_main
	199	do i=1,klon
	200	alp_bl(i)=alp_bl_stat(i)
	201	enddo
	202
	203	else
	204
	205	alp_bl_stat(:)=0.
	206
	207	endif !(iflag_clos_bl)
	208
[2558]	209	!
	210	!---------------------------------------
	211	ELSEIF (iflag_clos_bl .EQ. 3) THEN ! (iflag_clos_bl .LT. 3)
	212	!---------------------------------------
	213	!
	214	! New code with Effective Lifting Power
	215	! -------------------------------------
	216	IF (first) THEN
	217	cv_feed_area = 1.e10 ! m2
	218	CALL getin_p('cv_feed_area', cv_feed_area)
	219	first=.FALSE.
	220	ENDIF
	221
	222	!-----------Stochastic triggering-----------
	223	if (iflag_trig_bl.ge.1) then
	224	!
	225	IF (prt_level .GE. 10) THEN
	226	print *,'cin, ale_bl_stat, alp_bl_stat ', &
	227	cin, ale_bl_stat, alp_bl_stat
	228	ENDIF
	229
	230	! Use ale_bl_stat (Rochetin's code) or ale_bl (old code) according to
	231	! iflag_trig_bl value.
	232	IF (iflag_trig_bl.eq.1) then ! use ale_bl_stat (Rochetin computation)
	233	do i=1,klon
	234	ale_bl_ref(i)=ale_bl_stat(i)
	235	enddo
	236	ELSE IF (iflag_trig_bl.ge.2) then ! use ale_bl (old computation)
	237	do i=1,klon
	238	ale_bl_ref(i)=Ale_bl(i)
	239	enddo
	240	ENDIF ! (iflag_trig_bl.eq.1)
	241
	242
	243	!----Initializations and random number generation
	244	do i=1,klon
[2730]	245	!!jyg proba_notrig(i)=1.
[2558]	246	random_notrig(i)=1e6ale_bl_stat(i)-int(1e6ale_bl_stat(i))
	247	if ( ale_bl_trig(i) .lt. abs(cin(i))+1.e-10 ) then
	248	tau_trig(i)=tau_trig_shallow
	249	else
	250	tau_trig(i)=tau_trig_deep
	251	endif
	252	enddo
	253	!
	254	IF (prt_level .GE. 10) THEN
	255	print *,'random_notrig, tau_trig ', &
	256	random_notrig, tau_trig
	257	print *,'s_trig,s2,n2 ', &
	258	s_trig,s2,n2
	259	ENDIF
	260
	261	!----alp_bl computation
	262	do i=1,klon
	263	if ( (ale_bl_ref(i) .gt. abs(cin(i))+1.e-10) ) then
	264	birth_number = n2(i)*exp(-s_trig/s2(i))
	265	birth_rate(i) = birth_number/(tau_trig(i)*cell_area(i))
[3208]	266	!!! birth_rate(i) = max(birth_rate(i),1.e-18)
[2730]	267	proba_notrig(i)=proba_notrig(i)exp(-birth_numberdtime/tau_trig(i))
[2558]	268	Alp_bl(i) = Alp_bl(i)* &
	269	umexp(-birth_number*cv_feed_area/cell_area(i))/ &
	270	umexp(-birth_numberdtime/tau_trig(i)) &
	271	tau_trig(i)cv_feed_area/(dtimecell_area(i))
	272	else
[2730]	273	!!jyg proba_notrig(i)=1.
[3208]	274	birth_rate(i)=0.
[2558]	275	random_notrig(i)=0.
	276	alp_bl(i)=0.
	277	endif
	278	enddo
	279
	280	!----ale_bl_trig computation
	281	do i=1,klon
	282	if (random_notrig(i) .ge. proba_notrig(i)) then
	283	ale_bl_trig(i)=ale_bl_ref(i)
	284	else
	285	ale_bl_trig(i)=0.
	286	endif
	287	enddo
	288
	289	!
	290	IF (prt_level .GE. 10) THEN
	291	print *,'proba_notrig, ale_bl_trig ', &
	292	proba_notrig, ale_bl_trig
	293	ENDIF
	294
	295	endif !(iflag_trig_bl .ge. 1)
	296
	297	!---------------------------------------
	298	ENDIF ! (iflag_clos_bl .LT. 3)
	299	!---------------------------------------
	300
[2513]	301	IF (prt_level .GE. 10) THEN
[3209]	302	print *,'alpale_th: ale_bl_trig, alp_bl_stat, birth_rate ', &
	303	ale_bl_trig(1), alp_bl_stat(1), birth_rate(1)
[2513]	304	ENDIF
	305
	306	!cc fin nrlmd le 10/04/2012
[2554]	307	!
[2513]	308	!IM/FH: 2011/02/23
	309	! Couplage Thermiques/Emanuel seulement si T<0
	310	if (iflag_coupl==2) then
	311	IF (prt_level .GE. 10) THEN
	312	print*,'Couplage Thermiques/Emanuel seulement si T<0'
	313	ENDIF
	314	do i=1,klon
	315	if (t_seri(i,lmax_th(i))>273.) then
	316	Ale_bl(i)=0.
	317	endif
	318	enddo
[2554]	319	print *,'In order to run with iflag_coupl=2, you have to comment out the following stop'
	320	STOP
[2513]	321	endif
	322	RETURN
	323	END
	324

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

Download in other formats:

Original Format