da_moist_phys_adj.inc
References to this file elsewhere.
1 subroutine da_moist_phys_adj(xb, xa, xp)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Partition of the hydrometeors via the moist explicit scheme.
5 ! A warm rain process is used in this subroutine.
6 ! This is the adjoint code of the scheme.
7 !
8 ! Method: The warm rain process is according to Hsie and Anthes (1984)
9 ! and Dudhia (1989)
10 !
11 ! Assumptions: 1) Model level stored top down.
12 !---------------------------------------------------------------------------
13
14 implicit none
15
16 type (xb_type), intent(in) :: xb ! First guess structure.
17 type (x_type), intent(inout) :: xa ! increment structure.
18 type (xpose_type), intent(inout) :: xp ! Dimensions and xpose buffers.
19
20 real, dimension(ims:ime,jms:jme,kms:kme) :: T_OLD,T_NEW
21 real, dimension(ims:ime,jms:jme,kms:kme) :: Q_OLD,Q_NEW
22 real, dimension(ims:ime,jms:jme,kms:kme) :: QCW_OLD,QCW_NEW
23 real, dimension(ims:ime,jms:jme,kms:kme) :: QRN_OLD,QRN_NEW
24
25 real, dimension(kms:kme) :: EES, QVSS
26 real, dimension(kms:kme) :: EES9, QVSS9
27
28 real, dimension(kms:kme) :: DT
29 real, dimension(kms:kme) :: QVT,QCT,QRT,TTT
30 real, dimension(kms:kme) :: QVT9,QCT9,QRT9,TTT9
31 real, dimension(kms:kme) :: SCR2,SCR3,SCR4,SCR5,SCR6
32 real, dimension(kms:kme) :: DUM31
33 real, dimension(kms:kme) :: PRA,PRC,PRD,PRE
34 real, dimension(kms:kme) :: SCR31,SCR42,SCR71
35 real, dimension(kms:kme) :: DUM112,DUM113,DUM211,DUM411
36 real, dimension(kms:kme) :: PRA2,PRC2
37
38 real, dimension(kms:kme) :: SCR29,SCR39,SCR49,SCR59,SCR69
39 real, dimension(kms:kme) :: DUM319
40 real, dimension(kms:kme) :: PRA9,PRC9,PRD9,PRE9
41 real, dimension(kms:kme) :: SCR319,SCR429,SCR719
42 real, dimension(kms:kme) :: DUM1129,DUM1139,DUM2119,DUM4119
43 real, dimension(kms:kme) :: TMP
44
45 integer :: i, j, k
46
47 do k=kts,kte
48 do j=jts,jte
49 do i=its,ite
50 QRN_NEW(i,j,k) = xa % qrn (i,j,k)
51 QRN_OLD(i,j,k) = xa % qrn (i,j,k)
52 QCW_NEW(i,j,k) = xa % qcw (i,j,k)
53 QCW_OLD(i,j,k) = xa % qcw (i,j,k)
54 Q_NEW(i,j,k) = xa % q (i,j,k)
55 Q_OLD(i,j,k) = xa % q (i,j,k)
56 T_NEW(i,j,k) = xa % t (i,j,k)
57 T_OLD(i,j,k) = xa % t (i,j,k)
58 end do
59 end do
60 end do
61
62 call da_filter_adj(t_new, xp)
63 call da_filter_adj(q_new, xp)
64 call da_filter_adj(qcw_new, xp)
65 call da_filter_adj(qrn_new, xp)
66
67 do k=kts,kte
68 do j=jts,jte
69 do i=its,ite
70 xa % qrn (i,j,k) = QRN_NEW(i,j,k)
71 QRN_OLD(i,j,k) = QRN_OLD(i,j,k) - QRN_NEW(i,j,k)
72 xa % qcw (i,j,k) = QCW_NEW(i,j,k)
73 QCW_OLD(i,j,k) = QCW_OLD(i,j,k) - QCW_NEW(i,j,k)
74 xa % q (i,j,k) = Q_NEW(i,j,k)
75 Q_OLD(i,j,k) = Q_OLD(i,j,k) - Q_NEW(i,j,k)
76 xa % t (i,j,k) = T_NEW(i,j,k)
77 T_OLD(i,j,k) = T_OLD(i,j,k) - T_NEW(i,j,k)
78 end do
79 end do
80 end do
81
82 do j=jts,jte
83 do i=its,ite
84
85 ! Preparation
86
87 do K=kts,kte
88 DT(k) = xb%delt(i,j,k)
89 end do
90
91 do K=kts,kte
92
93 if (DT(k) <= 0.) cycle
94
95 if( xb%t(I,J,K) > TO )then
96 EES(K)=SVP1*EXP(SVP2*(xb%t(I,J,K)-SVPT0)/(xb%t(I,J,K)-SVP3))
97 else
98 EES(K)=.611*EXP(22.514-6.15E3/xb%t(I,J,K))
99 end if
100
101 QVSS(K)=622.*EES(K)/(xb%p(I,J,K)-EES(K))
102
103 SCR4(K)=xb%q(I,J,K)/QVSS(K)
104
105 end do
106
107 do K=kts,kte
108
109 if (DT(k) <= 0.) cycle
110
111 if(xb%qcw(I,J,K) > 0.) then
112 SCR2(K)=xb%qcw(I,J,K)
113 else
114 SCR2(K)=0.
115 end if
116 if(xb%qrn(I,J,K) > 1.E-25) then
117 SCR3(K)=xb%qrn(I,J,K)
118 else
119 SCR3(K)=1.E-25
120 end if
121 SCR5(K)=xb%q(I,J,K)/SCR4(K)
122
123 SCR6(K)=xb%p(I,J,K)/(gas_constant*xb%t(I,J,K))
124
125 DUM31(K)=3.1484E6-XLV1*xb%t(I,J,K)
126
127 end do
128
129 call da_autoc(DT,SCR2,PRC,kts,kte,kms,kme)
130
131 do k=kts,kte
132 PRC2(K)=PRC(K)
133 end do
134
135 call da_accre(DT,SCR2,SCR3,PRA,kts,kte,kms,kme)
136
137 do k=kts,kte
138 PRA2(K)=PRA(K)
139 end do
140
141 call da_evapo(DT,SCR3,SCR5,xb%q(I,J,:),PRE,SCR6,kts,kte,kms,kme)
142
143 do K=kts, kte
144
145 if (DT(k) <= 0.) cycle
146
147 ! Readjust
148
149 DUM112(K)=(PRC(k)+PRA(k))*DT(k)
150 if (DUM112(K) > SCR2(k)) then
151 PRC(k)=SCR2(K)*PRC(K)/DUM112(K)
152 PRA(k)=SCR2(K)*PRA(K)/DUM112(K)
153 end if
154 QVT(K)=-PRE(K)
155 QCT(K)=-PRC(K)-PRA(K)
156 QRT(K)=PRC(K)+PRA(K)+PRE(K)
157 if(xb%t(I,J,K).GT.TO)then
158 DUM411(K)=DUM31(K)
159 else
160 DUM411(K)=XLS
161 end if
162 PRD(K)=cp*(1.+0.887*xb%q(I,J,K))
163 TTT(K)=-DUM411(K)*QVT(K)/PRD(K)
164
165 DUM113(K)=xb%q(I,J,K)+DT(K)*QVT(K)
166 if(DUM113(K) > 1.E-12 ) then
167 SCR42(K)=DUM113(K)
168 else
169 SCR42(K)=1.E-12
170 end if
171 DUM211(K)=xb%qcw(I,J,K)+QCT(K)*DT(k)
172 if(DUM211(K) > 0.) then
173 SCR31(K)=DUM211(K)
174 else
175 SCR31(K)=0.
176 end if
177 SCR71(K)=xb%t(I,J,K)+TTT(K)*DT(k)
178 end do
179
180 call da_condens_adj(DT,SCR31,SCR42,SCR71,DUM31,PRD, &
181 QVT,QCT,QRT,TTT, &
182 xb%p(I,J,:),xb%t(I,J,:),xb%q(I,J,:), &
183 xb%qcw(I,J,:),xb%qrn(I,J,:), &
184 SCR319,SCR429,SCR719,DUM319,PRD9, &
185 QVT9,QCT9,QRT9,TTT9, &
186 xa%p(I,J,:),xa%t(I,J,:),xa%q(I,J,:), &
187 xa%qcw(I,J,:),xa%qrn(I,J,:),kts,kte)
188
189 do K=kts, kte
190 if (DT(k) <= 0.) cycle
191
192 ! Readjust
193
194 xa%t(I,J,K)=xa%t(I,J,K)+SCR719(K)
195 TTT9(K)=TTT9(K)+DT(k)*SCR719(K)
196 DUM2119(K)=0.
197 if(DUM211(K) > 0.) then
198 DUM2119(K)=SCR319(K)
199 end if
200 xa%qcw(I,J,K)=xa%qcw(I,J,K)+DUM2119(K)
201 QCT9(K)=QCT9(K)+DT(K)*DUM2119(K)
202 DUM1139(K)=0.
203 if(DUM113(K) > 1.E-12 ) then
204 DUM1139(K)=SCR429(K)
205 end if
206 xa%q(I,J,K)=xa%q(I,J,K)+DUM1139(K)
207 QVT9(K)=QVT9(K)+DT(K)*DUM1139(K)
208 PRD9(K)=PRD9(K)+DUM411(K)*QVT(K)/(PRD(K)*PRD(K))*TTT9(K)
209 QVT9(K)=QVT9(K)-DUM411(K)/PRD(K)*TTT9(K)
210 DUM4119(K)=-QVT(K)/PRD(K)*TTT9(K)
211 xa%q(I,J,K)=xa%q(I,J,K)+cp*0.887*PRD9(K)
212
213 if(xb%t(I,J,K).GT.TO)then
214 DUM319(K)=DUM319(K)+DUM4119(K)
215 end if
216 PRC9(K)=QRT9(K)
217 PRA9(K)=QRT9(K)
218 PRE9(K)=QRT9(K)
219 PRC9(K)=PRC9(K)-QCT9(K)
220 PRA9(K)=PRA9(K)-QCT9(K)
221 PRE9(K)=PRE9(K)-QVT9(K)
222
223 SCR29(K)=0.
224 DUM1129(K)=0.
225 if (DUM112(K) > SCR2(k)) then
226 DUM1129(K)=-SCR2(K)*PRA2(K)/(DUM112(K)*DUM112(K))*PRA9(K)
227 SCR29(K)=PRA2(K)/DUM112(K)*PRA9(K)
228 PRA9(K)=PRA9(K)*SCR2(K)/DUM112(K)
229 DUM1129(K)=DUM1129(K)-SCR2(K)*PRC2(K)/(DUM112(K)*DUM112(K))*PRC9(K)
230 SCR29(K)=SCR29(K)+PRC2(K)/DUM112(K)*PRC9(K)
231 PRC9(K)=PRC9(K)*SCR2(K)/DUM112(K)
232 end if
233 PRC9(K)=PRC9(K)+DT(K)*DUM1129(K)
234 PRA9(K)=PRA9(K)+DT(K)*DUM1129(K)
235
236 end do
237
238 call da_evapo_adj(DT,SCR3,SCR5,xb%q(I,J,:),PRE,SCR6, &
239 SCR39,SCR59,xa%q(I,J,:),PRE9,SCR69, &
240 kts,kte,kms,kme)
241
242 call da_accre_adj(DT,SCR2,SCR3,PRA,SCR29,SCR39,PRA9,kts,kte,kms,kme)
243
244 call da_autoc_adj(DT,SCR2,PRC,SCR29,PRC9,kts,kte,kms,kme)
245
246 ! Preparation
247
248 do K=kts,kte
249
250 if (DT(k) <= 0.) cycle
251
252 xa%t(I,J,K)=xa%t(I,J,K)-XLV1*DUM319(K)
253
254 xa%p(I,J,K)=xa%p(I,J,K)+SCR69(K)/(gas_constant*xb%t(I,J,K))
255 xa%t(I,J,K)=xa%t(I,J,K)-xb%p(I,J,K)/ &
256 (gas_constant*xb%t(I,J,K)**2)*SCR69(K)
257 xa%q(I,J,K)=xa%q(I,J,K)+SCR59(K)/SCR4(K)
258 SCR49(K)=-xb%q(I,J,K)/SCR4(K)**2*SCR59(K)
259
260 if(xb%qrn(I,J,K) > 1.E-25) then
261 xa%qrn(I,J,K)=xa%qrn(I,J,K)+SCR39(K)
262 end if
263 if(xb%qcw(I,J,K) > 0.) then
264 xa%qcw(I,J,K)=xa%qcw(I,J,K)+SCR29(K)
265 end if
266
267 end do
268
269 do K=kts,kte
270
271 if (DT(k) <= 0.) cycle
272
273 xa%q(I,J,K)=xa%q(I,J,K)+SCR49(K)/QVSS(K)
274 QVSS9(K)=-xb%q(I,J,K)/QVSS(K)**2*SCR49(K)
275 TMP(K)=622./((xb%p(I,J,K)-EES(K))**2)
276 EES9(K)=TMP(K)*xb%p(I,J,K)*QVSS9(K)
277 xa%p(I,J,K)=xa%p(I,J,K)-TMP(K)*EES(K)*QVSS9(K)
278 if( xb%t(I,J,K) > TO )then
279 xa%t(I,J,K)=xa%t(I,J,K)+EES(K)*SVP2*(SVPT0-SVP3)/ &
280 ((xb%t(I,J,K)-SVP3)*(xb%t(I,J,K)-SVP3))*EES9(K)
281 else
282 xa%t(I,J,K)=xa%t(I,J,K)+EES(K)*6.15E3/(xb%t(I,J,K)* &
283 xb%t(I,J,K))*EES9(K)
284 end if
285
286 end do
287
288 do K=kts,kte
289 xa % qt (i,j,k) = xa % qt (i,j,k) + xa % q(i,j,k)
290 xa % qcw(i,j,k) = xa % qcw(i,j,k) - xa % q(i,j,k)
291 xa % qrn(i,j,k) = xa % qrn(i,j,k) - xa % q(i,j,k)
292 end do
293 end do
294 end do
295
296 do k=kts,kte
297 do j=jts,jte
298 do i=its,ite
299 xa % qrn (i,j,k) = xa % qrn (i,j,k) + QRN_OLD(i,j,k)
300 xa % qcw (i,j,k) = xa % qcw (i,j,k) + QCW_OLD(i,j,k)
301 xa % q (i,j,k) = xa % q (i,j,k) + Q_OLD(i,j,k)
302 xa % t (i,j,k) = xa % t (i,j,k) + T_OLD(i,j,k)
303 end do
304 end do
305 end do
306
307 end subroutine da_moist_phys_adj
308
309