Context Navigation

alpale_th.F90 @ 3072

Last change on this file since 3072 was 3072, checked in by jyg, 7 years ago
bug correction in alpale_th for iflag_clos_bl=3
File size: 11.2 KB

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

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

Download in other formats:

Original Format