| blob | 62e1623de0625c73b97a8b5c378d00b298d00427 |
1 ! Generated by TAPENADE (INRIA, Tropics team)
2 ! Tapenade 3.7 (r4786) - 21 Feb 2013 15:53
3 !
4 ! Differentiation of ducu in reverse (adjoint) mode (with options r8):
5 ! gradient of useful results: th raincv t rthcuten pcps qv
6 ! z pratec rqvcuten rho dz8w
7 ! with respect to varying inputs: th raincv t rthcuten pcps qv
8 ! z pratec rqvcuten rho dz8w
9 MODULE a_module_cu_du
10 USE module_wrf_error
11 REAL , PARAMETER :: cincap = -10.
12 REAL , PARAMETER :: capemin = 10.
13 REAL , PARAMETER :: dpthmin = 1000.
14 REAL , PARAMETER :: alpha = 0.00002
15 REAL , PARAMETER :: eps = 0.5
16 REAL , PARAMETER :: Vfall = 5.
17 !--------------------------------------------------------------------
18 CONTAINS
19 SUBROUTINE DUCU_B(ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms&
20 & , kme, its, ite, jts, jte, kts, kte, dt, ktau, dx, rho, rhob, raincv, &
21 & raincvb, nca, pratec, pratecb, u, v, th, thb, t, tb, w, dz8w, dz8wb, z&
22 & , zb, pcps, pcpsb, pi, w0avg, cp, rd, rv, g, xlv, ep2, svp1, svp2, &
23 & svp3, svpt0, stepcu, cu_act_flag, warm_rain, cutop, cubot, qv, qvb, &
24 & rthcuten, rthcutenb, rqvcuten, rqvcutenb)
25 IMPLICIT NONE
26 !
27 !-------------------------------------------------------------
28 INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
29 & jme, kms, kme, its, ite, jts, jte, kts, kte
30 INTEGER, INTENT(IN) :: stepcu
31 LOGICAL, INTENT(IN) :: warm_rain
32 REAL, INTENT(IN) :: xlv
33 REAL, INTENT(IN) :: cp, rd, rv, g, ep2
34 REAL, INTENT(IN) :: svp1, svp2, svp3, svpt0
35 INTEGER, INTENT(IN) :: ktau
36 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: u, v, w, th&
37 & , t, qv, dz8w, z, pcps, rho, pi
38 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: thb, tb, qvb, dz8wb, zb&
39 & , pcpsb, rhob
40 !
41 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: w0avg
42 REAL, INTENT(IN) :: dt, dx
43 !
44 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: raincv, pratec
45 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: raincvb
46 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
47 REAL, DIMENSION(ims:ime, jms:jme), INTENT(OUT) :: cubot, cutop
48 LOGICAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: cu_act_flag
49 !
50 ! Optional arguments
51 !
52 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL, INTENT(INOUT) ::&
53 & rthcuten, rqvcuten
54 REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL :: rthcutenb, &
55 & rqvcutenb
56 !
57 ! LOCAL VARS
58 LOGICAL :: flag_qr, flag_qi, flag_qs
59 REAL, DIMENSION(kts:kte) :: u1d, v1d, t1d, th1d, dz1d, z1d, qv1d, p1d&
60 & , rho1d, w0avg1d
61 REAL, DIMENSION(kts:kte) :: t1db, th1db, dz1db, z1db, qv1db, p1db, &
62 & rho1db
63 REAL, DIMENSION(kts:kte) :: dqvdt, dthdt
64 REAL, DIMENSION(kts:kte) :: dqvdtb, dthdtb
65 REAL :: pprate, tst, tv, prs, rhoe, w0, scr1, dxsq, rthcumax
66 REAL :: pprateb
67 INTEGER :: i, j, k, i_start, i_end, j_start, j_end, sz, ntst, icldck
68 INTEGER :: branch
69 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: pratecb
70 !
71 ntst = stepcu
72 icldck = MOD(ktau, ntst)
73 IF (icldck .EQ. 0 .OR. ktau .EQ. 1) THEN
74 !
75 ! Keep away from specified and relaxation zone (should be for just specified and nested bc)
76 sz = 1
77 IF (ids + sz .LT. its) THEN
78 i_start = its
79 ELSE
80 i_start = ids + sz
81 END IF
82 IF (ide - 1 - sz .GT. ite) THEN
83 i_end = ite
84 ELSE
85 i_end = ide - 1 - sz
86 END IF
87 IF (jds + sz .LT. jts) THEN
88 j_start = jts
89 ELSE
90 j_start = jds + sz
91 END IF
92 IF (jde - 1 - sz .GT. jte) THEN
93 j_end = jte
94 ELSE
95 j_end = jde - 1 - sz
96 END IF
97 !
98 DO j=j_start,j_end
99 DO i=i_start,i_end
100 !
101 ! assign vars from 3D to 1D
102 DO k=kts,kte
103 CALL PUSHREAL8(t1d(k))
104 t1d(k) = t(i, k, j)
105 CALL PUSHREAL8(th1d(k))
106 th1d(k) = th(i, k, j)
107 CALL PUSHREAL8(rho1d(k))
108 rho1d(k) = rho(i, k, j)
109 CALL PUSHREAL8(qv1d(k))
110 qv1d(k) = qv(i, k, j)
111 IF (qv1d(k) .LT. 1.e-08) THEN
112 qv1d(k) = 1.e-08
113 CALL PUSHCONTROL1B(0)
114 ELSE
115 CALL PUSHCONTROL1B(1)
116 END IF
117 CALL PUSHREAL8(p1d(k))
118 p1d(k) = pcps(i, k, j)
119 CALL PUSHREAL8(dz1d(k))
120 dz1d(k) = dz8w(i, k, j)
121 CALL PUSHREAL8(z1d(k))
122 z1d(k) = z(i, k, j)
123 END DO
124 IF (PRESENT(rthcuten) .AND. PRESENT(rqvcuten)) THEN
125 CALL PUSHCONTROL1B(1)
126 ELSE
127 CALL PUSHCONTROL1B(0)
128 END IF
129 END DO
130 END DO
131 qv1db = 0.0_8
132 dthdtb = 0.0_8
133 th1db = 0.0_8
134 rho1db = 0.0_8
135 p1db = 0.0_8
136 z1db = 0.0_8
137 dz1db = 0.0_8
138 t1db = 0.0_8
139 dqvdtb = 0.0_8
140 DO j=j_end,j_start,-1
141 DO i=i_end,i_start,-1
142 CALL POPCONTROL1B(branch)
143 IF (branch .EQ. 0) THEN
144 pprateb = 0.0_8
145 ELSE
146 pprateb = pratecb(i, j) + dt*raincvb(i, j)
147 raincvb(i, j) = 0.0_8
148 pratecb(i, j) = 0.0_8
149 DO k=kte,kts,-1
150 dqvdtb(k) = dqvdtb(k) + rqvcutenb(i, k, j)
151 rqvcutenb(i, k, j) = 0.0_8
152 dthdtb(k) = dthdtb(k) + rthcutenb(i, k, j)
153 rthcutenb(i, k, j) = 0.0_8
154 END DO
155 END IF
156 CALL DUCU1D_B(i, j, u1d, v1d, t1d, t1db, qv1d, qv1db, p1d, p1db&
157 & , dz1d, dz1db, z1d, z1db, w0avg1d, dt, dx, dxsq, rho1d, &
158 & rho1db, th1d, th1db, xlv, cp, rd, rv, g, ep2, svp1, svp2&
159 & , svp3, svpt0, dqvdt, dqvdtb, dthdt, dthdtb, pprate, &
160 & pprateb, nca, ntst, cutop, cubot, ids, ide, jds, jde, &
161 & kds, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, &
162 & jte, kts, kte)
163 DO k=kte,kts,-1
164 CALL POPREAL8(z1d(k))
165 zb(i, k, j) = zb(i, k, j) + z1db(k)
166 z1db(k) = 0.0_8
167 CALL POPREAL8(dz1d(k))
168 dz8wb(i, k, j) = dz8wb(i, k, j) + dz1db(k)
169 dz1db(k) = 0.0_8
170 CALL POPREAL8(p1d(k))
171 pcpsb(i, k, j) = pcpsb(i, k, j) + p1db(k)
172 p1db(k) = 0.0_8
173 CALL POPCONTROL1B(branch)
174 IF (branch .EQ. 0) qv1db(k) = 0.0_8
175 CALL POPREAL8(qv1d(k))
176 qvb(i, k, j) = qvb(i, k, j) + qv1db(k)
177 qv1db(k) = 0.0_8
178 CALL POPREAL8(rho1d(k))
179 rhob(i, k, j) = rhob(i, k, j) + rho1db(k)
180 rho1db(k) = 0.0_8
181 CALL POPREAL8(th1d(k))
182 thb(i, k, j) = thb(i, k, j) + th1db(k)
183 th1db(k) = 0.0_8
184 CALL POPREAL8(t1d(k))
185 tb(i, k, j) = tb(i, k, j) + t1db(k)
186 t1db(k) = 0.0_8
187 END DO
188 pratecb(i, j) = 0.0_8
189 raincvb(i, j) = 0.0_8
190 DO k=kte,kts,-1
191 dthdtb(k) = 0.0_8
192 dqvdtb(k) = 0.0_8
193 END DO
194 END DO
195 END DO
196 END IF
197 END SUBROUTINE DUCU_B
199 ! Generated by TAPENADE (INRIA, Tropics team)
200 ! Tapenade 3.7 (r4786) - 21 Feb 2013 15:53
201 !
202 ! Differentiation of ducu1d in reverse (adjoint) mode (with options r8):
203 ! gradient of useful results: dthdt p0 z t0 th0 pprate rhoe
204 ! dzq dqvdt qv0
205 ! with respect to varying inputs: dthdt p0 z t0 th0 rhoe dzq
206 ! dqvdt qv0
207 ! ****************************************************************************
208 !-----------------------------------------------------------
209 SUBROUTINE DUCU1D_B(i, j, u0, v0, t0, t0b, qv0, qv0b, p0, p0b, dzq, dzqb&
210 & , z, zb0, w0avg1d, delt, dx, dxsq, rhoe, rhoeb, th0, th0b, xlv, cp, rd&
211 & , rv, g, ep2, svp1, svp2, svp3, svpt0, dqvdt, dqvdtb, dthdt, dthdtb, &
212 & pprate, pprateb, nca, ntst, cutop, cubot, ids, ide, jds, jde, kds, kde&
213 & , ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
214 IMPLICIT NONE
215 !-----------------------------------------------------------------------
216 !-----------------------------------------------------------
217 INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
218 & jme, kms, kme, its, ite, jts, jte, kts, kte, i, j, ntst
219 !
220 REAL, DIMENSION(kts:kte), INTENT(IN) :: u0, v0, t0, th0, qv0, p0, rhoe&
221 & , dzq, z, w0avg1d
222 REAL, DIMENSION(kts:kte) :: t0b, th0b, qv0b, p0b, rhoeb, dzqb, zb0
223 !
224 REAL, INTENT(IN) :: delt, dx, dxsq
225 !
226 REAL, INTENT(IN) :: xlv, cp, rd, rv, g
227 REAL, INTENT(IN) :: ep2, svp1, svp2, svp3, svpt0
228 !
229 REAL, DIMENSION(kts:kte), INTENT(INOUT) :: dqvdt, dthdt
230 REAL, DIMENSION(kts:kte), INTENT(INOUT) :: dqvdtb
231 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
232 REAL, DIMENSION(ims:ime, jms:jme), INTENT(OUT) :: cubot, cutop
233 REAL :: pprate
234 REAL :: pprateb
235 !
236 !...DEFINE LOCAL VARIABLES...
237 !
238 REAL, DIMENSION(kts:kte) :: cond, h, hs, qs, x
239 REAL, DIMENSION(kts:kte) :: condb, hb, hsb, qsb, xb
240 REAL :: buoy, cape, cin, dh, dq, dt, dtm, ep, es, evap, hp, mp, qp, &
241 & qsp, rrk, rrkp, tadp, tdp, zb, zg, zi, zt
242 REAL :: dhb, dqb, dtb, esb, evapb, hpb, mpb, qpb, qspb, rrkb, rrkpb
243 INTEGER :: ipos, isat, k, kb, ki, kt
244 INTEGER :: branch
245 INTEGER :: ad_count
246 INTEGER :: i0
247 INTEGER :: ad_from
248 INTEGER :: ad_to
249 INTEGER :: ad_from0
250 REAL, DIMENSION(kts:kte), INTENT(INOUT) :: dthdtb
251 REAL :: temp3
252 REAL :: temp2
253 REAL :: temp1
254 REAL :: temp0
255 REAL :: temp7b
256 REAL :: temp9b1
257 REAL :: temp9b0
258 REAL :: temp7b0
259 REAL :: temp6b
260 REAL :: temp12
261 REAL :: temp11
262 REAL :: temp10
263 REAL :: temp9b
264 REAL :: temp5b0
265 REAL :: tempb
266 REAL :: temp2b
267 REAL :: temp5b
268 REAL :: temp8b0
269 REAL :: temp8b
270 REAL :: abs1
271 REAL :: abs0
272 REAL :: temp1b
273 REAL :: temp
274 REAL :: temp6b0
275 REAL :: temp9
276 REAL :: temp1b1
277 REAL :: temp8
278 REAL :: temp1b0
279 REAL :: temp7
280 REAL :: temp10b
281 REAL :: temp6
282 REAL :: temp4b
283 REAL :: temp5
284 REAL :: temp10b0
285 REAL :: temp4
286 !
287 !...DEFINE PROFILES
288 DO k=kts,kte
289 h(k) = cp*t0(k) + g*z(k) + xlv*qv0(k)
290 CALL PUSHREAL8(es)
291 es = 1000.*svp1*EXP(svp2*(t0(k)-svpt0)/(t0(k)-svp3))
292 qs(k) = ep2*es/(p0(k)-es)
293 hs(k) = cp*t0(k) + g*z(k) + xlv*qs(k)
294 x(k) = xlv*xlv*qs(k)/(cp*rv*t0(k)*t0(k))
295 END DO
296 !
297 !...LOOP OVER PARCELS
298 loop_origin:DO ki=kts,kte
299 hp = h(ki)
300 qp = qv0(ki)
301 CALL PUSHREAL8(mp)
302 mp = alpha*rhoe(ki)*dzq(ki)
303 buoy = 0.
304 cape = 0.
305 dtm = 0.
306 isat = 0
307 ipos = 0
308 CALL PUSHINTEGER4(kt)
309 kt = 0
310 kb = 0
311 zt = 0.
312 zb = 0.
313 cond = 0.
314 CALL PUSHINTEGER4(k)
315 ad_count = 1
316 !
317 !...LIFT PARCEL
318 loop_lift:DO k=ki+1,kte
319 tadp = t0(ki) + g/cp*(z(ki)-z(k))
320 ep = p0(k)*qv0(ki)/(ep2+qv0(ki))
321 tdp = (svpt0-svp3/svp2*ALOG(0.001*ep/svp1))/(1.-1./svp2*ALOG(0.001&
322 & *ep/svp1))
323 IF (tadp .GE. tdp) THEN
324 dt = tadp - t0(k)
325 cond(k) = 0.
326 CALL PUSHCONTROL1B(1)
327 ELSE
328 ! saturated
329 IF (isat .EQ. 0) THEN
330 kb = k
331 zb = z(k) - 0.5*dzq(k)
332 END IF
333 isat = 1
334 CALL PUSHREAL8(dh)
335 dh = hp - hs(k)
336 dt = dh/cp/(1.+x(k))
337 CALL PUSHREAL8(qsp)
338 qsp = qs(k) + dh/xlv*x(k)/(1.+x(k))
339 !...CONDENSATE PRODUCED
340 cond(k) = mp*(qp-qsp)
341 CALL PUSHREAL8(qp)
342 qp = qsp
343 CALL PUSHCONTROL1B(0)
344 END IF
345 buoy = buoy + g*dt*dzq(k)/t0(k)
346 IF (cape .LT. buoy) THEN
347 cape = buoy
348 ELSE
349 cape = cape
350 END IF
351 IF (dt .GE. 0.) THEN
352 kt = k
353 zt = z(k) + 0.5*dzq(k)
354 ELSE IF (dt .LT. 0. .AND. dtm .GE. 0.) THEN
355 IF (dt .GE. 0.) THEN
356 abs0 = dt
357 ELSE
358 abs0 = -dt
359 END IF
360 IF (dtm .GE. 0.) THEN
361 abs1 = dtm
362 ELSE
363 abs1 = -dtm
364 END IF
365 ! cloud top is level closest to parcel temperature
366 IF (abs0 .LT. abs1) THEN
367 kt = k
368 zt = z(k) + 0.5*dzq(k)
369 END IF
370 END IF
371 dtm = dt
372 ! continue lifting until buoyancy is gone
373 IF (buoy .LT. cincap) THEN
374 GOTO 100
375 ELSE
376 IF (buoy .GT. 0.) ipos = 1
377 ! positive area detected
378 IF (k .EQ. 1) THEN
379 kt = k
380 zt = z(k) + 0.5*dzq(k)
381 END IF
382 CALL PUSHINTEGER4(k)
383 ad_count = ad_count + 1
384 END IF
385 END DO loop_lift
386 CALL PUSHCONTROL1B(0)
387 CALL PUSHINTEGER4(ad_count)
388 GOTO 110
389 100 CALL PUSHCONTROL1B(1)
390 CALL PUSHINTEGER4(ad_count)
391 !
392 !...CHECK FOR CLOUD
393 110 IF (isat .EQ. 0) THEN
394 CALL PUSHCONTROL3B(3)
395 ELSE IF (zt - zb .LE. dpthmin) THEN
396 ! no cloud from lifting - no convection
397 CALL PUSHCONTROL3B(2)
398 ELSE IF (ipos .EQ. 0) THEN
399 ! not more than one cloud level - no convection
400 CALL PUSHCONTROL3B(1)
401 ELSE IF (cape .LE. capemin) THEN
402 ! no buoyancy in cloud - no convection
403 CALL PUSHCONTROL3B(0)
404 ELSE
405 CALL PUSHREAL8(dh)
406 ! not enough cape
407 !
408 !...IF CHECK FOR CLOUD SUCCESSFUL
409 !
410 !...DETRAINMENT
411 ad_from = kt - 1
412 CALL PUSHINTEGER4(k)
413 !
414 !...SUBSIDENCE
415 k = ki - 1
416 CALL PUSHINTEGER4(k + 1)
417 CALL PUSHINTEGER4(ad_from)
418 !
419 !...RAINFALL AND EVAPORATION
420 rrkp = 0.
421 ad_from0 = kt
422 loop_rainfall:DO k=ad_from0,1,-1
423 rrk = rrkp + cond(k)
424 CALL PUSHREAL8(evap)
425 evap = dzq(k)*rrkp/vfall*eps*(qs(k)-qv0(k))
426 ! restrict evap to below cloud base
427 IF (k .GE. kb) THEN
428 evap = 0.
429 CALL PUSHCONTROL1B(1)
430 ELSE
431 CALL PUSHCONTROL1B(0)
432 END IF
433 IF (rrk .GT. evap) THEN
434 CALL PUSHCONTROL1B(0)
435 evap = evap
436 ELSE
437 evap = rrk
438 CALL PUSHCONTROL1B(1)
439 END IF
440 rrk = rrk - evap
441 CALL PUSHREAL8(rrkp)
442 rrkp = rrk
443 END DO loop_rainfall
444 CALL PUSHINTEGER4(ad_from0)
445 CALL PUSHCONTROL3B(4)
446 END IF
447 END DO loop_origin
448 hb = 0.0_8
449 qsb = 0.0_8
450 xb = 0.0_8
451 hsb = 0.0_8
452 DO ki=kte,kts,-1
453 CALL POPCONTROL3B(branch)
454 IF (branch .LT. 2) THEN
455 IF (branch .EQ. 0) THEN
456 hpb = 0.0_8
457 condb = 0.0_8
458 mpb = 0.0_8
459 ELSE
460 hpb = 0.0_8
461 condb = 0.0_8
462 mpb = 0.0_8
463 END IF
464 ELSE IF (branch .EQ. 2) THEN
465 hpb = 0.0_8
466 condb = 0.0_8
467 mpb = 0.0_8
468 ELSE IF (branch .EQ. 3) THEN
469 hpb = 0.0_8
470 condb = 0.0_8
471 mpb = 0.0_8
472 ELSE
473 rrkpb = pprateb
474 condb = 0.0_8
475 CALL POPINTEGER4(ad_from0)
476 DO k=1,ad_from0,1
477 temp9 = rhoe(k)*dzq(k)
478 temp9b1 = -(evap*dqvdtb(k)/temp9**2)
479 CALL POPREAL8(rrkp)
480 rrkb = rrkpb
481 temp12 = rhoe(k)*dzq(k)
482 temp11 = cp*t0(kt)
483 temp10 = temp11*temp12
484 temp10b = -(xlv*dthdtb(k)/temp10)
485 temp10b0 = -(th0(kt)*evap*temp10b/temp10)
486 th0b(kt) = th0b(kt) + evap*temp10b
487 evapb = dqvdtb(k)/temp9 - rrkb + th0(kt)*temp10b
488 t0b(kt) = t0b(kt) + temp12*cp*temp10b0
489 rhoeb(k) = rhoeb(k) + dzq(k)*temp9b1 + temp11*dzq(k)*temp10b0
490 dzqb(k) = dzqb(k) + rhoe(k)*temp9b1 + temp11*rhoe(k)*temp10b0
491 CALL POPCONTROL1B(branch)
492 IF (branch .NE. 0) THEN
493 rrkb = rrkb + evapb
494 evapb = 0.0_8
495 END IF
496 CALL POPCONTROL1B(branch)
497 IF (branch .NE. 0) evapb = 0.0_8
498 CALL POPREAL8(evap)
499 temp9b = eps*dzq(k)*rrkp*evapb/vfall
500 temp9b0 = eps*(qs(k)-qv0(k))*evapb
501 qsb(k) = qsb(k) + temp9b
502 qv0b(k) = qv0b(k) - temp9b
503 dzqb(k) = dzqb(k) + rrkp*temp9b0/vfall
504 rrkpb = rrkb + dzq(k)*temp9b0/vfall
505 condb(k) = condb(k) + rrkb
506 END DO
507 mp = alpha*rhoe(ki)*dzq(ki)
508 mpb = 0.0_8
509 CALL POPINTEGER4(ad_from)
510 CALL POPINTEGER4(ad_to)
511 DO k=ad_to,ad_from,1
512 temp7 = rhoe(k)*dzq(k)
513 temp7b = dthdtb(k)/temp7
514 temp7b0 = -(mp*(th0(k+1)-th0(k))*temp7b/temp7)
515 temp8 = rhoe(k)*dzq(k)
516 temp8b = dqvdtb(k)/temp8
517 temp8b0 = -(mp*(qv0(k+1)-qv0(k))*temp8b/temp8)
518 mpb = mpb + (th0(k+1)-th0(k))*temp7b + (qv0(k+1)-qv0(k))*temp8b
519 qv0b(k+1) = qv0b(k+1) + mp*temp8b
520 qv0b(k) = qv0b(k) - mp*temp8b
521 rhoeb(k) = rhoeb(k) + dzq(k)*temp7b0 + dzq(k)*temp8b0
522 dzqb(k) = dzqb(k) + rhoe(k)*temp7b0 + rhoe(k)*temp8b0
523 th0b(k+1) = th0b(k+1) + mp*temp7b
524 th0b(k) = th0b(k) - mp*temp7b
525 END DO
526 CALL POPINTEGER4(k)
527 temp3 = cp*(x(kt)+1.)
528 temp5 = t0(kt)*rhoe(kt)
529 temp5b = dthdtb(kt)/(temp5*dzq(kt))
530 hp = h(ki)
531 dh = hp - hs(kt)
532 dq = qs(kt) + dh/xlv*x(kt)/(1.+x(kt)) - qv0(kt)
533 temp6 = rhoe(kt)*dzq(kt)
534 temp6b = dqvdtb(kt)/temp6
535 temp6b0 = -(mp*dq*temp6b/temp6)
536 dqb = mp*temp6b
537 dt = dh/cp/(1.+x(kt))
538 mpb = mpb + th0(kt)*dt*temp5b + dq*temp6b
539 temp5b0 = -(th0(kt)*mp*dt*temp5b/(temp5*dzq(kt)))
540 rhoeb(kt) = rhoeb(kt) + dzq(kt)*t0(kt)*temp5b0 + dzq(kt)*temp6b0
541 dzqb(kt) = dzqb(kt) + temp5*temp5b0 + rhoe(kt)*temp6b0
542 th0b(kt) = th0b(kt) + mp*dt*temp5b
543 dtb = th0(kt)*mp*temp5b
544 t0b(kt) = t0b(kt) + dzq(kt)*rhoe(kt)*temp5b0
545 temp4 = xlv*(x(kt)+1.)
546 temp4b = dqb/temp4
547 qsb(kt) = qsb(kt) + dqb
548 dhb = dtb/temp3 + x(kt)*temp4b
549 xb(kt) = xb(kt) + (dh-dh*x(kt)*xlv/temp4)*temp4b - dh*cp*dtb/temp3&
550 & **2
551 qv0b(kt) = qv0b(kt) - dqb
552 CALL POPREAL8(dh)
553 hpb = dhb
554 hsb(kt) = hsb(kt) - dhb
555 END IF
556 CALL POPINTEGER4(ad_count)
557 DO i0=1,ad_count
558 IF (i0 .EQ. 1) THEN
559 CALL POPCONTROL1B(branch)
560 IF (branch .EQ. 0) THEN
561 qpb = 0.0_8
562 GOTO 120
563 ELSE
564 qpb = 0.0_8
565 END IF
566 END IF
567 CALL POPCONTROL1B(branch)
568 IF (branch .EQ. 0) THEN
569 CALL POPREAL8(qp)
570 qspb = qpb - mp*condb(k)
571 mpb = mpb + (qp-qsp)*condb(k)
572 qpb = mp*condb(k)
573 condb(k) = 0.0_8
574 CALL POPREAL8(qsp)
575 temp2 = xlv*(x(k)+1.)
576 temp2b = qspb/temp2
577 qsb(k) = qsb(k) + qspb
578 dhb = x(k)*temp2b
579 xb(k) = xb(k) + (dh-dh*x(k)*xlv/temp2)*temp2b
580 CALL POPREAL8(dh)
581 hpb = hpb + dhb
582 hsb(k) = hsb(k) - dhb
583 ELSE
584 condb(k) = 0.0_8
585 END IF
586 120 CALL POPINTEGER4(k)
587 END DO
588 CALL POPINTEGER4(kt)
589 CALL POPREAL8(mp)
590 rhoeb(ki) = rhoeb(ki) + alpha*dzq(ki)*mpb
591 dzqb(ki) = dzqb(ki) + alpha*rhoe(ki)*mpb
592 qv0b(ki) = qv0b(ki) + qpb
593 hb(ki) = hb(ki) + hpb
594 END DO
595 DO k=kte,kts,-1
596 dqvdtb(k) = 0.0_8
597 dthdtb(k) = 0.0_8
598 temp1 = cp*rv*t0(k)**2
599 temp1b = xlv**2*xb(k)/temp1
600 qsb(k) = qsb(k) + xlv*hsb(k) + temp1b
601 xb(k) = 0.0_8
602 zb0(k) = zb0(k) + g*hb(k) + g*hsb(k)
603 temp1b0 = ep2*qsb(k)/(p0(k)-es)
604 temp1b1 = -(es*temp1b0/(p0(k)-es))
605 esb = temp1b0 - temp1b1
606 p0b(k) = p0b(k) + temp1b1
607 qsb(k) = 0.0_8
608 CALL POPREAL8(es)
609 temp0 = t0(k) - svp3
610 temp = (t0(k)-svpt0)/temp0
611 tempb = svp2*EXP(svp2*temp)*svp1*1000.*esb/temp0
612 t0b(k) = t0b(k) + cp*hsb(k) + cp*hb(k) + (1.0-temp)*tempb - cp*rv*qs&
613 & (k)*2*t0(k)*temp1b/temp1
614 hsb(k) = 0.0_8
615 qv0b(k) = qv0b(k) + xlv*hb(k)
616 hb(k) = 0.0_8
617 END DO
618 END SUBROUTINE DUCU1D_B
620 ! Generated by TAPENADE (INRIA, Tropics team)
621 ! Tapenade 3.6 (r4343) - 10 Feb 2012 10:52
622 !
623 ! Differentiation of ducuinit in reverse (adjoint) mode:
624 ! gradient of useful results: rqccuten rthcuten w0avg rqvcuten
625 ! with respect to varying inputs: rqccuten rthcuten w0avg rqvcuten
626 SUBROUTINE DUCUINIT_B(rthcuten, rthcutenb, rqvcuten, rqvcutenb, rqccuten&
627 & , rqccutenb, rqrcuten, rqicuten, rqscuten, nca, w0avg, w0avgb, p_qc, &
628 & p_qr, svp1, svp2, svp3, svpt0, p_first_scalar, restart, &
629 & allowed_to_read, ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms&
630 & , kme, its, ite, jts, jte, kts, kte)
631 IMPLICIT NONE
632 !--------------------------------------------------------------------
633 LOGICAL, INTENT(IN) :: restart, allowed_to_read
634 INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
635 & jme, kms, kme, its, ite, jts, jte, kts, kte
636 INTEGER, INTENT(IN) :: p_qc, p_qr, p_first_scalar
637 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: rthcuten, rqvcuten, &
638 & rqccuten, rqrcuten, rqicuten, rqscuten
639 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: rthcutenb, rqvcutenb, &
640 & rqccutenb
641 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: w0avg
642 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: w0avgb
643 REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
644 INTEGER :: i, j, k, itf, jtf, ktf
645 REAL, INTENT(IN) :: svp1, svp2, svp3, svpt0
646 INTEGER :: branch
647 INTRINSIC MIN0
648 IF (jte .GT. jde - 1) THEN
649 jtf = jde - 1
650 ELSE
651 jtf = jte
652 END IF
653 IF (kte .GT. kde - 1) THEN
654 ktf = kde - 1
655 ELSE
656 ktf = kte
657 END IF
658 IF (ite .GT. ide - 1) THEN
659 itf = ide - 1
660 ELSE
661 itf = ite
662 END IF
663 IF (.NOT.restart) THEN
664 IF (p_qc .GE. p_first_scalar) THEN
665 CALL PUSHCONTROL1B(1)
666 ELSE
667 CALL PUSHCONTROL1B(0)
668 END IF
669 DO j=jtf,jts,-1
670 DO k=ktf,kts,-1
671 DO i=itf,its,-1
672 w0avgb(i, k, j) = 0.0
673 END DO
674 END DO
675 END DO
676 CALL POPCONTROL1B(branch)
677 IF (branch .NE. 0) THEN
678 DO j=jtf,jts,-1
679 DO k=ktf,kts,-1
680 DO i=itf,its,-1
681 rqccutenb(i, k, j) = 0.0
682 END DO
683 END DO
684 END DO
685 END IF
686 DO j=jtf,jts,-1
687 DO k=ktf,kts,-1
688 DO i=itf,its,-1
689 rqvcutenb(i, k, j) = 0.0
690 rthcutenb(i, k, j) = 0.0
691 END DO
692 END DO
693 END DO
694 END IF
695 END SUBROUTINE DUCUINIT_B
697 END MODULE a_module_cu_du