Context Navigation

alpale_th.F90 @ 3412

Last change on this file since 3412 was 3209, checked in by jyg, 7 years ago
error in the prvious commit
File size: 11.9 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, birth_rate)
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	REAL, DIMENSION(klon), INTENT(OUT) :: birth_rate
45
46	include "thermcell.h"
47
48	! Local variables
49	!----------------
50	INTEGER :: i
51	LOGICAL, SAVE :: first = .TRUE.
52	REAL, SAVE :: random_notrig_max=1.
53	REAL, SAVE :: cv_feed_area
54	REAL :: birth_number
55	REAL, DIMENSION(klon) :: ale_bl_ref
56	REAL, DIMENSION(klon) :: tau_trig
57	!
58	!$OMP THREADPRIVATE(random_notrig_max)
59	!$OMP THREADPRIVATE(cv_feed_area)
60	!$OMP THREADPRIVATE(first)
61	!
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) + &
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)
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	!
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	!!
89	!---------------------------------------
90	IF (iflag_clos_bl .LT. 3) THEN
91	!---------------------------------------
92	!
93	! Original code (Nicolas Rochetin)
94	! --------------------------------
95
96	IF (first) THEN
97	random_notrig_max=1.
98	CALL getin_p('random_notrig_max',random_notrig_max)
99	first=.FALSE.
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
106	print *,'cin, ale_bl_stat, alp_bl, alp_bl_stat ', &
107	cin, ale_bl_stat, alp_bl, alp_bl_stat
108	ENDIF
109
110
111	!----Initialisations
112	do i=1,klon
113	!!jyg proba_notrig(i)=1.
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
136	proba_notrig(i)=proba_notrig(i)* &
137	(1.-exp(-s_trig/s2(i)))*(n2(i)dtime/tau_trig(i))
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
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)
146	else
147	!!jyg proba_notrig(i)=1.
148	birth_rate(i) = 0.
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
158	proba_notrig(i)=proba_notrig(i)* &
159	(1.-exp(-s_trig/s2(i)))*(n2(i)dtime/tau_trig(i))
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
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)
168	else
169	!!jyg proba_notrig(i)=1.
170	birth_rate(i) = 0.
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
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
245	!!jyg proba_notrig(i)=1.
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))
266	!!! birth_rate(i) = max(birth_rate(i),1.e-18)
267	proba_notrig(i)=proba_notrig(i)exp(-birth_numberdtime/tau_trig(i))
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
273	!!jyg proba_notrig(i)=1.
274	birth_rate(i)=0.
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
301	IF (prt_level .GE. 10) THEN
302	print *,'alpale_th: ale_bl_trig, alp_bl_stat, birth_rate ', &
303	ale_bl_trig(1), alp_bl_stat(1), birth_rate(1)
304	ENDIF
305
306	!cc fin nrlmd le 10/04/2012
307	!
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
319	print *,'In order to run with iflag_coupl=2, you have to comment out the following stop'
320	STOP
321	endif
322	RETURN
323	END
324

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

Download in other formats:

Original Format