GCC Code Coverage Report


Directory: ./
File: phys/thermcell_alim.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 18 50 36.0%
Branches: 19 60 31.7%

Line Branch Exec Source
1 !
2 ! $Id: thermcell_plume.F90 2311 2015-06-25 07:45:24Z emillour $
3 !
4 480 SUBROUTINE thermcell_alim(flag,ngrid,klev,ztv,d_temp,zlev,alim_star,lalim)
5 IMPLICIT NONE
6
7 !--------------------------------------------------------------------------
8 ! FH : 2015/11/06
9 ! thermcell_alim: calcule la distribution verticale de l'alimentation
10 ! laterale a la base des panaches thermiques
11 !--------------------------------------------------------------------------
12
13 !
14 ! $Header$
15 !
16 ! ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
17 ! veillez � n'utiliser que des ! pour les commentaires
18 ! et � bien positionner les & des lignes de continuation
19 ! (les placer en colonne 6 et en colonne 73)
20 !
21 !
22 ! A1.0 Fundamental constants
23 REAL RPI,RCLUM,RHPLA,RKBOL,RNAVO
24 ! A1.1 Astronomical constants
25 REAL RDAY,REA,REPSM,RSIYEA,RSIDAY,ROMEGA
26 ! A1.1.bis Constantes concernant l'orbite de la Terre:
27 REAL R_ecc, R_peri, R_incl
28 ! A1.2 Geoide
29 REAL RA,RG,R1SA
30 ! A1.3 Radiation
31 ! REAL RSIGMA,RI0
32 REAL RSIGMA
33 ! A1.4 Thermodynamic gas phase
34 REAL RMO3,RMCO2,RMC,RMCH4,RMN2O,RMCFC11,RMCFC12
35 REAL R,RMD,RMV,RD,RV,RCPD,RCPV,RCVD,RCVV
36 REAL RKAPPA,RETV, eps_w
37 ! A1.5,6 Thermodynamic liquid,solid phases
38 REAL RCW,RCS
39 ! A1.7 Thermodynamic transition of phase
40 REAL RLVTT,RLSTT,RLMLT,RTT,RATM
41 ! A1.8 Curve of saturation
42 REAL RESTT,RALPW,RBETW,RGAMW,RALPS,RBETS,RGAMS
43 REAL RALPD,RBETD,RGAMD
44 !
45 COMMON/YOMCST/RPI ,RCLUM ,RHPLA ,RKBOL ,RNAVO &
46 & ,RDAY ,REA ,REPSM ,RSIYEA,RSIDAY,ROMEGA &
47 & ,R_ecc, R_peri, R_incl &
48 & ,RA ,RG ,R1SA &
49 & ,RSIGMA &
50 & ,R ,RMD ,RMV ,RD ,RV ,RCPD &
51 & ,RMO3 ,RMCO2 ,RMC ,RMCH4 ,RMN2O ,RMCFC11 ,RMCFC12 &
52 & ,RCPV ,RCVD ,RCVV ,RKAPPA,RETV, eps_w &
53 & ,RCW ,RCS &
54 & ,RLVTT ,RLSTT ,RLMLT ,RTT ,RATM &
55 & ,RESTT ,RALPW ,RBETW ,RGAMW ,RALPS ,RBETS ,RGAMS &
56 & ,RALPD ,RBETD ,RGAMD
57 ! ------------------------------------------------------------------
58 !$OMP THREADPRIVATE(/YOMCST/)
59 !
60 ! $Id: YOETHF.h 2799 2017-02-24 18:50:33Z jyg $
61 !
62 ! ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
63 ! veillez n'utiliser que des ! pour les commentaires
64 ! et bien positionner les & des lignes de continuation
65 ! (les placer en colonne 6 et en colonne 73)
66 !
67 !* COMMON *YOETHF* DERIVED CONSTANTS SPECIFIC TO ECMWF THERMODYNAMICS
68 !
69 ! *R__ES* *CONSTANTS USED FOR COMPUTATION OF SATURATION
70 ! MIXING RATIO OVER LIQUID WATER(*R_LES*) OR
71 ! ICE(*R_IES*).
72 ! *RVTMP2* *RVTMP2=RCPV/RCPD-1.
73 ! *RHOH2O* *DENSITY OF LIQUID WATER. (RATM/100.)
74 !
75 REAL R2ES, R3LES, R3IES, R4LES, R4IES, R5LES, R5IES
76 REAL RVTMP2, RHOH2O
77 REAL R5ALVCP,R5ALSCP,RALVDCP,RALSDCP,RALFDCP,RTWAT,RTBER,RTBERCU
78 REAL RTICE,RTICECU,RTWAT_RTICE_R,RTWAT_RTICECU_R,RKOOP1,RKOOP2
79 LOGICAL OK_BAD_ECMWF_THERMO ! If TRUE, then variables set by rrtm/suphec.F90
80 ! If FALSE, then variables set by suphel.F90
81 COMMON /YOETHF/R2ES, R3LES, R3IES, R4LES, R4IES, R5LES, R5IES, &
82 & RVTMP2, RHOH2O, &
83 & R5ALVCP,R5ALSCP,RALVDCP,RALSDCP, &
84 & RALFDCP,RTWAT,RTBER,RTBERCU, &
85 & RTICE,RTICECU,RTWAT_RTICE_R,RTWAT_RTICECU_R,RKOOP1,&
86 & RKOOP2, &
87 & OK_BAD_ECMWF_THERMO
88
89 !$OMP THREADPRIVATE(/YOETHF/)
90 !
91 ! $Header$
92 !
93 !
94 ! ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
95 ! veillez n'utiliser que des ! pour les commentaires
96 ! et bien positionner les & des lignes de continuation
97 ! (les placer en colonne 6 et en colonne 73)
98 !
99 ! ------------------------------------------------------------------
100 ! This COMDECK includes the Thermodynamical functions for the cy39
101 ! ECMWF Physics package.
102 ! Consistent with YOMCST Basic physics constants, assuming the
103 ! partial pressure of water vapour is given by a first order
104 ! Taylor expansion of Qs(T) w.r.t. to Temperature, using constants
105 ! in YOETHF
106 ! ------------------------------------------------------------------
107 REAL PTARG, PDELARG, P5ARG, PQSARG, PCOARG
108 REAL FOEEW, FOEDE, qsats, qsatl, dqsats, dqsatl
109 LOGICAL thermcep
110 PARAMETER (thermcep=.TRUE.)
111 !
112 FOEEW ( PTARG,PDELARG ) = EXP ( &
113 & (R3LES*(1.-PDELARG)+R3IES*PDELARG) * (PTARG-RTT) &
114 & / (PTARG-(R4LES*(1.-PDELARG)+R4IES*PDELARG)) )
115 !
116 FOEDE ( PTARG,PDELARG,P5ARG,PQSARG,PCOARG ) = PQSARG*PCOARG*P5ARG &
117 & / (PTARG-(R4LES*(1.-PDELARG)+R4IES*PDELARG))**2
118 !
119 qsats(ptarg) = 100.0 * 0.622 * 10.0 &
120 & ** (2.07023 - 0.00320991 * ptarg &
121 & - 2484.896 / ptarg + 3.56654 * LOG10(ptarg))
122 qsatl(ptarg) = 100.0 * 0.622 * 10.0 &
123 & ** (23.8319 - 2948.964 / ptarg &
124 & - 5.028 * LOG10(ptarg) &
125 & - 29810.16 * EXP( - 0.0699382 * ptarg) &
126 & + 25.21935 * EXP( - 2999.924 / ptarg))
127 !
128 dqsats(ptarg,pqsarg) = RLVTT/RCPD*pqsarg * (3.56654/ptarg &
129 & +2484.896*LOG(10.)/ptarg**2 &
130 & -0.00320991*LOG(10.))
131 dqsatl(ptarg,pqsarg) = RLVTT/RCPD*pqsarg*LOG(10.)* &
132 & (2948.964/ptarg**2-5.028/LOG(10.)/ptarg &
133 & +25.21935*2999.924/ptarg**2*EXP(-2999.924/ptarg) &
134 & +29810.16*0.0699382*EXP(-0.0699382*ptarg))
135 integer :: iflag_thermals,nsplit_thermals
136
137 !!! nrlmd le 10/04/2012
138 integer :: iflag_trig_bl,iflag_clos_bl
139 integer :: tau_trig_shallow,tau_trig_deep
140 real :: s_trig
141 !!! fin nrlmd le 10/04/2012
142
143 real,parameter :: r_aspect_thermals=2.,l_mix_thermals=30.
144 real :: alp_bl_k
145 real :: tau_thermals,fact_thermals_ed_dz
146 integer,parameter :: w2di_thermals=0
147 integer :: isplit
148
149 integer :: iflag_coupl,iflag_clos,iflag_wake
150 integer :: iflag_thermals_ed,iflag_thermals_optflux,iflag_thermals_closure
151
152 common/ctherm1/iflag_thermals,nsplit_thermals,iflag_thermals_closure
153 common/ctherm2/tau_thermals,alp_bl_k,fact_thermals_ed_dz
154 common/ctherm4/iflag_coupl,iflag_clos,iflag_wake
155 common/ctherm5/iflag_thermals_ed,iflag_thermals_optflux
156
157 !!! nrlmd le 10/04/2012
158 common/ctherm6/iflag_trig_bl,iflag_clos_bl
159 common/ctherm7/tau_trig_shallow,tau_trig_deep
160 common/ctherm8/s_trig
161 !!! fin nrlmd le 10/04/2012
162
163 !$OMP THREADPRIVATE(/ctherm1/,/ctherm2/,/ctherm4/,/ctherm5/)
164 !$OMP THREADPRIVATE(/ctherm6/,/ctherm7/,/ctherm8/)
165
166 ! fort(10) ptimestep,ztv,zthl,po,zl,rhobarz,zlev,pplev,pphi,zpspsk,f0
167 INTEGER, INTENT(IN) :: ngrid,klev
168 REAL, INTENT(IN) :: ztv(ngrid,klev)
169 REAL, INTENT(IN) :: d_temp(ngrid)
170 REAL, INTENT(IN) :: zlev(ngrid,klev+1)
171 REAL, INTENT(OUT) :: alim_star(ngrid,klev)
172 INTEGER, INTENT(OUT) :: lalim(ngrid)
173 INTEGER, INTENT(IN) :: flag
174
175 960 REAL :: alim_star_tot(ngrid),zi(ngrid),zh(ngrid)
176 480 REAL :: zlay(ngrid,klev)
177 REAL ztv_parcel
178
179 INTEGER ig,l
180
181 REAL h,z,falim
182 falim(h,z)=0.2*((z-h)**5+h**5)
183
184
185 !===================================================================
186
187
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 lalim(:)=1
188
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 alim_star_tot(:)=0.
189
190 !-------------------------------------------------------------------------
191 ! Definition de l'alimentation a l'origine dans thermcell_init
192 !-------------------------------------------------------------------------
193
1/2
✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
480 IF (flag==0) THEN ! CMIP5 version
194
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 do l=1,klev-1
195
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 do ig=1,ngrid
196
4/4
✓ Branch 0 taken 645327 times.
✓ Branch 1 taken 17485233 times.
✓ Branch 2 taken 606965 times.
✓ Branch 3 taken 38362 times.
18148800 if (ztv(ig,l)> ztv(ig,l+1) .and. ztv(ig,1)>=ztv(ig,l) ) then
197 alim_star(ig,l)=MAX((ztv(ig,l)-ztv(ig,l+1)),0.) &
198 606965 & *sqrt(zlev(ig,l+1))
199 606965 lalim(ig)=l+1
200 606965 alim_star_tot(ig)=alim_star_tot(ig)+alim_star(ig,l)
201 endif
202 enddo
203 enddo
204
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 do l=1,klev
205
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 do ig=1,ngrid
206
2/2
✓ Branch 0 taken 9688770 times.
✓ Branch 1 taken 8918910 times.
18626400 if (alim_star_tot(ig) > 1.e-10 ) then
207 9688770 alim_star(ig,l)=alim_star(ig,l)/alim_star_tot(ig)
208 endif
209 enddo
210 enddo
211 480 alim_star_tot(:)=1.
212
213 !-------------------------------------------------------------------------
214 ! Nouvelle definition avec possibilite d'introduire un DT en surface
215 ! On suppose que la forme du profile d'alimentation scale avec la hauteur
216 ! d'inversion calculée avec une particule partant de la premieere couche
217
218 ! Fonction f(z) = z ( h - z ) , avec h = zi/3
219 ! On utilise l'integralle
220 ! Int_0^z f(z') dz' = z^2 ( h/2 - z/3 ) = falim(h,z)
221 ! Pour calculer l'alimentation des couches
222 !-------------------------------------------------------------------------
223 ELSE
224 ! Computing inversion height zi and zh=zi/3.
225 zi(:)=0.
226 ! Il faut recalculer zlay qui n'est pas dispo dans thermcell_plume
227 ! A changer eventuellement.
228 do l=1,klev
229 zlay(:,l)=0.5*(zlev(:,l)+zlev(:,l+1))
230 enddo
231
232 do l=klev-1,1,-1
233 do ig=1,ngrid
234 ztv_parcel=ztv(ig,1)+d_temp(ig)
235 if (ztv_parcel<ztv(ig,l+1)) lalim(ig)=l
236 enddo
237 enddo
238
239 do ig=1,ngrid
240 l=lalim(ig)
241 IF (l==1) THEN
242 zi(ig)=0.
243 ELSE
244 ztv_parcel=ztv(ig,1)+d_temp(ig)
245 zi(ig)=zlay(ig,l)+(zlay(ig,l+1)-zlay(ig,l))/(ztv(ig,l+1)-ztv(ig,l))*(ztv_parcel-ztv(ig,l))
246 ENDIF
247 enddo
248
249 zh(:)=zi(:)/2.
250 alim_star_tot(:)=0.
251 alim_star(:,:)=0.
252 lalim(:)=0
253 do l=1,klev-1
254 do ig=1,ngrid
255 IF (zh(ig)==0.) THEN
256 alim_star(ig,l)=0.
257 lalim(ig)=1
258 ELSE IF (zlev(ig,l+1)<=zh(ig)) THEN
259 alim_star(ig,l)=(falim(zh(ig),zlev(ig,l+1))-falim(zh(ig),zlev(ig,l)))/falim(zh(ig),zh(ig))
260 lalim(ig)=l
261 ELSE IF (zlev(ig,l)<=zh(ig)) THEN
262 alim_star(ig,l)=(falim(zh(ig),zh(ig))-falim(zh(ig),zlev(ig,l)))/falim(zh(ig),zh(ig))
263 lalim(ig)=l
264 ELSE
265 alim_star(ig,l)=0.
266 ENDIF
267 ENDDO
268 alim_star_tot(:)=alim_star_tot(:)+alim_star(:,l)
269 ENDDO
270 IF (ngrid==1) print*,'NEW ALIM CALCUL DE ZI ',alim_star_tot,lalim,zi,zh
271 alim_star_tot(:)=1.
272
273 ENDIF
274
275
276 480 RETURN
277 END
278