CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
progsigma_calc.f90
1
6 module progsigma
7
8 implicit none
9
10 public progsigma_calc
11
12 contains
13
19 subroutine progsigma_calc (im,km,flag_init,flag_restart,flag_shallow,&
20 flag_mid,del,tmf,qmicro,dbyo1,zdqca,omega_u,zeta,hvap, &
21 delt,qadv,kbcon1,ktcon,cnvflg,betascu,betamcu,betadcu, &
22 sigmind,sigminm,sigmins,sigmain,sigmaout,sigmab)
23!
24!
25 use machine, only : kind_phys
26 use funcphys, only : fpvs
27
28 implicit none
29
30! intent in
31 integer, intent(in) :: im,km,kbcon1(im),ktcon(im)
32 real(kind=kind_phys), intent(in) :: hvap,delt,betascu,betamcu,betadcu, &
33 sigmind,sigminm,sigmins
34 real(kind=kind_phys), intent(in) :: qadv(im,km),del(im,km), &
35 qmicro(im,km),tmf(im,km),dbyo1(im,km),zdqca(im,km), &
36 omega_u(im,km),zeta(im,km)
37 logical, intent(in) :: flag_init,flag_restart,cnvflg(im),flag_shallow,flag_mid
38 real(kind=kind_phys), intent(in) :: sigmain(im,km)
39
40! intent out
41 real(kind=kind_phys), intent(out) :: sigmaout(im,km)
42 real(kind=kind_phys), intent(out) :: sigmab(im)
43
44
45! Local variables
46 integer :: i,k,km1
47 real(kind=kind_phys) :: terma(im),termb(im),termc(im),termd(im)
48 real(kind=kind_phys) :: mcons(im),fdqa(im),form(im,km), &
49 dp(im,km),inbu(im,km)
50
51
52 real(kind=kind_phys) :: gcvalmx,epsilon,zz,cvg,mcon,buy2, &
53 fdqb,dtdyn,dxlim,rmulacvg,tem, &
54 den,dp1,invdelt
55
56 !Parameters
57 gcvalmx = 0.1
58 rmulacvg=10.
59 epsilon=1.e-11
60 km1=km-1
61 invdelt = 1./delt
62
63 !Initialization 2D
64 do k = 1,km
65 do i = 1,im
66 inbu(i,k)=0.
67 form(i,k)=0.
68 dp(i,k)=0.
69 enddo
70 enddo
71
72 !Initialization 1D
73 do i=1,im
74 sigmab(i)=0.
75 terma(i)=0.
76 termb(i)=0.
77 termc(i)=0.
78 termd(i)=0.
79 fdqa(i)=0.
80 mcons(i)=0.
81 enddo
82
83 do k = 2,km1
84 do i = 1,im
85 if(cnvflg(i))then
86 dp(i,k) = 1000. * del(i,k)
87 endif
88 enddo
89 enddo
90
91 !Initial computations, place maximum sigmain in sigmab
92 do i=1,im
93 if(cnvflg(i))then
94 do k=2,km
95 if(sigmain(i,k)>sigmab(i))then
96 sigmab(i)=sigmain(i,k)
97 endif
98 enddo
99 endif
100 enddo
101
102 do i=1,im
103 if(cnvflg(i))then
104 if(sigmab(i) < 1.e-5)then !after advection
105 sigmab(i)=0.
106 endif
107 endif
108 enddo
109
110 !compute termD "The vertical integral of the latent heat convergence is limited to the
111 !buoyant layers with positive moisture convergence (accumulated from the surface).
112 !Lowest level:
113 do i = 1,im
114 dp1 = 1000. * del(i,1)
115 mcons(i)=(hvap*(qadv(i,1)+tmf(i,1)+qmicro(i,1))*dp1)
116 enddo
117 !Levels above:
118 do k = 2,km1
119 do i = 1,im
120 if(cnvflg(i))then
121 mcon = (hvap*(qadv(i,k)+tmf(i,k)+qmicro(i,k))*dp(i,k))
122 buy2 = termd(i)+mcon+mcons(i)
123! Do the integral over buoyant layers with positive mcon acc from surface
124 if(dbyo1(i,k)>0 .and. buy2 > 0.)then
125 inbu(i,k)=1.
126 endif
127 inbu(i,k-1)=max(inbu(i,k-1),inbu(i,k))
128 termd(i) = termd(i) + inbu(i,k-1)*mcons(i)
129 mcons(i)=mcon
130 endif
131 enddo
132 enddo
133
134 !termA
135 do k = 2,km1
136 do i = 1,im
137 if(cnvflg(i))then
138 tem=sigmab(i)*zeta(i,k)*inbu(i,k)*dbyo1(i,k)*dp(i,k)
139 terma(i)=terma(i)+tem
140 endif
141 enddo
142 enddo
143
144 !termB
145 do k = 2,km1
146 do i = 1,im
147 if(cnvflg(i))then
148 tem=zeta(i,k)*dbyo1(i,k)*inbu(i,k)*dp(i,k)
149 termb(i)=termb(i)+tem
150 endif
151 enddo
152 enddo
153
154 !termC
155 do k = 2,km1
156 do i = 1,im
157 if(cnvflg(i))then
158 form(i,k)=-1.0*inbu(i,k)*(omega_u(i,k)*delt)
159 fdqb=0.5*((form(i,k)*zdqca(i,k)))
160 termc(i)=termc(i)+inbu(i,k)* &
161 (fdqb+fdqa(i))*hvap*zeta(i,k)
162 fdqa(i)=fdqb
163 endif
164 enddo
165 enddo
166
167 !sigmab
168 if(flag_init .and. .not. flag_restart)then
169 do i = 1,im
170 if(cnvflg(i))then
171 sigmab(i)=0.03
172 endif
173 enddo
174 else
175 do i = 1,im
176 if(cnvflg(i))then
177 den=min(termc(i)+termb(i),1.e8)
178 cvg=termd(i)*delt
179 zz=max(0.0,sign(1.0,terma(i))) &
180 *max(0.0,sign(1.0,termb(i))) &
181 *max(0.0,sign(1.0,termc(i)-epsilon))
182 cvg=max(0.0,cvg)
183 sigmab(i)=(zz*(terma(i)+cvg))/(den+(1.0-zz))
184 if(sigmab(i)>0.)then
185 sigmab(i)=min(sigmab(i),0.95)
186 sigmab(i)=max(sigmab(i),0.01)
187 endif
188 endif!cnvflg
189 enddo
190 endif
191
192 do k=1,km
193 do i=1,im
194 if(cnvflg(i))then
195 sigmaout(i,k)=sigmab(i)
196 endif
197 enddo
198 enddo
199
200 !Reduce area fraction before coupling back to mass-flux computation.
201 if(flag_shallow)then
202 do i= 1, im
203 if(cnvflg(i)) then
204 sigmab(i)=sigmab(i)/betascu
205 sigmab(i)=max(sigmins,sigmab(i))
206 endif
207 enddo
208 elseif(flag_mid)then
209 do i= 1, im
210 if(cnvflg(i)) then
211 sigmab(i)=sigmab(i)/betamcu
212 sigmab(i)=max(sigminm,sigmab(i))
213 endif
214 enddo
215 else
216 do i= 1, im
217 if(cnvflg(i)) then
218 sigmab(i)=sigmab(i)/betadcu
219 sigmab(i)=max(sigmind,sigmab(i))
220 endif
221 enddo
222 endif
223 do i= 1, im
224 sigmab(i) = min(0.95,sigmab(i))
225 enddo
226
227 end subroutine progsigma_calc
228
229end module progsigma