Line data Source code
1 : module zm_conv_convtran
2 :
3 : use ccpp_kinds, only: kind_phys
4 :
5 : implicit none
6 :
7 : save
8 : private ! Make default type private to the module
9 : !
10 : ! PUBLIC: interfaces
11 : !
12 : public zm_conv_convtran_run ! convective transport
13 :
14 :
15 : contains
16 :
17 : !===============================================================================
18 : !> \section arg_table_zm_conv_convtran_run Argument Table
19 : !! \htmlinclude zm_conv_convtran_run.html
20 : !!
21 1495368 : subroutine zm_conv_convtran_run(ncol, pver, &
22 1495368 : doconvtran,q ,ncnst ,mu ,md , &
23 1495368 : du ,eu ,ed ,dp ,dsubcld , &
24 2990736 : jt ,mx ,ideep ,il1g ,il2g , &
25 4486104 : nstep ,fracis ,dqdt ,dpdry)
26 : !-----------------------------------------------------------------------
27 : !
28 : ! Purpose:
29 : ! Convective transport of trace species
30 : !
31 : ! Mixing ratios may be with respect to either dry or moist air
32 : !
33 : ! Method:
34 : ! <Describe the algorithm(s) used in the routine.>
35 : ! <Also include any applicable external references.>
36 : !
37 : ! Author: P. Rasch
38 : !
39 : !-----------------------------------------------------------------------
40 : ! CACNOTE - replace with CCPP constituents
41 : use constituents, only: cnst_get_type_byind
42 :
43 : implicit none
44 : !-----------------------------------------------------------------------
45 : !
46 : ! Input arguments
47 : !
48 : integer, intent(in) :: ncol
49 : integer, intent(in) :: pver
50 : integer, intent(in) :: ncnst ! number of tracers to transport
51 : logical, intent(in) :: doconvtran(:) ! flag for doing convective transport (ncnst)
52 : real(kind_phys), intent(in) :: q(:,:,:) ! Tracer array including moisture (ncol,pver,ncnst)
53 : real(kind_phys), intent(in) :: mu(:,:) ! Mass flux up (ncol,pver)
54 : real(kind_phys), intent(in) :: md(:,:) ! Mass flux down (ncol,pver)
55 : real(kind_phys), intent(in) :: du(:,:) ! Mass detraining from updraft (ncol,pver)
56 : real(kind_phys), intent(in) :: eu(:,:) ! Mass entraining from updraft (ncol,pver)
57 : real(kind_phys), intent(in) :: ed(:,:) ! Mass entraining from downdraft (ncol,pver)
58 : real(kind_phys), intent(in) :: dp(:,:) ! Delta pressure between interfaces (ncol,pver)
59 : real(kind_phys), intent(in) :: dsubcld(:) ! Delta pressure from cloud base to sfc (ncol)
60 : real(kind_phys), intent(in) :: fracis(:,:,:) ! fraction of tracer that is insoluble (ncol,pver,ncnst)
61 :
62 : integer, intent(in) :: jt(:) ! Index of cloud top for each column (ncol)
63 : integer, intent(in) :: mx(:) ! Index of cloud top for each column (ncol)
64 : integer, intent(in) :: ideep(:) ! Gathering array (ncol)
65 : integer, intent(in) :: il1g ! Gathered min lon indices over which to operate
66 : integer, intent(in) :: il2g ! Gathered max lon indices over which to operate
67 : integer, intent(in) :: nstep ! Time step index
68 :
69 : real(kind_phys), intent(in) :: dpdry(:,:) ! Delta pressure between interfaces (ncol,pver)
70 :
71 : ! input/output
72 :
73 : real(kind_phys), intent(out) :: dqdt(:,:,:) ! Tracer tendency array (ncol,pver,ncnst)
74 :
75 : !--------------------------Local Variables------------------------------
76 :
77 : integer i ! Work index
78 : integer k ! Work index
79 : integer kbm ! Highest altitude index of cloud base
80 : integer kk ! Work index
81 : integer kkp1 ! Work index
82 : integer km1 ! Work index
83 : integer kp1 ! Work index
84 : integer ktm ! Highest altitude index of cloud top
85 : integer m ! Work index
86 :
87 : real(kind_phys) cabv ! Mix ratio of constituent above
88 : real(kind_phys) cbel ! Mix ratio of constituent below
89 : real(kind_phys) cdifr ! Normalized diff between cabv and cbel
90 2990736 : real(kind_phys) chat(ncol,pver) ! Mix ratio in env at interfaces
91 2990736 : real(kind_phys) cond(ncol,pver) ! Mix ratio in downdraft at interfaces
92 2990736 : real(kind_phys) const(ncol,pver) ! Gathered tracer array
93 2990736 : real(kind_phys) fisg(ncol,pver) ! gathered insoluble fraction of tracer
94 2990736 : real(kind_phys) conu(ncol,pver) ! Mix ratio in updraft at interfaces
95 2990736 : real(kind_phys) dcondt(ncol,pver) ! Gathered tend array
96 : real(kind_phys) small ! A small number
97 : real(kind_phys) mbsth ! Threshold for mass fluxes
98 : real(kind_phys) mupdudp ! A work variable
99 : real(kind_phys) minc ! A work variable
100 : real(kind_phys) maxc ! A work variable
101 : real(kind_phys) fluxin ! A work variable
102 : real(kind_phys) fluxout ! A work variable
103 : real(kind_phys) netflux ! A work variable
104 :
105 2990736 : real(kind_phys) dutmp(ncol,pver) ! Mass detraining from updraft
106 2990736 : real(kind_phys) eutmp(ncol,pver) ! Mass entraining from updraft
107 2990736 : real(kind_phys) edtmp(ncol,pver) ! Mass entraining from downdraft
108 2990736 : real(kind_phys) dptmp(ncol,pver) ! Delta pressure between interfaces
109 : real(kind_phys) total(ncol)
110 : real(kind_phys) negadt,qtmp
111 :
112 : !-----------------------------------------------------------------------
113 : !
114 1495368 : small = 1.e-36_kind_phys
115 : ! mbsth is the threshold below which we treat the mass fluxes as zero (in mb/s)
116 1495368 : mbsth = 1.e-15_kind_phys
117 :
118 : ! Find the highest level top and bottom levels of convection
119 1495368 : ktm = pver
120 1495368 : kbm = pver
121 13077368 : do i = il1g, il2g
122 11582000 : ktm = min(ktm,jt(i))
123 13077368 : kbm = min(kbm,mx(i))
124 : end do
125 :
126 : ! Loop ever each constituent
127 4486104 : do m = 2, ncnst
128 4486104 : if (doconvtran(m)) then
129 :
130 2990736 : if (cnst_get_type_byind(m).eq.'dry') then
131 0 : do k = 1,pver
132 0 : do i =il1g,il2g
133 0 : dptmp(i,k) = dpdry(i,k)
134 0 : dutmp(i,k) = du(i,k)*dp(i,k)/dpdry(i,k)
135 0 : eutmp(i,k) = eu(i,k)*dp(i,k)/dpdry(i,k)
136 0 : edtmp(i,k) = ed(i,k)*dp(i,k)/dpdry(i,k)
137 : end do
138 : end do
139 : else
140 83740608 : do k = 1,pver
141 683013872 : do i =il1g,il2g
142 602264000 : dptmp(i,k) = dp(i,k)
143 602264000 : dutmp(i,k) = du(i,k)
144 602264000 : eutmp(i,k) = eu(i,k)
145 680023136 : edtmp(i,k) = ed(i,k)
146 : end do
147 : end do
148 : endif
149 : ! dptmp = dp
150 :
151 : ! Gather up the constituent and set tend to zero
152 80749872 : do k = 1,pver
153 683013872 : do i =il1g,il2g
154 602264000 : const(i,k) = q(ideep(i),k,m)
155 680023136 : fisg(i,k) = fracis(ideep(i),k,m)
156 : end do
157 : end do
158 :
159 : ! From now on work only with gathered data
160 :
161 : ! Interpolate environment tracer values to interfaces
162 80749872 : do k = 1,pver
163 77759136 : km1 = max(1,k-1)
164 683013872 : do i = il1g, il2g
165 602264000 : minc = min(const(i,km1),const(i,k))
166 602264000 : maxc = max(const(i,km1),const(i,k))
167 602264000 : if (minc < 0) then
168 : cdifr = 0._kind_phys
169 : else
170 602264000 : cdifr = abs(const(i,k)-const(i,km1))/max(maxc,small)
171 : endif
172 :
173 : ! If the two layers differ significantly use a geometric averaging
174 : ! procedure
175 602264000 : if (cdifr > 1.E-6_kind_phys) then
176 414024028 : cabv = max(const(i,km1),maxc*1.e-12_kind_phys)
177 414024028 : cbel = max(const(i,k),maxc*1.e-12_kind_phys)
178 414024028 : chat(i,k) = log(cabv/cbel)/(cabv-cbel)*cabv*cbel
179 :
180 : else ! Small diff, so just arithmetic mean
181 188239972 : chat(i,k) = 0.5_kind_phys* (const(i,k)+const(i,km1))
182 : end if
183 :
184 : ! Provisional up and down draft values
185 602264000 : conu(i,k) = chat(i,k)
186 602264000 : cond(i,k) = chat(i,k)
187 :
188 : ! provisional tends
189 680023136 : dcondt(i,k) = 0._kind_phys
190 :
191 : end do
192 : end do
193 :
194 : ! Do levels adjacent to top and bottom
195 26154736 : k = 2
196 26154736 : km1 = 1
197 26154736 : kk = pver
198 26154736 : do i = il1g,il2g
199 23164000 : mupdudp = mu(i,kk) + dutmp(i,kk)*dptmp(i,kk)
200 23164000 : if (mupdudp > mbsth) then
201 22223084 : conu(i,kk) = (+eutmp(i,kk)*fisg(i,kk)*const(i,kk)*dptmp(i,kk))/mupdudp
202 : endif
203 26154736 : if (md(i,k) < -mbsth) then
204 0 : cond(i,k) = (-edtmp(i,km1)*fisg(i,km1)*const(i,km1)*dptmp(i,km1))/md(i,k)
205 : endif
206 : end do
207 :
208 : ! Updraft from bottom to top
209 77759136 : do kk = pver-1,1,-1
210 74768400 : kkp1 = min(pver,kk+1)
211 656859136 : do i = il1g,il2g
212 579100000 : mupdudp = mu(i,kk) + dutmp(i,kk)*dptmp(i,kk)
213 653868400 : if (mupdudp > mbsth) then
214 118024902 : conu(i,kk) = ( mu(i,kkp1)*conu(i,kkp1)+eutmp(i,kk)*fisg(i,kk)* &
215 236049804 : const(i,kk)*dptmp(i,kk) )/mupdudp
216 : endif
217 : end do
218 : end do
219 :
220 : ! Downdraft from top to bottom
221 74768400 : do k = 3,pver
222 71777664 : km1 = max(1,k-1)
223 630704400 : do i = il1g,il2g
224 627713664 : if (md(i,k) < -mbsth) then
225 196768404 : cond(i,k) = ( md(i,km1)*cond(i,km1)-edtmp(i,km1)*fisg(i,km1)*const(i,km1) &
226 196768404 : *dptmp(i,km1) )/md(i,k)
227 : endif
228 : end do
229 : end do
230 :
231 :
232 18397672 : do k = ktm,pver
233 15406936 : km1 = max(1,k-1)
234 15406936 : kp1 = min(pver,k+1)
235 186633220 : do i = il1g,il2g
236 :
237 : ! limit fluxes outside convection to mass in appropriate layer
238 : ! these limiters are probably only safe for positive definite quantitities
239 : ! it assumes that mu and md already satify a courant number limit of 1
240 672942192 : fluxin = mu(i,kp1)*conu(i,kp1)+ mu(i,k)*min(chat(i,k),const(i,km1)) &
241 841177740 : -(md(i,k) *cond(i,k) + md(i,kp1)*min(chat(i,kp1),const(i,kp1)))
242 : fluxout = mu(i,k)*conu(i,k) + mu(i,kp1)*min(chat(i,kp1),const(i,k)) &
243 168235548 : -(md(i,kp1)*cond(i,kp1) + md(i,k)*min(chat(i,k),const(i,k)))
244 :
245 168235548 : netflux = fluxin - fluxout
246 168235548 : if (abs(netflux) < max(fluxin,fluxout)*1.e-12_kind_phys) then
247 12612158 : netflux = 0._kind_phys
248 : endif
249 183642484 : dcondt(i,k) = netflux/dptmp(i,k)
250 : end do
251 : end do
252 : ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
253 : !
254 6246608 : do k = kbm,pver
255 3255872 : km1 = max(1,k-1)
256 31958224 : do i = il1g,il2g
257 28967488 : if (k == mx(i)) then
258 :
259 23164000 : fluxin = mu(i,k)*min(chat(i,k),const(i,km1)) - md(i,k)*cond(i,k)
260 23164000 : fluxout = mu(i,k)*conu(i,k) - md(i,k)*min(chat(i,k),const(i,k))
261 :
262 23164000 : netflux = fluxin - fluxout
263 23164000 : if (abs(netflux) < max(fluxin,fluxout)*1.e-12_kind_phys) then
264 0 : netflux = 0._kind_phys
265 : endif
266 23164000 : dcondt(i,k) = netflux/dptmp(i,k)
267 2547616 : else if (k > mx(i)) then
268 758318 : dcondt(i,k) = 0._kind_phys
269 : end if
270 : end do
271 : end do
272 :
273 : ! Initialize to zero everywhere, then scatter tendency back to full array
274 1301387472 : dqdt(:,:,m) = 0._kind_phys
275 80749872 : do k = 1,pver
276 77759136 : kp1 = min(pver,k+1)
277 683013872 : do i = il1g,il2g
278 680023136 : dqdt(ideep(i),k,m) = dcondt(i,k)
279 : end do
280 : end do
281 :
282 : end if ! for doconvtran
283 :
284 : end do
285 :
286 1495368 : return
287 : end subroutine zm_conv_convtran_run
288 :
289 :
290 : end module zm_conv_convtran
|
