FV3 Bundle
fv_sg_adm.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU General Public License *
3 !* This file is a part of fvGFS. *
4 !* *
5 !* fvGFS is free software; you can redistribute it and/or modify it *
6 !* and are expected to follow the terms of the GNU General Public *
7 !* License as published by the Free Software Foundation; either *
8 !* version 2 of the License, or (at your option) any later version. *
9 !* *
10 !* fvGFS is distributed in the hope that it will be useful, but *
11 !* WITHOUT ANY WARRANTY; without even the implied warranty of *
12 !* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
13 !* General Public License for more details. *
14 !* *
15 !* For the full text of the GNU General Public License, *
16 !* write to: Free Software Foundation, Inc., *
17 !* 675 Mass Ave, Cambridge, MA 02139, USA. *
18 !* or see: http://www.gnu.org/licenses/gpl.html *
19 !***********************************************************************
20 module fv_sg_adm_mod
21 
22 !-----------------------------------------------------------------------
23 ! FV sub-grid mixing
24 !-----------------------------------------------------------------------
25  use constants_mod, only: rdgas, rvgas, cp_air, cp_vapor, hlv, hlf, kappa, grav
26  use tracer_manager_mod, only: get_tracer_index
27  use field_manager_mod, only: model_atmos
28  use lin_cld_microphys_mod, only: wqs2, wqsat2_moist
29  use fv_mp_nlm_mod, only: mp_reduce_min, is_master
30  use tapenade_iter, only: pushcontrol, popcontrol, pushinteger, popinteger, &
31  pushrealarray, poprealarray, pushrealarray_adm, poprealarray_adm
32 
33 implicit none
34 private
35 
36 public fv_subgrid_z_fwd, fv_subgrid_z_bwd
37 public fv_subgrid_z, neg_adj3
38 
39  real, parameter:: esl = 0.621971831
40  real, parameter:: tice = 273.16
41 ! real, parameter:: c_ice = 2106. ! Emanuel table, page 566
42  real, parameter:: c_ice = 1972. ! -15 C
43  real, parameter:: c_liq = 4.1855e+3 ! GFS
44 ! real, parameter:: c_liq = 4218. ! ECMWF-IFS
45  real, parameter:: cv_vap = cp_vapor - rvgas ! 1384.5
46  real, parameter:: c_con = c_ice
47 
48 ! real, parameter:: dc_vap = cp_vapor - c_liq ! = -2368.
49  real, parameter:: dc_vap = cv_vap - c_liq ! = -2368.
50  real, parameter:: dc_ice = c_liq - c_ice ! = 2112.
51 ! Values at 0 Deg C
52  real, parameter:: hlv0 = 2.5e6
53  real, parameter:: hlf0 = 3.3358e5
54 ! real, parameter:: hlv0 = 2.501e6 ! Emanual Appendix-2
55 ! real, parameter:: hlf0 = 3.337e5 ! Emanual
56  real, parameter:: t_ice = 273.16
57  real, parameter:: ri_max = 1.
58  real, parameter:: ri_min = 0.25
59  real, parameter:: t1_min = 160.
60  real, parameter:: t2_min = 165.
61  real, parameter:: t2_max = 315.
62  real, parameter:: t3_max = 325.
63  real, parameter:: lv0 = hlv0 - dc_vap*t_ice ! = 3.147782e6
64  real, parameter:: li0 = hlf0 - dc_ice*t_ice ! = -2.431928e5
65 
66  real, parameter:: zvir = rvgas/rdgas - 1. ! = 0.607789855
67  real, allocatable:: table(:),des(:)
68  real:: lv00, d0_vap
69 
70 CONTAINS
71 ! Differentiation of fv_subgrid_z in reverse (adjoint) mode, forward sweep (with options split(a2b_edge_mod.a2b_ord4 a2b_edge
72 !_mod.a2b_ord2 dyn_core_mod.dyn_core dyn_core_mod.pk3_halo dyn_core_mod.pln_halo dyn_core_mod.pe_halo dyn_core_mod.adv_pe dyn_cor
73 !e_mod.p_grad_c dyn_core_mod.nh_p_grad dyn_core_mod.split_p_grad dyn_core_mod.one_grad_p dyn_core_mod.grad1_p_update dyn_core_mod
74 !.mix_dp dyn_core_mod.geopk dyn_core_mod.del2_cubed dyn_core_mod.Rayleigh_fast fv_dynamics_mod.fv_dynamics fv_dynamics_mod.Raylei
75 !gh_Super fv_dynamics_mod.Rayleigh_Friction fv_dynamics_mod.compute_aam fv_grid_utils_mod.cubed_to_latlon fv_grid_utils_mod.c2l_o
76 !rd4 fv_grid_utils_mod.c2l_ord2 fv_mapz_mod.Lagrangian_to_Eulerian fv_mapz_mod.compute_total_energy fv_mapz_mod.pkez fv_mapz_mod.
77 !remap_z fv_mapz_mod.map_scalar fv_mapz_mod.map1_ppm fv_mapz_mod.mapn_tracer fv_mapz_mod.map1_q2 fv_mapz_mod.remap_2d
78 ! fv_mapz_mod.scalar_profile fv_mapz_mod.cs_profile fv_mapz_mod.cs_limiters fv_mapz_mod.ppm_profile fv_mapz_mod.ppm_limiter
79 !s fv_mapz_mod.steepz fv_mapz_mod.rst_remap fv_mapz_mod.mappm fv_mapz_mod.moist_cv fv_mapz_mod.moist_cp fv_mapz_mod.map1_cubic fv
80 !_restart_mod.d2c_setup fv_tracer2d_mod.tracer_2d_1L fv_tracer2d_mod.tracer_2d fv_tracer2d_mod.tracer_2d_nested fv_sg_mod.fv_subg
81 !rid_z main_mod.compute_pressures main_mod.run nh_core_mod.Riem_Solver3 nh_utils_mod.update_dz_c nh_utils_mod.update_dz_d nh_util
82 !s_mod.Riem_Solver_c nh_utils_mod.Riem_Solver3test nh_utils_mod.imp_diff_w nh_utils_mod.RIM_2D nh_utils_mod.SIM3_solver nh_utils_
83 !mod.SIM3p0_solver nh_utils_mod.SIM1_solver nh_utils_mod.SIM_solver nh_utils_mod.edge_scalar nh_utils_mod.edge_profile nh_utils_m
84 !od.nest_halo_nh sw_core_mod.c_sw sw_core_mod.d_sw sw_core_mod.divergence_corner sw_core_mod.divergence_corner_nest sw_core_mod.
85 !d2a2c_vect sw_core_mod.fill3_4corners sw_core_mod.fill2_4corners sw_core_mod.fill_4corners sw_core_mod.xtp_u sw_core_mod.ytp_
86 !v_fb sw_core_mod.compute_divergence_damping sw_core_mod.smag_corner tp_core_mod.mp_ghost_ew tp_core_mod.fv_tp_2d tp_cor
87 !e_mod.copy_corners tp_core_mod.xppm tp_core_mod.yppm tp_core_mod.deln_flux a2b_edge_mod.extrap_corner fv_grid_uti
88 !ls_mod.great_circle_dist sw_core_mod.edge_interpolate4)):
89 ! gradient of useful results: qa w delp delz pkz ta
90 ! with respect to varying inputs: qa peln w delp ua delz va pkz
91 ! pe ta
92  SUBROUTINE fv_subgrid_z_fwd(isd, ied, jsd, jed, is, ie, js, je, km, nq&
93 & , dt, tau, nwat, delp, pe, peln, pkz, ta, qa, ua, va, hydrostatic, w&
94 & , delz, u_dt, v_dt, t_dt, k_bot)
95  IMPLICIT NONE
96 ! Dry convective adjustment-mixing
97 !-------------------------------------------
98  INTEGER, INTENT(IN) :: is, ie, js, je, km, nq, nwat
99  INTEGER, INTENT(IN) :: isd, ied, jsd, jed
100 ! Relaxation time scale
101  INTEGER, INTENT(IN) :: tau
102 ! model time step
103  REAL, INTENT(IN) :: dt
104  REAL, INTENT(IN) :: pe(is-1:ie+1, km+1, js-1:je+1)
105  REAL, INTENT(IN) :: peln(is:ie, km+1, js:je)
106 ! Delta p at each model level
107  REAL, INTENT(IN) :: delp(isd:ied, jsd:jed, km)
108 ! Delta z at each model level
109  REAL, INTENT(IN) :: delz(isd:, jsd:, :)
110  REAL, INTENT(IN) :: pkz(is:ie, js:je, km)
111  LOGICAL, INTENT(IN) :: hydrostatic
112  INTEGER, INTENT(IN), OPTIONAL :: k_bot
113 !
114  REAL, INTENT(INOUT) :: ua(isd:ied, jsd:jed, km)
115  REAL, INTENT(INOUT) :: va(isd:ied, jsd:jed, km)
116  REAL, INTENT(INOUT) :: w(isd:, jsd:, :)
117 ! Temperature
118  REAL, INTENT(INOUT) :: ta(isd:ied, jsd:jed, km)
119 ! Specific humidity & tracers
120  REAL, INTENT(INOUT) :: qa(isd:ied, jsd:jed, km, nq)
121  REAL, INTENT(INOUT) :: u_dt(isd:ied, jsd:jed, km)
122  REAL, INTENT(INOUT) :: v_dt(isd:ied, jsd:jed, km)
123  REAL, INTENT(INOUT) :: t_dt(is:ie, js:je, km)
124 !---------------------------Local variables-----------------------------
125  REAL, DIMENSION(is:ie, km) :: u0, v0, w0, t0, hd, te, gz, tvm, pm, &
126 & den
127  REAL :: q0(is:ie, km, nq), qcon(is:ie, km)
128  REAL, DIMENSION(is:ie) :: gzh, lcp2, icp2, cvm, cpm, qs
129  REAL :: ri_ref, ri, pt1, pt2, ratio, tv, cv, tmp, q_liq, q_sol
130  REAL :: tv1, tv2, g2, h0, mc, fra, rk, rz, rdt, tvd, tv_surf
131  REAL :: dh, dq, qsw, dqsdt, tcp3, t_max, t_min
132  INTEGER :: i, j, k, kk, n, m, iq, km1, im, kbot
133  REAL, PARAMETER :: ustar2=1.e-4
134  REAL :: cv_air, xvir
135  INTEGER :: sphum, liq_wat, rainwat, snowwat, graupel, ice_wat, &
136 & cld_amt
137  INTRINSIC present
138  INTRINSIC min
139  INTRINSIC real
140  INTRINSIC max
141  INTEGER :: min1
142  REAL :: max1
143  REAL :: max2
144  INTEGER :: ad_from
145  REAL :: y1
146 ! = rdgas * (7/2-1) = 2.5*rdgas=717.68
147  cv_air = cp_air - rdgas
148  rk = cp_air/rdgas + 1.
149  g2 = 0.5*grav
150  IF (PRESENT(k_bot)) THEN
151  IF (k_bot .LT. 3) THEN
152  CALL pushcontrol(1,0)
153  GOTO 100
154  ELSE
155  CALL pushcontrol(1,0)
156  kbot = k_bot
157  END IF
158  ELSE
159  CALL pushcontrol(1,1)
160  kbot = km
161  END IF
162  IF (pe(is, 1, js) .LT. 2.) THEN
163  t_min = t1_min
164  ELSE
165  t_min = t2_min
166  END IF
167  IF (km .GT. 24) THEN
168  min1 = 24
169  ELSE
170  min1 = km
171  END IF
172  IF (k_bot .LT. min1) THEN
173  t_max = t2_max
174  ELSE
175  t_max = t3_max
176  END IF
177  sphum = 1
178  rainwat = -1
179  snowwat = -1
180  graupel = -1
181  IF (nwat .EQ. 0) THEN
182  CALL pushcontrol(1,0)
183  xvir = 0.
184  rz = 0.
185  ELSE
186  xvir = zvir
187 ! rz = zvir * rdgas
188  rz = rvgas - rdgas
189  IF (nwat .EQ. 3) THEN
190  CALL pushcontrol(1,1)
191  liq_wat = 2
192  ice_wat = 3
193  ELSE
194  CALL pushcontrol(1,1)
195  END IF
196  END IF
197 !------------------------------------------------------------------------
198 ! The nonhydrostatic pressure changes if there is heating (under constant
199 ! volume and mass is locally conserved).
200 !------------------------------------------------------------------------
201  m = 3
202  fra = dt/REAL(tau)
203 !$OMP parallel do default(none) shared(im,is,ie,js,je,nq,kbot,qa,ta,sphum,ua,va,delp,peln, &
204 !$OMP hydrostatic,pe,delz,g2,w,liq_wat,rainwat,ice_wat, &
205 !$OMP snowwat,cv_air,m,graupel,pkz,rk,rz,fra, t_max, t_min, &
206 !$OMP u_dt,rdt,v_dt,xvir,nwat) &
207 !$OMP private(kk,lcp2,icp2,tcp3,dh,dq,den,qs,qsw,dqsdt,qcon,q0, &
208 !$OMP t0,u0,v0,w0,h0,pm,gzh,tvm,tmp,cpm,cvm,q_liq,q_sol, &
209 !$OMP tv,gz,hd,te,ratio,pt1,pt2,tv1,tv2,ri_ref, ri,mc,km1)
210  DO j=js,je
211  DO iq=1,nq
212  DO k=1,kbot
213  DO i=is,ie
214  CALL pushrealarray(q0(i, k, iq))
215  q0(i, k, iq) = qa(i, j, k, iq)
216  END DO
217  END DO
218  END DO
219  DO k=1,kbot
220  DO i=is,ie
221  CALL pushrealarray(t0(i, k))
222  t0(i, k) = ta(i, j, k)
223  tvm(i, k) = t0(i, k)*(1.+xvir*q0(i, k, sphum))
224  CALL pushrealarray(u0(i, k))
225  u0(i, k) = ua(i, j, k)
226  CALL pushrealarray(v0(i, k))
227  v0(i, k) = va(i, j, k)
228  CALL pushrealarray(pm(i, k))
229  pm(i, k) = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
230  END DO
231  END DO
232  DO i=is,ie
233  CALL pushrealarray(gzh(i))
234  gzh(i) = 0.
235  END DO
236  IF (hydrostatic) THEN
237  DO k=kbot,1,-1
238  DO i=is,ie
239  CALL pushrealarray(tv)
240  tv = rdgas*tvm(i, k)
241  CALL pushrealarray(gz(i, k))
242  gz(i, k) = gzh(i) + tv*(1.-pe(i, k, j)/pm(i, k))
243  CALL pushrealarray(hd(i, k))
244  hd(i, k) = cp_air*tvm(i, k) + gz(i, k) + 0.5*(u0(i, k)**2+v0&
245 & (i, k)**2)
246  CALL pushrealarray(gzh(i))
247  gzh(i) = gzh(i) + tv*(peln(i, k+1, j)-peln(i, k, j))
248  END DO
249  END DO
250  CALL pushcontrol(1,1)
251  ELSE
252  DO k=kbot,1,-1
253  IF (nwat .EQ. 0) THEN
254  DO i=is,ie
255  CALL pushrealarray(cpm(i))
256  cpm(i) = cp_air
257  CALL pushrealarray(cvm(i))
258  cvm(i) = cv_air
259  END DO
260  CALL pushcontrol(3,5)
261  ELSE IF (nwat .EQ. 1) THEN
262  DO i=is,ie
263  CALL pushrealarray(cpm(i))
264  cpm(i) = (1.-q0(i, k, sphum))*cp_air + q0(i, k, sphum)*&
265 & cp_vapor
266  CALL pushrealarray(cvm(i))
267  cvm(i) = (1.-q0(i, k, sphum))*cv_air + q0(i, k, sphum)*&
268 & cv_vap
269  END DO
270  CALL pushcontrol(3,4)
271  ELSE IF (nwat .EQ. 2) THEN
272 ! GFS
273  DO i=is,ie
274  CALL pushrealarray(cpm(i))
275  cpm(i) = (1.-q0(i, k, sphum))*cp_air + q0(i, k, sphum)*&
276 & cp_vapor
277  CALL pushrealarray(cvm(i))
278  cvm(i) = (1.-q0(i, k, sphum))*cv_air + q0(i, k, sphum)*&
279 & cv_vap
280  END DO
281  CALL pushcontrol(3,3)
282  ELSE IF (nwat .EQ. 3) THEN
283  DO i=is,ie
284  q_liq = q0(i, k, liq_wat)
285  q_sol = q0(i, k, ice_wat)
286  CALL pushrealarray(cpm(i))
287  cpm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cp_air + q0(i&
288 & , k, sphum)*cp_vapor + q_liq*c_liq + q_sol*c_ice
289  CALL pushrealarray(cvm(i))
290  cvm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cv_air + q0(i&
291 & , k, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
292  END DO
293  CALL pushcontrol(3,2)
294  ELSE IF (nwat .EQ. 4) THEN
295  DO i=is,ie
296  q_liq = q0(i, k, liq_wat) + q0(i, k, rainwat)
297  CALL pushrealarray(cpm(i))
298  cpm(i) = (1.-(q0(i, k, sphum)+q_liq))*cp_air + q0(i, k, &
299 & sphum)*cp_vapor + q_liq*c_liq
300  CALL pushrealarray(cvm(i))
301  cvm(i) = (1.-(q0(i, k, sphum)+q_liq))*cv_air + q0(i, k, &
302 & sphum)*cv_vap + q_liq*c_liq
303  END DO
304  CALL pushcontrol(3,1)
305  ELSE
306  DO i=is,ie
307  q_liq = q0(i, k, liq_wat) + q0(i, k, rainwat)
308  q_sol = q0(i, k, ice_wat) + q0(i, k, snowwat) + q0(i, k, &
309 & graupel)
310  CALL pushrealarray(cpm(i))
311  cpm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cp_air + q0(i&
312 & , k, sphum)*cp_vapor + q_liq*c_liq + q_sol*c_ice
313  CALL pushrealarray(cvm(i))
314  cvm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cv_air + q0(i&
315 & , k, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
316  END DO
317  CALL pushcontrol(3,0)
318  END IF
319  DO i=is,ie
320  CALL pushrealarray(w0(i, k))
321  w0(i, k) = w(i, j, k)
322  CALL pushrealarray(gz(i, k))
323  gz(i, k) = gzh(i) - g2*delz(i, j, k)
324  tmp = gz(i, k) + 0.5*(u0(i, k)**2+v0(i, k)**2+w0(i, k)**2)
325  CALL pushrealarray(hd(i, k))
326  hd(i, k) = cpm(i)*t0(i, k) + tmp
327  CALL pushrealarray(te(i, k))
328  te(i, k) = cvm(i)*t0(i, k) + tmp
329  CALL pushrealarray(gzh(i))
330  gzh(i) = gzh(i) - grav*delz(i, j, k)
331  END DO
332  END DO
333  CALL pushcontrol(1,0)
334  END IF
335  DO n=1,m
336  IF (m .EQ. 3) THEN
337  IF (n .EQ. 1) THEN
338  CALL pushrealarray(ratio)
339  ratio = 0.25
340  CALL pushcontrol(1,0)
341  ELSE
342  CALL pushcontrol(1,1)
343  END IF
344  IF (n .EQ. 2) THEN
345  CALL pushrealarray(ratio)
346  ratio = 0.5
347  CALL pushcontrol(1,0)
348  ELSE
349  CALL pushcontrol(1,1)
350  END IF
351  IF (n .EQ. 3) THEN
352  CALL pushrealarray(ratio)
353  ratio = 0.999
354  CALL pushcontrol(2,2)
355  ELSE
356  CALL pushcontrol(2,1)
357  END IF
358  ELSE
359  CALL pushrealarray(ratio)
360  ratio = REAL(n)/REAL(m)
361  CALL pushcontrol(2,0)
362  END IF
363  DO i=is,ie
364  CALL pushrealarray(gzh(i))
365  gzh(i) = 0.
366  END DO
367 ! Compute total condensate
368  IF (nwat .LT. 2) THEN
369  DO k=1,kbot
370  DO i=is,ie
371  CALL pushrealarray(qcon(i, k))
372  qcon(i, k) = 0.
373  END DO
374  END DO
375  CALL pushcontrol(3,4)
376  ELSE IF (nwat .EQ. 2) THEN
377 ! GFS_2015
378  DO k=1,kbot
379  DO i=is,ie
380  CALL pushrealarray(qcon(i, k))
381  qcon(i, k) = q0(i, k, liq_wat)
382  END DO
383  END DO
384  CALL pushcontrol(3,3)
385  ELSE IF (nwat .EQ. 3) THEN
386  DO k=1,kbot
387  DO i=is,ie
388  CALL pushrealarray(qcon(i, k))
389  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, ice_wat)
390  END DO
391  END DO
392  CALL pushcontrol(3,2)
393  ELSE IF (nwat .EQ. 4) THEN
394  DO k=1,kbot
395  DO i=is,ie
396  CALL pushrealarray(qcon(i, k))
397  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, rainwat)
398  END DO
399  END DO
400  CALL pushcontrol(3,1)
401  ELSE
402  DO k=1,kbot
403  DO i=is,ie
404  CALL pushrealarray(qcon(i, k))
405  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, ice_wat) + q0(i&
406 & , k, snowwat) + q0(i, k, rainwat) + q0(i, k, graupel)
407  END DO
408  END DO
409  CALL pushcontrol(3,0)
410  END IF
411  DO k=kbot,2,-1
412  km1 = k - 1
413  DO i=is,ie
414 ! Richardson number = g*delz * del_theta/theta / (del_u**2 + del_v**2)
415 ! Use exact form for "density temperature"
416  tv1 = t0(i, km1)*(1.+xvir*q0(i, km1, sphum)-qcon(i, km1))
417  tv2 = t0(i, k)*(1.+xvir*q0(i, k, sphum)-qcon(i, k))
418  pt1 = tv1/pkz(i, j, km1)
419  pt2 = tv2/pkz(i, j, k)
420 !
421  CALL pushrealarray(ri)
422  ri = (gz(i, km1)-gz(i, k))*(pt1-pt2)/(0.5*(pt1+pt2)*((u0(i, &
423 & km1)-u0(i, k))**2+(v0(i, km1)-v0(i, k))**2+ustar2))
424  IF (tv1 .GT. t_max .AND. tv1 .GT. tv2) THEN
425 ! top layer unphysically warm
426  ri = 0.
427  CALL pushcontrol(2,0)
428  ELSE IF (tv2 .LT. t_min) THEN
429  IF (ri .GT. 0.2) THEN
430  ri = 0.2
431  CALL pushcontrol(2,3)
432  ELSE
433  CALL pushcontrol(2,2)
434  ri = ri
435  END IF
436  ELSE
437  CALL pushcontrol(2,1)
438  END IF
439  IF (400.e2 - pm(i, k) .LT. 0.) THEN
440  CALL pushcontrol(1,0)
441  max2 = 0.
442  ELSE
443  max2 = 400.e2 - pm(i, k)
444  CALL pushcontrol(1,1)
445  END IF
446  y1 = ri_min + (ri_max-ri_min)*max2/200.e2
447  IF (ri_max .GT. y1) THEN
448  CALL pushrealarray(ri_ref)
449  ri_ref = y1
450  CALL pushcontrol(1,0)
451  ELSE
452  CALL pushrealarray(ri_ref)
453  ri_ref = ri_max
454  CALL pushcontrol(1,1)
455  END IF
456 ! Enhancing mixing at the model top
457  IF (k .EQ. 2) THEN
458  ri_ref = 4.*ri_ref
459  CALL pushcontrol(2,0)
460  ELSE IF (k .EQ. 3) THEN
461  ri_ref = 2.*ri_ref
462  CALL pushcontrol(2,1)
463  ELSE IF (k .EQ. 4) THEN
464  ri_ref = 1.5*ri_ref
465  CALL pushcontrol(2,2)
466  ELSE
467  CALL pushcontrol(2,3)
468  END IF
469  IF (ri .LT. ri_ref) THEN
470  IF (0.0 .LT. ri/ri_ref) THEN
471  CALL pushrealarray(max1)
472  max1 = ri/ri_ref
473  CALL pushcontrol(1,0)
474  ELSE
475  CALL pushrealarray(max1)
476  max1 = 0.0
477  CALL pushcontrol(1,1)
478  END IF
479  CALL pushrealarray(mc)
480  mc = ratio*delp(i, j, km1)*delp(i, j, k)/(delp(i, j, km1)+&
481 & delp(i, j, k))*(1.-max1)**2
482  DO iq=1,nq
483  CALL pushrealarray(h0)
484  h0 = mc*(q0(i, k, iq)-q0(i, km1, iq))
485  CALL pushrealarray(q0(i, km1, iq))
486  q0(i, km1, iq) = q0(i, km1, iq) + h0/delp(i, j, km1)
487  CALL pushrealarray(q0(i, k, iq))
488  q0(i, k, iq) = q0(i, k, iq) - h0/delp(i, j, k)
489  END DO
490 ! Recompute qcon
491  IF (nwat .LT. 2) THEN
492  CALL pushrealarray(qcon(i, km1))
493  qcon(i, km1) = 0.
494  CALL pushcontrol(3,0)
495  ELSE IF (nwat .EQ. 2) THEN
496 ! GFS_2015
497  CALL pushrealarray(qcon(i, km1))
498  qcon(i, km1) = q0(i, km1, liq_wat)
499  CALL pushcontrol(3,1)
500  ELSE IF (nwat .EQ. 3) THEN
501 ! AM3/AM4
502  CALL pushrealarray(qcon(i, km1))
503  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, ice_wat)
504  CALL pushcontrol(3,2)
505  ELSE IF (nwat .EQ. 4) THEN
506 ! K_warm_rain scheme with fake ice
507  CALL pushrealarray(qcon(i, km1))
508  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, rainwat)
509  CALL pushcontrol(3,3)
510  ELSE
511  CALL pushrealarray(qcon(i, km1))
512  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, ice_wat)&
513 & + q0(i, km1, snowwat) + q0(i, km1, rainwat) + q0(i, &
514 & km1, graupel)
515  CALL pushcontrol(3,4)
516  END IF
517 ! u:
518  CALL pushrealarray(h0)
519  h0 = mc*(u0(i, k)-u0(i, k-1))
520  CALL pushrealarray(u0(i, k-1))
521  u0(i, k-1) = u0(i, k-1) + h0/delp(i, j, k-1)
522  CALL pushrealarray(u0(i, k))
523  u0(i, k) = u0(i, k) - h0/delp(i, j, k)
524 ! v:
525  CALL pushrealarray(h0)
526  h0 = mc*(v0(i, k)-v0(i, k-1))
527  CALL pushrealarray(v0(i, k-1))
528  v0(i, k-1) = v0(i, k-1) + h0/delp(i, j, k-1)
529  CALL pushrealarray(v0(i, k))
530  v0(i, k) = v0(i, k) - h0/delp(i, j, k)
531  IF (hydrostatic) THEN
532 ! Static energy
533  CALL pushrealarray(h0)
534  h0 = mc*(hd(i, k)-hd(i, k-1))
535  CALL pushrealarray(hd(i, k-1))
536  hd(i, k-1) = hd(i, k-1) + h0/delp(i, j, k-1)
537  CALL pushrealarray(hd(i, k))
538  hd(i, k) = hd(i, k) - h0/delp(i, j, k)
539  CALL pushcontrol(2,2)
540  ELSE
541 ! Total energy
542  CALL pushrealarray(h0)
543  h0 = mc*(hd(i, k)-hd(i, k-1))
544  CALL pushrealarray(te(i, k-1))
545  te(i, k-1) = te(i, k-1) + h0/delp(i, j, k-1)
546  CALL pushrealarray(te(i, k))
547  te(i, k) = te(i, k) - h0/delp(i, j, k)
548 ! w:
549  h0 = mc*(w0(i, k)-w0(i, k-1))
550  CALL pushrealarray(w0(i, k-1))
551  w0(i, k-1) = w0(i, k-1) + h0/delp(i, j, k-1)
552  CALL pushrealarray(w0(i, k))
553  w0(i, k) = w0(i, k) - h0/delp(i, j, k)
554  CALL pushcontrol(2,1)
555  END IF
556  ELSE
557  CALL pushcontrol(2,0)
558  END IF
559  END DO
560 !--------------
561 ! Retrive Temp:
562 !--------------
563  IF (hydrostatic) THEN
564  CALL pushinteger(kk)
565  kk = k
566  DO i=is,ie
567  CALL pushrealarray(t0(i, kk))
568  t0(i, kk) = (hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, kk)&
569 & **2))/(rk-pe(i, kk, j)/pm(i, kk))
570  CALL pushrealarray(gzh(i))
571  gzh(i) = gzh(i) + t0(i, kk)*(peln(i, kk+1, j)-peln(i, kk, &
572 & j))
573  CALL pushrealarray(t0(i, kk))
574  t0(i, kk) = t0(i, kk)/(rdgas+rz*q0(i, kk, sphum))
575  END DO
576  kk = k - 1
577  DO i=is,ie
578  CALL pushrealarray(t0(i, kk))
579  t0(i, kk) = (hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, kk)&
580 & **2))/((rk-pe(i, kk, j)/pm(i, kk))*(rdgas+rz*q0(i, kk, &
581 & sphum)))
582  END DO
583  CALL pushcontrol(1,1)
584  ELSE
585  ad_from = k - 1
586  CALL pushinteger(kk)
587 ! Non-hydrostatic under constant volume heating/cooling
588  DO kk=ad_from,k
589  IF (nwat .EQ. 0) THEN
590  DO i=is,ie
591  CALL pushrealarray(cpm(i))
592  cpm(i) = cp_air
593  CALL pushrealarray(cvm(i))
594  cvm(i) = cv_air
595  END DO
596  CALL pushcontrol(3,5)
597  ELSE IF (nwat .EQ. 1) THEN
598  DO i=is,ie
599  CALL pushrealarray(cpm(i))
600  cpm(i) = (1.-q0(i, kk, sphum))*cp_air + q0(i, kk, &
601 & sphum)*cp_vapor
602  CALL pushrealarray(cvm(i))
603  cvm(i) = (1.-q0(i, kk, sphum))*cv_air + q0(i, kk, &
604 & sphum)*cv_vap
605  END DO
606  CALL pushcontrol(3,4)
607  ELSE IF (nwat .EQ. 2) THEN
608  DO i=is,ie
609  CALL pushrealarray(cpm(i))
610  cpm(i) = (1.-q0(i, kk, sphum))*cp_air + q0(i, kk, &
611 & sphum)*cp_vapor
612  CALL pushrealarray(cvm(i))
613  cvm(i) = (1.-q0(i, kk, sphum))*cv_air + q0(i, kk, &
614 & sphum)*cv_vap
615  END DO
616  CALL pushcontrol(3,3)
617  ELSE IF (nwat .EQ. 3) THEN
618  DO i=is,ie
619  q_liq = q0(i, kk, liq_wat)
620  q_sol = q0(i, kk, ice_wat)
621  CALL pushrealarray(cpm(i))
622  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cp_air + &
623 & q0(i, kk, sphum)*cp_vapor + q_liq*c_liq + q_sol*&
624 & c_ice
625  CALL pushrealarray(cvm(i))
626  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cv_air + &
627 & q0(i, kk, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
628  END DO
629  CALL pushcontrol(3,2)
630  ELSE IF (nwat .EQ. 4) THEN
631  DO i=is,ie
632  q_liq = q0(i, kk, liq_wat) + q0(i, kk, rainwat)
633  CALL pushrealarray(cpm(i))
634  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq))*cp_air + q0(i, &
635 & kk, sphum)*cp_vapor + q_liq*c_liq
636  CALL pushrealarray(cvm(i))
637  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq))*cv_air + q0(i, &
638 & kk, sphum)*cv_vap + q_liq*c_liq
639  END DO
640  CALL pushcontrol(3,1)
641  ELSE
642  DO i=is,ie
643  q_liq = q0(i, kk, liq_wat) + q0(i, kk, rainwat)
644  q_sol = q0(i, kk, ice_wat) + q0(i, kk, snowwat) + q0(i&
645 & , kk, graupel)
646  CALL pushrealarray(cpm(i))
647  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cp_air + &
648 & q0(i, kk, sphum)*cp_vapor + q_liq*c_liq + q_sol*&
649 & c_ice
650  CALL pushrealarray(cvm(i))
651  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cv_air + &
652 & q0(i, kk, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
653  END DO
654  CALL pushcontrol(3,0)
655  END IF
656  DO i=is,ie
657  CALL pushrealarray(tv)
658  tv = gz(i, kk) + 0.5*(u0(i, kk)**2+v0(i, kk)**2+w0(i, kk&
659 & )**2)
660  CALL pushrealarray(t0(i, kk))
661  t0(i, kk) = (te(i, kk)-tv)/cvm(i)
662  CALL pushrealarray(hd(i, kk))
663  hd(i, kk) = cpm(i)*t0(i, kk) + tv
664  END DO
665  END DO
666  CALL pushinteger(kk - 1)
667  CALL pushinteger(ad_from)
668  CALL pushcontrol(1,0)
669  END IF
670  END DO
671  END DO
672 ! k-loop
673 ! n-loop
674 !--------------------
675  IF (fra .LT. 1.) THEN
676  DO k=1,kbot
677  DO i=is,ie
678  CALL pushrealarray(t0(i, k))
679  t0(i, k) = ta(i, j, k) + (t0(i, k)-ta(i, j, k))*fra
680  CALL pushrealarray(u0(i, k))
681  u0(i, k) = ua(i, j, k) + (u0(i, k)-ua(i, j, k))*fra
682  CALL pushrealarray(v0(i, k))
683  v0(i, k) = va(i, j, k) + (v0(i, k)-va(i, j, k))*fra
684  END DO
685  END DO
686  IF (.NOT.hydrostatic) THEN
687  DO k=1,kbot
688  DO i=is,ie
689  CALL pushrealarray(w0(i, k))
690  w0(i, k) = w(i, j, k) + (w0(i, k)-w(i, j, k))*fra
691  END DO
692  END DO
693  CALL pushcontrol(1,1)
694  ELSE
695  CALL pushcontrol(1,0)
696  END IF
697  DO iq=1,nq
698  DO k=1,kbot
699  DO i=is,ie
700  CALL pushrealarray(q0(i, k, iq))
701  q0(i, k, iq) = qa(i, j, k, iq) + (q0(i, k, iq)-qa(i, j, k&
702 & , iq))*fra
703  END DO
704  END DO
705  END DO
706  CALL pushcontrol(1,1)
707  ELSE
708  CALL pushcontrol(1,0)
709  END IF
710  DO k=1,kbot
711  DO i=is,ie
712 ! *** temperature updated ***
713  CALL pushrealarray(ta(i, j, k))
714  ta(i, j, k) = t0(i, k)
715  END DO
716  END DO
717  IF (.NOT.hydrostatic) THEN
718  CALL pushcontrol(1,0)
719  ELSE
720  CALL pushcontrol(1,1)
721  END IF
722  END DO
723  CALL pushrealarray(te, (ie-is+1)*km)
724  CALL pushinteger(liq_wat)
725  CALL pushinteger(ice_wat)
726  CALL pushrealarray(xvir)
727  CALL pushinteger(sphum)
728  CALL pushrealarray(pm, (ie-is+1)*km)
729  CALL pushrealarray(mc)
730  CALL pushrealarray(u0, (ie-is+1)*km)
731  CALL pushinteger(rainwat)
732  CALL pushrealarray(qcon, (ie-is+1)*km)
733  CALL pushrealarray(h0)
734  CALL pushinteger(kbot)
735  CALL pushrealarray(gzh, ie - is + 1)
736  CALL pushinteger(snowwat)
737  CALL pushrealarray(ratio)
738  CALL pushrealarray(cv_air)
739  CALL pushrealarray(rz)
740  CALL pushrealarray(ri_ref)
741  CALL pushrealarray(fra)
742  CALL pushrealarray(q0, (ie-is+1)*km*nq)
743  CALL pushrealarray(t0, (ie-is+1)*km)
744  CALL pushrealarray(g2)
745  CALL pushrealarray(rk)
746  CALL pushrealarray(ri)
747  CALL pushrealarray(w0, (ie-is+1)*km)
748  CALL pushrealarray(hd, (ie-is+1)*km)
749  CALL pushinteger(kk)
750  CALL pushrealarray(cpm, ie - is + 1)
751  CALL pushrealarray(gz, (ie-is+1)*km)
752  CALL pushrealarray(tv)
753  CALL pushinteger(m)
754  CALL pushrealarray(v0, (ie-is+1)*km)
755  CALL pushrealarray(cvm, ie - is + 1)
756  CALL pushrealarray(max1)
757  CALL pushinteger(graupel)
758  CALL pushcontrol(1,1)
759  100 CONTINUE
760  END SUBROUTINE fv_subgrid_z_fwd
761 ! Differentiation of fv_subgrid_z in reverse (adjoint) mode, backward sweep (with options split(a2b_edge_mod.a2b_ord4 a2b_edg
762 !e_mod.a2b_ord2 dyn_core_mod.dyn_core dyn_core_mod.pk3_halo dyn_core_mod.pln_halo dyn_core_mod.pe_halo dyn_core_mod.adv_pe dyn_co
763 !re_mod.p_grad_c dyn_core_mod.nh_p_grad dyn_core_mod.split_p_grad dyn_core_mod.one_grad_p dyn_core_mod.grad1_p_update dyn_core_mo
764 !d.mix_dp dyn_core_mod.geopk dyn_core_mod.del2_cubed dyn_core_mod.Rayleigh_fast fv_dynamics_mod.fv_dynamics fv_dynamics_mod.Rayle
765 !igh_Super fv_dynamics_mod.Rayleigh_Friction fv_dynamics_mod.compute_aam fv_grid_utils_mod.cubed_to_latlon fv_grid_utils_mod.c2l_
766 !ord4 fv_grid_utils_mod.c2l_ord2 fv_mapz_mod.Lagrangian_to_Eulerian fv_mapz_mod.compute_total_energy fv_mapz_mod.pkez fv_mapz_mod
767 !.remap_z fv_mapz_mod.map_scalar fv_mapz_mod.map1_ppm fv_mapz_mod.mapn_tracer fv_mapz_mod.map1_q2 fv_mapz_mod.remap_2
768 !d fv_mapz_mod.scalar_profile fv_mapz_mod.cs_profile fv_mapz_mod.cs_limiters fv_mapz_mod.ppm_profile fv_mapz_mod.ppm_limite
769 !rs fv_mapz_mod.steepz fv_mapz_mod.rst_remap fv_mapz_mod.mappm fv_mapz_mod.moist_cv fv_mapz_mod.moist_cp fv_mapz_mod.map1_cubic f
770 !v_restart_mod.d2c_setup fv_tracer2d_mod.tracer_2d_1L fv_tracer2d_mod.tracer_2d fv_tracer2d_mod.tracer_2d_nested fv_sg_mod.fv_sub
771 !grid_z main_mod.compute_pressures main_mod.run nh_core_mod.Riem_Solver3 nh_utils_mod.update_dz_c nh_utils_mod.update_dz_d nh_uti
772 !ls_mod.Riem_Solver_c nh_utils_mod.Riem_Solver3test nh_utils_mod.imp_diff_w nh_utils_mod.RIM_2D nh_utils_mod.SIM3_solver nh_utils
773 !_mod.SIM3p0_solver nh_utils_mod.SIM1_solver nh_utils_mod.SIM_solver nh_utils_mod.edge_scalar nh_utils_mod.edge_profile nh_utils_
774 !mod.nest_halo_nh sw_core_mod.c_sw sw_core_mod.d_sw sw_core_mod.divergence_corner sw_core_mod.divergence_corner_nest sw_core_mod
775 !.d2a2c_vect sw_core_mod.fill3_4corners sw_core_mod.fill2_4corners sw_core_mod.fill_4corners sw_core_mod.xtp_u sw_core_mod.ytp
776 !_v_fb sw_core_mod.compute_divergence_damping sw_core_mod.smag_corner tp_core_mod.mp_ghost_ew tp_core_mod.fv_tp_2d tp_co
777 !re_mod.copy_corners tp_core_mod.xppm tp_core_mod.yppm tp_core_mod.deln_flux a2b_edge_mod.extrap_corner fv_grid_ut
778 !ils_mod.great_circle_dist sw_core_mod.edge_interpolate4)):
779 ! gradient of useful results: qa w delp delz pkz ta
780 ! with respect to varying inputs: qa peln w delp ua delz va pkz
781 ! pe ta
782  SUBROUTINE fv_subgrid_z_bwd(isd, ied, jsd, jed, is, ie, js, je, km, nq&
783 & , dt, tau, nwat, delp, delp_ad, pe, pe_ad, peln, peln_ad, pkz, &
784 & pkz_ad, ta, ta_ad, qa, qa_ad, ua, ua_ad, va, va_ad, hydrostatic, w, &
785 & w_ad, delz, delz_ad, u_dt, v_dt, t_dt, k_bot)
786  IMPLICIT NONE
787  INTEGER, INTENT(IN) :: is, ie, js, je, km, nq, nwat
788  INTEGER, INTENT(IN) :: isd, ied, jsd, jed
789  INTEGER, INTENT(IN) :: tau
790  REAL, INTENT(IN) :: dt
791  REAL, INTENT(IN) :: pe(is-1:ie+1, km+1, js-1:je+1)
792  REAL :: pe_ad(is-1:ie+1, km+1, js-1:je+1)
793  REAL, INTENT(IN) :: peln(is:ie, km+1, js:je)
794  REAL :: peln_ad(is:ie, km+1, js:je)
795  REAL, INTENT(IN) :: delp(isd:ied, jsd:jed, km)
796  REAL :: delp_ad(isd:ied, jsd:jed, km)
797  REAL, INTENT(IN) :: delz(isd:, jsd:, :)
798  REAL :: delz_ad(isd:, jsd:, :)
799  REAL, INTENT(IN) :: pkz(is:ie, js:je, km)
800  REAL :: pkz_ad(is:ie, js:je, km)
801  LOGICAL, INTENT(IN) :: hydrostatic
802  INTEGER, INTENT(IN), OPTIONAL :: k_bot
803  REAL, INTENT(INOUT) :: ua(isd:ied, jsd:jed, km)
804  REAL, INTENT(INOUT) :: ua_ad(isd:ied, jsd:jed, km)
805  REAL, INTENT(INOUT) :: va(isd:ied, jsd:jed, km)
806  REAL, INTENT(INOUT) :: va_ad(isd:ied, jsd:jed, km)
807  REAL, INTENT(INOUT) :: w(isd:, jsd:, :)
808  REAL, INTENT(INOUT) :: w_ad(isd:, jsd:, :)
809  REAL, INTENT(INOUT) :: ta(isd:ied, jsd:jed, km)
810  REAL, INTENT(INOUT) :: ta_ad(isd:ied, jsd:jed, km)
811  REAL, INTENT(INOUT) :: qa(isd:ied, jsd:jed, km, nq)
812  REAL, INTENT(INOUT) :: qa_ad(isd:ied, jsd:jed, km, nq)
813  REAL, INTENT(INOUT) :: u_dt(isd:ied, jsd:jed, km)
814  REAL, INTENT(INOUT) :: v_dt(isd:ied, jsd:jed, km)
815  REAL, INTENT(INOUT) :: t_dt(is:ie, js:je, km)
816  REAL, DIMENSION(is:ie, km) :: u0, v0, w0, t0, hd, te, gz, tvm, pm, &
817 & den
818  REAL, DIMENSION(is:ie, km) :: u0_ad, v0_ad, w0_ad, t0_ad, hd_ad, &
819 & te_ad, gz_ad, tvm_ad, pm_ad
820  REAL :: q0(is:ie, km, nq), qcon(is:ie, km)
821  REAL :: q0_ad(is:ie, km, nq), qcon_ad(is:ie, km)
822  REAL, DIMENSION(is:ie) :: gzh, lcp2, icp2, cvm, cpm, qs
823  REAL, DIMENSION(is:ie) :: gzh_ad, cvm_ad, cpm_ad
824  REAL :: ri_ref, ri, pt1, pt2, ratio, tv, cv, tmp, q_liq, q_sol
825  REAL :: ri_ref_ad, ri_ad, pt1_ad, pt2_ad, tv_ad, tmp_ad, q_liq_ad, &
826 & q_sol_ad
827  REAL :: tv1, tv2, g2, h0, mc, fra, rk, rz, rdt, tvd, tv_surf
828  REAL :: tv1_ad, tv2_ad, h0_ad, mc_ad
829  REAL :: dh, dq, qsw, dqsdt, tcp3, t_max, t_min
830  INTEGER :: i, j, k, kk, n, m, iq, km1, im, kbot
831  REAL, PARAMETER :: ustar2=1.e-4
832  REAL :: cv_air, xvir
833  INTEGER :: sphum, liq_wat, rainwat, snowwat, graupel, ice_wat, &
834 & cld_amt
835  INTRINSIC present
836  INTRINSIC min
837  INTRINSIC real
838  INTRINSIC max
839  INTEGER :: min1
840  REAL :: max1
841  REAL :: max1_ad
842  REAL :: max2
843  REAL :: max2_ad
844  REAL :: temp
845  REAL :: temp_ad
846  REAL :: temp0
847  REAL :: temp1
848  REAL :: temp_ad0
849  REAL :: temp_ad1
850  REAL :: temp_ad2
851  REAL :: temp_ad3
852  REAL :: temp_ad4
853  REAL :: temp_ad5
854  REAL :: temp2
855  REAL :: temp3
856  REAL :: temp4
857  REAL :: temp5
858  REAL :: temp6
859  REAL :: temp_ad6
860  REAL :: temp_ad7
861  REAL :: temp_ad8
862  REAL :: temp_ad9
863  REAL :: temp_ad10
864  REAL :: temp_ad11
865  REAL :: temp_ad12
866  REAL :: temp_ad13
867  REAL :: temp_ad14
868  REAL :: temp_ad15
869  REAL :: y1_ad
870  REAL :: temp7
871  REAL :: temp8
872  REAL :: temp_ad16
873  REAL :: temp_ad17
874  REAL :: temp_ad18
875  REAL :: temp_ad19
876  REAL :: temp_ad20
877  REAL :: temp_ad21
878  REAL :: temp_ad22
879  REAL :: temp_ad23
880  REAL :: temp_ad24
881  REAL :: temp_ad25
882  REAL :: temp_ad26
883  REAL :: temp_ad27
884  REAL :: temp_ad28
885  REAL :: temp_ad29
886  REAL :: temp_ad30
887  REAL :: temp9
888  REAL :: temp10
889  REAL :: temp11
890  REAL :: temp_ad31
891  REAL :: temp_ad32
892  REAL :: temp_ad33
893  REAL :: temp12
894  REAL :: temp13
895  REAL :: temp14
896  REAL :: temp15
897  REAL :: temp_ad34
898  REAL :: temp_ad35
899  REAL :: temp_ad36
900  REAL :: temp_ad37
901  REAL :: temp_ad38
902  REAL :: temp_ad39
903  REAL :: temp_ad40
904  REAL :: temp_ad41
905  REAL :: temp_ad42
906  INTEGER :: branch
907  INTEGER :: ad_from
908  INTEGER :: ad_to
909  REAL :: y1
910  CALL popcontrol(1,branch)
911  IF (branch .EQ. 0) THEN
912  peln_ad = 0.0
913  ua_ad = 0.0
914  va_ad = 0.0
915  pe_ad = 0.0
916  ELSE
917  CALL popinteger(graupel)
918  CALL poprealarray(max1)
919  CALL poprealarray(cvm, ie - is + 1)
920  CALL poprealarray(v0, (ie-is+1)*km)
921  CALL popinteger(m)
922  CALL poprealarray(tv)
923  CALL poprealarray(gz, (ie-is+1)*km)
924  CALL poprealarray(cpm, ie - is + 1)
925  CALL popinteger(kk)
926  CALL poprealarray(hd, (ie-is+1)*km)
927  CALL poprealarray(w0, (ie-is+1)*km)
928  CALL poprealarray(ri)
929  CALL poprealarray(rk)
930  CALL poprealarray(g2)
931  CALL poprealarray(t0, (ie-is+1)*km)
932  CALL poprealarray(q0, (ie-is+1)*km*nq)
933  CALL poprealarray(fra)
934  CALL poprealarray(ri_ref)
935  CALL poprealarray(rz)
936  CALL poprealarray(cv_air)
937  CALL poprealarray(ratio)
938  CALL popinteger(snowwat)
939  CALL poprealarray(gzh, ie - is + 1)
940  CALL popinteger(kbot)
941  CALL poprealarray(h0)
942  CALL poprealarray(qcon, (ie-is+1)*km)
943  CALL popinteger(rainwat)
944  CALL poprealarray(u0, (ie-is+1)*km)
945  CALL poprealarray(mc)
946  CALL poprealarray(pm, (ie-is+1)*km)
947  CALL popinteger(sphum)
948  CALL poprealarray(xvir)
949  CALL popinteger(ice_wat)
950  CALL popinteger(liq_wat)
951  CALL poprealarray(te, (ie-is+1)*km)
952  peln_ad = 0.0
953  ua_ad = 0.0
954  va_ad = 0.0
955  pe_ad = 0.0
956  cvm_ad = 0.0
957  v0_ad = 0.0
958  gz_ad = 0.0
959  cpm_ad = 0.0
960  hd_ad = 0.0
961  w0_ad = 0.0
962  tvm_ad = 0.0
963  t0_ad = 0.0
964  q0_ad = 0.0
965  gzh_ad = 0.0
966  qcon_ad = 0.0
967  u0_ad = 0.0
968  pm_ad = 0.0
969  te_ad = 0.0
970  DO j=je,js,-1
971  CALL popcontrol(1,branch)
972  IF (branch .EQ. 0) THEN
973  DO k=kbot,1,-1
974  DO i=ie,is,-1
975  w0_ad(i, k) = w0_ad(i, k) + w_ad(i, j, k)
976  w_ad(i, j, k) = 0.0
977  END DO
978  END DO
979  END IF
980  DO k=kbot,1,-1
981  DO iq=nq,1,-1
982  DO i=ie,is,-1
983  q0_ad(i, k, iq) = q0_ad(i, k, iq) + qa_ad(i, j, k, iq)
984  qa_ad(i, j, k, iq) = 0.0
985  END DO
986  END DO
987  DO i=ie,is,-1
988  v0_ad(i, k) = v0_ad(i, k) + va_ad(i, j, k)
989  va_ad(i, j, k) = 0.0
990  u0_ad(i, k) = u0_ad(i, k) + ua_ad(i, j, k)
991  ua_ad(i, j, k) = 0.0
992  CALL poprealarray(ta(i, j, k))
993  t0_ad(i, k) = t0_ad(i, k) + ta_ad(i, j, k)
994  ta_ad(i, j, k) = 0.0
995  END DO
996  END DO
997  CALL popcontrol(1,branch)
998  IF (branch .NE. 0) THEN
999  DO iq=nq,1,-1
1000  DO k=kbot,1,-1
1001  DO i=ie,is,-1
1002  CALL poprealarray(q0(i, k, iq))
1003  qa_ad(i, j, k, iq) = qa_ad(i, j, k, iq) + (1.0-fra)*&
1004 & q0_ad(i, k, iq)
1005  q0_ad(i, k, iq) = fra*q0_ad(i, k, iq)
1006  END DO
1007  END DO
1008  END DO
1009  CALL popcontrol(1,branch)
1010  IF (branch .NE. 0) THEN
1011  DO k=kbot,1,-1
1012  DO i=ie,is,-1
1013  CALL poprealarray(w0(i, k))
1014  w_ad(i, j, k) = w_ad(i, j, k) + (1.0-fra)*w0_ad(i, k)
1015  w0_ad(i, k) = fra*w0_ad(i, k)
1016  END DO
1017  END DO
1018  END IF
1019  DO k=kbot,1,-1
1020  DO i=ie,is,-1
1021  CALL poprealarray(v0(i, k))
1022  va_ad(i, j, k) = va_ad(i, j, k) + (1.0-fra)*v0_ad(i, k)
1023  v0_ad(i, k) = fra*v0_ad(i, k)
1024  CALL poprealarray(u0(i, k))
1025  ua_ad(i, j, k) = ua_ad(i, j, k) + (1.0-fra)*u0_ad(i, k)
1026  u0_ad(i, k) = fra*u0_ad(i, k)
1027  CALL poprealarray(t0(i, k))
1028  ta_ad(i, j, k) = ta_ad(i, j, k) + (1.0-fra)*t0_ad(i, k)
1029  t0_ad(i, k) = fra*t0_ad(i, k)
1030  END DO
1031  END DO
1032  END IF
1033  DO n=m,1,-1
1034  DO k=2,kbot,1
1035  CALL popcontrol(1,branch)
1036  IF (branch .EQ. 0) THEN
1037  CALL popinteger(ad_from)
1038  CALL popinteger(ad_to)
1039  DO kk=ad_to,ad_from,-1
1040  DO i=ie,is,-1
1041  CALL poprealarray(hd(i, kk))
1042  cpm_ad(i) = cpm_ad(i) + t0(i, kk)*hd_ad(i, kk)
1043  t0_ad(i, kk) = t0_ad(i, kk) + cpm(i)*hd_ad(i, kk)
1044  CALL poprealarray(t0(i, kk))
1045  temp_ad41 = t0_ad(i, kk)/cvm(i)
1046  tv_ad = hd_ad(i, kk) - temp_ad41
1047  hd_ad(i, kk) = 0.0
1048  te_ad(i, kk) = te_ad(i, kk) + temp_ad41
1049  cvm_ad(i) = cvm_ad(i) - (te(i, kk)-tv)*temp_ad41/cvm(i&
1050 & )
1051  t0_ad(i, kk) = 0.0
1052  CALL poprealarray(tv)
1053  temp_ad42 = 0.5*tv_ad
1054  gz_ad(i, kk) = gz_ad(i, kk) + tv_ad
1055  u0_ad(i, kk) = u0_ad(i, kk) + 2*u0(i, kk)*temp_ad42
1056  v0_ad(i, kk) = v0_ad(i, kk) + 2*v0(i, kk)*temp_ad42
1057  w0_ad(i, kk) = w0_ad(i, kk) + 2*w0(i, kk)*temp_ad42
1058  END DO
1059  CALL popcontrol(3,branch)
1060  IF (branch .LT. 3) THEN
1061  IF (branch .EQ. 0) THEN
1062  DO i=ie,is,-1
1063  temp_ad40 = cp_air*cpm_ad(i)
1064  CALL poprealarray(cvm(i))
1065  temp_ad39 = cv_air*cvm_ad(i)
1066  q_liq_ad = c_liq*cpm_ad(i) - temp_ad40 + c_liq*&
1067 & cvm_ad(i) - temp_ad39
1068  q_sol_ad = c_ice*cpm_ad(i) - temp_ad40 + c_ice*&
1069 & cvm_ad(i) - temp_ad39
1070  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + &
1071 & cp_vapor*cpm_ad(i) - temp_ad40 + cv_vap*cvm_ad(i&
1072 & ) - temp_ad39
1073  cvm_ad(i) = 0.0
1074  CALL poprealarray(cpm(i))
1075  cpm_ad(i) = 0.0
1076  q0_ad(i, kk, ice_wat) = q0_ad(i, kk, ice_wat) + &
1077 & q_sol_ad
1078  q0_ad(i, kk, snowwat) = q0_ad(i, kk, snowwat) + &
1079 & q_sol_ad
1080  q0_ad(i, kk, graupel) = q0_ad(i, kk, graupel) + &
1081 & q_sol_ad
1082  q0_ad(i, kk, liq_wat) = q0_ad(i, kk, liq_wat) + &
1083 & q_liq_ad
1084  q0_ad(i, kk, rainwat) = q0_ad(i, kk, rainwat) + &
1085 & q_liq_ad
1086  END DO
1087  ELSE IF (branch .EQ. 1) THEN
1088  DO i=ie,is,-1
1089  CALL poprealarray(cvm(i))
1090  q_liq_ad = (c_liq-cp_air)*cpm_ad(i) + (c_liq-&
1091 & cv_air)*cvm_ad(i)
1092  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + (&
1093 & cp_vapor-cp_air)*cpm_ad(i) + (cv_vap-cv_air)*&
1094 & cvm_ad(i)
1095  cvm_ad(i) = 0.0
1096  CALL poprealarray(cpm(i))
1097  cpm_ad(i) = 0.0
1098  q0_ad(i, kk, liq_wat) = q0_ad(i, kk, liq_wat) + &
1099 & q_liq_ad
1100  q0_ad(i, kk, rainwat) = q0_ad(i, kk, rainwat) + &
1101 & q_liq_ad
1102  END DO
1103  ELSE
1104  DO i=ie,is,-1
1105  temp_ad38 = cp_air*cpm_ad(i)
1106  CALL poprealarray(cvm(i))
1107  temp_ad37 = cv_air*cvm_ad(i)
1108  q_liq_ad = c_liq*cpm_ad(i) - temp_ad38 + c_liq*&
1109 & cvm_ad(i) - temp_ad37
1110  q_sol_ad = c_ice*cpm_ad(i) - temp_ad38 + c_ice*&
1111 & cvm_ad(i) - temp_ad37
1112  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + &
1113 & cp_vapor*cpm_ad(i) - temp_ad38 + cv_vap*cvm_ad(i&
1114 & ) - temp_ad37
1115  cvm_ad(i) = 0.0
1116  CALL poprealarray(cpm(i))
1117  cpm_ad(i) = 0.0
1118  q0_ad(i, kk, ice_wat) = q0_ad(i, kk, ice_wat) + &
1119 & q_sol_ad
1120  q0_ad(i, kk, liq_wat) = q0_ad(i, kk, liq_wat) + &
1121 & q_liq_ad
1122  END DO
1123  END IF
1124  ELSE IF (branch .EQ. 3) THEN
1125  DO i=ie,is,-1
1126  CALL poprealarray(cvm(i))
1127  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + (&
1128 & cp_vapor-cp_air)*cpm_ad(i) + (cv_vap-cv_air)*&
1129 & cvm_ad(i)
1130  cvm_ad(i) = 0.0
1131  CALL poprealarray(cpm(i))
1132  cpm_ad(i) = 0.0
1133  END DO
1134  ELSE IF (branch .EQ. 4) THEN
1135  DO i=ie,is,-1
1136  CALL poprealarray(cvm(i))
1137  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + (&
1138 & cp_vapor-cp_air)*cpm_ad(i) + (cv_vap-cv_air)*&
1139 & cvm_ad(i)
1140  cvm_ad(i) = 0.0
1141  CALL poprealarray(cpm(i))
1142  cpm_ad(i) = 0.0
1143  END DO
1144  ELSE
1145  DO i=ie,is,-1
1146  CALL poprealarray(cvm(i))
1147  cvm_ad(i) = 0.0
1148  CALL poprealarray(cpm(i))
1149  cpm_ad(i) = 0.0
1150  END DO
1151  END IF
1152  END DO
1153  CALL popinteger(kk)
1154  ELSE
1155  DO i=ie,is,-1
1156  CALL poprealarray(t0(i, kk))
1157  temp15 = rdgas + rz*q0(i, kk, sphum)
1158  temp14 = pm(i, kk)
1159  temp13 = pe(i, kk, j)/temp14
1160  temp12 = (rk-temp13)*temp15
1161  temp_ad34 = t0_ad(i, kk)/temp12
1162  temp_ad35 = -((hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, &
1163 & kk)**2))*temp_ad34/temp12)
1164  temp_ad36 = -(temp15*temp_ad35/temp14)
1165  hd_ad(i, kk) = hd_ad(i, kk) + temp_ad34
1166  gzh_ad(i) = gzh_ad(i) - temp_ad34
1167  u0_ad(i, kk) = u0_ad(i, kk) - 0.5*2*u0(i, kk)*temp_ad34
1168  v0_ad(i, kk) = v0_ad(i, kk) - 0.5*2*v0(i, kk)*temp_ad34
1169  pe_ad(i, kk, j) = pe_ad(i, kk, j) + temp_ad36
1170  pm_ad(i, kk) = pm_ad(i, kk) - temp13*temp_ad36
1171  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) + (rk-temp13)*&
1172 & rz*temp_ad35
1173  t0_ad(i, kk) = 0.0
1174  END DO
1175  kk = k
1176  DO i=ie,is,-1
1177  CALL poprealarray(t0(i, kk))
1178  temp11 = rdgas + rz*q0(i, kk, sphum)
1179  q0_ad(i, kk, sphum) = q0_ad(i, kk, sphum) - t0(i, kk)*rz&
1180 & *t0_ad(i, kk)/temp11**2
1181  t0_ad(i, kk) = (peln(i, kk+1, j)-peln(i, kk, j))*gzh_ad(&
1182 & i) + t0_ad(i, kk)/temp11
1183  CALL poprealarray(gzh(i))
1184  temp_ad31 = t0(i, kk)*gzh_ad(i)
1185  peln_ad(i, kk+1, j) = peln_ad(i, kk+1, j) + temp_ad31
1186  peln_ad(i, kk, j) = peln_ad(i, kk, j) - temp_ad31
1187  CALL poprealarray(t0(i, kk))
1188  temp10 = pm(i, kk)
1189  temp9 = pe(i, kk, j)/temp10
1190  temp_ad32 = t0_ad(i, kk)/(rk-temp9)
1191  temp_ad33 = (hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, kk&
1192 & )**2))*temp_ad32/((rk-temp9)*temp10)
1193  hd_ad(i, kk) = hd_ad(i, kk) + temp_ad32
1194  gzh_ad(i) = gzh_ad(i) - temp_ad32
1195  u0_ad(i, kk) = u0_ad(i, kk) - 0.5*2*u0(i, kk)*temp_ad32
1196  v0_ad(i, kk) = v0_ad(i, kk) - 0.5*2*v0(i, kk)*temp_ad32
1197  pe_ad(i, kk, j) = pe_ad(i, kk, j) + temp_ad33
1198  pm_ad(i, kk) = pm_ad(i, kk) - temp9*temp_ad33
1199  t0_ad(i, kk) = 0.0
1200  END DO
1201  CALL popinteger(kk)
1202  END IF
1203  km1 = k - 1
1204  DO i=ie,is,-1
1205  CALL popcontrol(2,branch)
1206  IF (branch .EQ. 0) THEN
1207  ri_ad = 0.0
1208  ri_ref_ad = 0.0
1209  ELSE
1210  IF (branch .EQ. 1) THEN
1211  temp_ad30 = te_ad(i, k-1)/delp(i, j, k-1)
1212  temp_ad28 = w0_ad(i, k-1)/delp(i, j, k-1)
1213  CALL poprealarray(w0(i, k))
1214  temp_ad27 = -(w0_ad(i, k)/delp(i, j, k))
1215  h0_ad = temp_ad28 + temp_ad27
1216  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad27/&
1217 & delp(i, j, k)
1218  CALL poprealarray(w0(i, k-1))
1219  delp_ad(i, j, k-1) = delp_ad(i, j, k-1) - h0*temp_ad28&
1220 & /delp(i, j, k-1)
1221  mc_ad = (w0(i, k)-w0(i, k-1))*h0_ad
1222  w0_ad(i, k) = w0_ad(i, k) + mc*h0_ad
1223  w0_ad(i, k-1) = w0_ad(i, k-1) - mc*h0_ad
1224  h0 = mc*(hd(i, k)-hd(i, k-1))
1225  CALL poprealarray(te(i, k))
1226  temp_ad29 = -(te_ad(i, k)/delp(i, j, k))
1227  h0_ad = temp_ad30 + temp_ad29
1228  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad29/&
1229 & delp(i, j, k)
1230  CALL poprealarray(te(i, k-1))
1231  delp_ad(i, j, k-1) = delp_ad(i, j, k-1) - h0*temp_ad30&
1232 & /delp(i, j, k-1)
1233  CALL poprealarray(h0)
1234  mc_ad = mc_ad + (hd(i, k)-hd(i, k-1))*h0_ad
1235  hd_ad(i, k) = hd_ad(i, k) + mc*h0_ad
1236  hd_ad(i, k-1) = hd_ad(i, k-1) - mc*h0_ad
1237  ELSE
1238  temp_ad26 = hd_ad(i, k-1)/delp(i, j, k-1)
1239  CALL poprealarray(hd(i, k))
1240  temp_ad25 = -(hd_ad(i, k)/delp(i, j, k))
1241  h0_ad = temp_ad26 + temp_ad25
1242  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad25/&
1243 & delp(i, j, k)
1244  CALL poprealarray(hd(i, k-1))
1245  delp_ad(i, j, k-1) = delp_ad(i, j, k-1) - h0*temp_ad26&
1246 & /delp(i, j, k-1)
1247  CALL poprealarray(h0)
1248  mc_ad = (hd(i, k)-hd(i, k-1))*h0_ad
1249  hd_ad(i, k) = hd_ad(i, k) + mc*h0_ad
1250  hd_ad(i, k-1) = hd_ad(i, k-1) - mc*h0_ad
1251  END IF
1252  temp_ad24 = u0_ad(i, k-1)/delp(i, j, k-1)
1253  temp_ad22 = v0_ad(i, k-1)/delp(i, j, k-1)
1254  CALL poprealarray(v0(i, k))
1255  temp_ad21 = -(v0_ad(i, k)/delp(i, j, k))
1256  h0_ad = temp_ad22 + temp_ad21
1257  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad21/delp(&
1258 & i, j, k)
1259  CALL poprealarray(v0(i, k-1))
1260  delp_ad(i, j, k-1) = delp_ad(i, j, k-1) - h0*temp_ad22/&
1261 & delp(i, j, k-1)
1262  CALL poprealarray(h0)
1263  mc_ad = mc_ad + (v0(i, k)-v0(i, k-1))*h0_ad
1264  v0_ad(i, k) = v0_ad(i, k) + mc*h0_ad
1265  v0_ad(i, k-1) = v0_ad(i, k-1) - mc*h0_ad
1266  CALL poprealarray(u0(i, k))
1267  temp_ad23 = -(u0_ad(i, k)/delp(i, j, k))
1268  h0_ad = temp_ad24 + temp_ad23
1269  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad23/delp(&
1270 & i, j, k)
1271  CALL poprealarray(u0(i, k-1))
1272  delp_ad(i, j, k-1) = delp_ad(i, j, k-1) - h0*temp_ad24/&
1273 & delp(i, j, k-1)
1274  CALL poprealarray(h0)
1275  mc_ad = mc_ad + (u0(i, k)-u0(i, k-1))*h0_ad
1276  u0_ad(i, k) = u0_ad(i, k) + mc*h0_ad
1277  u0_ad(i, k-1) = u0_ad(i, k-1) - mc*h0_ad
1278  CALL popcontrol(3,branch)
1279  IF (branch .LT. 2) THEN
1280  IF (branch .EQ. 0) THEN
1281  CALL poprealarray(qcon(i, km1))
1282  qcon_ad(i, km1) = 0.0
1283  ELSE
1284  CALL poprealarray(qcon(i, km1))
1285  q0_ad(i, km1, liq_wat) = q0_ad(i, km1, liq_wat) + &
1286 & qcon_ad(i, km1)
1287  qcon_ad(i, km1) = 0.0
1288  END IF
1289  ELSE IF (branch .EQ. 2) THEN
1290  CALL poprealarray(qcon(i, km1))
1291  q0_ad(i, km1, liq_wat) = q0_ad(i, km1, liq_wat) + &
1292 & qcon_ad(i, km1)
1293  q0_ad(i, km1, ice_wat) = q0_ad(i, km1, ice_wat) + &
1294 & qcon_ad(i, km1)
1295  qcon_ad(i, km1) = 0.0
1296  ELSE IF (branch .EQ. 3) THEN
1297  CALL poprealarray(qcon(i, km1))
1298  q0_ad(i, km1, liq_wat) = q0_ad(i, km1, liq_wat) + &
1299 & qcon_ad(i, km1)
1300  q0_ad(i, km1, rainwat) = q0_ad(i, km1, rainwat) + &
1301 & qcon_ad(i, km1)
1302  qcon_ad(i, km1) = 0.0
1303  ELSE
1304  CALL poprealarray(qcon(i, km1))
1305  q0_ad(i, km1, liq_wat) = q0_ad(i, km1, liq_wat) + &
1306 & qcon_ad(i, km1)
1307  q0_ad(i, km1, ice_wat) = q0_ad(i, km1, ice_wat) + &
1308 & qcon_ad(i, km1)
1309  q0_ad(i, km1, snowwat) = q0_ad(i, km1, snowwat) + &
1310 & qcon_ad(i, km1)
1311  q0_ad(i, km1, rainwat) = q0_ad(i, km1, rainwat) + &
1312 & qcon_ad(i, km1)
1313  q0_ad(i, km1, graupel) = q0_ad(i, km1, graupel) + &
1314 & qcon_ad(i, km1)
1315  qcon_ad(i, km1) = 0.0
1316  END IF
1317  DO iq=nq,1,-1
1318  temp_ad20 = q0_ad(i, km1, iq)/delp(i, j, km1)
1319  CALL poprealarray(q0(i, k, iq))
1320  temp_ad19 = -(q0_ad(i, k, iq)/delp(i, j, k))
1321  h0_ad = temp_ad20 + temp_ad19
1322  delp_ad(i, j, k) = delp_ad(i, j, k) - h0*temp_ad19/&
1323 & delp(i, j, k)
1324  CALL poprealarray(q0(i, km1, iq))
1325  delp_ad(i, j, km1) = delp_ad(i, j, km1) - h0*temp_ad20&
1326 & /delp(i, j, km1)
1327  CALL poprealarray(h0)
1328  mc_ad = mc_ad + (q0(i, k, iq)-q0(i, km1, iq))*h0_ad
1329  q0_ad(i, k, iq) = q0_ad(i, k, iq) + mc*h0_ad
1330  q0_ad(i, km1, iq) = q0_ad(i, km1, iq) - mc*h0_ad
1331  END DO
1332  CALL poprealarray(mc)
1333  temp8 = delp(i, j, km1) + delp(i, j, k)
1334  temp7 = delp(i, j, k)/temp8
1335  temp_ad16 = ratio*mc_ad
1336  temp_ad17 = delp(i, j, km1)*(1.-max1)**2*temp_ad16/temp8
1337  temp_ad18 = -(temp7*temp_ad17)
1338  delp_ad(i, j, km1) = delp_ad(i, j, km1) + temp_ad18 + &
1339 & temp7*(1.-max1)**2*temp_ad16
1340  max1_ad = -(2*(1.-max1)*delp(i, j, km1)*temp7*temp_ad16)
1341  delp_ad(i, j, k) = delp_ad(i, j, k) + temp_ad18 + &
1342 & temp_ad17
1343  CALL popcontrol(1,branch)
1344  IF (branch .EQ. 0) THEN
1345  CALL poprealarray(max1)
1346  ri_ad = max1_ad/ri_ref
1347  ri_ref_ad = -(ri*max1_ad/ri_ref**2)
1348  ELSE
1349  CALL poprealarray(max1)
1350  ri_ad = 0.0
1351  ri_ref_ad = 0.0
1352  END IF
1353  END IF
1354  CALL popcontrol(2,branch)
1355  IF (branch .LT. 2) THEN
1356  IF (branch .EQ. 0) THEN
1357  ri_ref_ad = 4.*ri_ref_ad
1358  ELSE
1359  ri_ref_ad = 2.*ri_ref_ad
1360  END IF
1361  ELSE IF (branch .EQ. 2) THEN
1362  ri_ref_ad = 1.5*ri_ref_ad
1363  END IF
1364  CALL popcontrol(1,branch)
1365  IF (branch .EQ. 0) THEN
1366  CALL poprealarray(ri_ref)
1367  y1_ad = ri_ref_ad
1368  ELSE
1369  CALL poprealarray(ri_ref)
1370  y1_ad = 0.0
1371  END IF
1372  max2_ad = (ri_max-ri_min)*y1_ad/200.e2
1373  CALL popcontrol(1,branch)
1374  IF (branch .NE. 0) pm_ad(i, k) = pm_ad(i, k) - max2_ad
1375  CALL popcontrol(2,branch)
1376  IF (branch .LT. 2) THEN
1377  IF (branch .EQ. 0) THEN
1378  tv2 = t0(i, k)*(1.+xvir*q0(i, k, sphum)-qcon(i, k))
1379  ri_ad = 0.0
1380  GOTO 100
1381  END IF
1382  ELSE IF (branch .NE. 2) THEN
1383  ri_ad = 0.0
1384  END IF
1385  tv2 = t0(i, k)*(1.+xvir*q0(i, k, sphum)-qcon(i, k))
1386  100 tv1 = t0(i, km1)*(1.+xvir*q0(i, km1, sphum)-qcon(i, km1))
1387  pt1 = tv1/pkz(i, j, km1)
1388  pt2 = tv2/pkz(i, j, k)
1389  CALL poprealarray(ri)
1390  temp5 = v0(i, km1) - v0(i, k)
1391  temp4 = u0(i, km1) - u0(i, k)
1392  temp3 = ustar2 + temp4**2 + temp5**2
1393  temp2 = 0.5*(pt1+pt2)*temp3
1394  temp_ad6 = ri_ad/temp2
1395  temp6 = gz(i, km1) - gz(i, k)
1396  temp_ad7 = -(temp6*(pt1-pt2)*temp_ad6/temp2)
1397  temp_ad8 = temp3*0.5*temp_ad7
1398  temp_ad9 = 0.5*(pt1+pt2)*temp_ad7
1399  temp_ad10 = 2*temp4*temp_ad9
1400  temp_ad11 = 2*temp5*temp_ad9
1401  gz_ad(i, km1) = gz_ad(i, km1) + (pt1-pt2)*temp_ad6
1402  gz_ad(i, k) = gz_ad(i, k) - (pt1-pt2)*temp_ad6
1403  pt1_ad = temp_ad8 + temp6*temp_ad6
1404  pt2_ad = temp_ad8 - temp6*temp_ad6
1405  u0_ad(i, km1) = u0_ad(i, km1) + temp_ad10
1406  u0_ad(i, k) = u0_ad(i, k) - temp_ad10
1407  v0_ad(i, km1) = v0_ad(i, km1) + temp_ad11
1408  v0_ad(i, k) = v0_ad(i, k) - temp_ad11
1409  temp_ad12 = pt2_ad/pkz(i, j, k)
1410  tv2_ad = temp_ad12
1411  pkz_ad(i, j, k) = pkz_ad(i, j, k) - tv2*temp_ad12/pkz(i, j&
1412 & , k)
1413  temp_ad13 = pt1_ad/pkz(i, j, km1)
1414  tv1_ad = temp_ad13
1415  pkz_ad(i, j, km1) = pkz_ad(i, j, km1) - tv1*temp_ad13/pkz(&
1416 & i, j, km1)
1417  temp_ad14 = t0(i, k)*tv2_ad
1418  t0_ad(i, k) = t0_ad(i, k) + (xvir*q0(i, k, sphum)-qcon(i, &
1419 & k)+1.)*tv2_ad
1420  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + xvir*temp_ad14
1421  qcon_ad(i, k) = qcon_ad(i, k) - temp_ad14
1422  temp_ad15 = t0(i, km1)*tv1_ad
1423  t0_ad(i, km1) = t0_ad(i, km1) + (xvir*q0(i, km1, sphum)-&
1424 & qcon(i, km1)+1.)*tv1_ad
1425  q0_ad(i, km1, sphum) = q0_ad(i, km1, sphum) + xvir*&
1426 & temp_ad15
1427  qcon_ad(i, km1) = qcon_ad(i, km1) - temp_ad15
1428  END DO
1429  END DO
1430  CALL popcontrol(3,branch)
1431  IF (branch .LT. 2) THEN
1432  IF (branch .EQ. 0) THEN
1433  DO k=kbot,1,-1
1434  DO i=ie,is,-1
1435  CALL poprealarray(qcon(i, k))
1436  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + qcon_ad(&
1437 & i, k)
1438  q0_ad(i, k, ice_wat) = q0_ad(i, k, ice_wat) + qcon_ad(&
1439 & i, k)
1440  q0_ad(i, k, snowwat) = q0_ad(i, k, snowwat) + qcon_ad(&
1441 & i, k)
1442  q0_ad(i, k, rainwat) = q0_ad(i, k, rainwat) + qcon_ad(&
1443 & i, k)
1444  q0_ad(i, k, graupel) = q0_ad(i, k, graupel) + qcon_ad(&
1445 & i, k)
1446  qcon_ad(i, k) = 0.0
1447  END DO
1448  END DO
1449  ELSE
1450  DO k=kbot,1,-1
1451  DO i=ie,is,-1
1452  CALL poprealarray(qcon(i, k))
1453  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + qcon_ad(&
1454 & i, k)
1455  q0_ad(i, k, rainwat) = q0_ad(i, k, rainwat) + qcon_ad(&
1456 & i, k)
1457  qcon_ad(i, k) = 0.0
1458  END DO
1459  END DO
1460  END IF
1461  ELSE IF (branch .EQ. 2) THEN
1462  DO k=kbot,1,-1
1463  DO i=ie,is,-1
1464  CALL poprealarray(qcon(i, k))
1465  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + qcon_ad(i&
1466 & , k)
1467  q0_ad(i, k, ice_wat) = q0_ad(i, k, ice_wat) + qcon_ad(i&
1468 & , k)
1469  qcon_ad(i, k) = 0.0
1470  END DO
1471  END DO
1472  ELSE IF (branch .EQ. 3) THEN
1473  DO k=kbot,1,-1
1474  DO i=ie,is,-1
1475  CALL poprealarray(qcon(i, k))
1476  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + qcon_ad(i&
1477 & , k)
1478  qcon_ad(i, k) = 0.0
1479  END DO
1480  END DO
1481  ELSE
1482  DO k=kbot,1,-1
1483  DO i=ie,is,-1
1484  CALL poprealarray(qcon(i, k))
1485  qcon_ad(i, k) = 0.0
1486  END DO
1487  END DO
1488  END IF
1489  DO i=ie,is,-1
1490  CALL poprealarray(gzh(i))
1491  gzh_ad(i) = 0.0
1492  END DO
1493  CALL popcontrol(2,branch)
1494  IF (branch .EQ. 0) THEN
1495  CALL poprealarray(ratio)
1496  ELSE
1497  IF (branch .NE. 1) CALL poprealarray(ratio)
1498  CALL popcontrol(1,branch)
1499  IF (branch .EQ. 0) CALL poprealarray(ratio)
1500  CALL popcontrol(1,branch)
1501  IF (branch .EQ. 0) CALL poprealarray(ratio)
1502  END IF
1503  END DO
1504  CALL popcontrol(1,branch)
1505  IF (branch .EQ. 0) THEN
1506  DO k=1,kbot,1
1507  DO i=ie,is,-1
1508  tmp_ad = hd_ad(i, k) + te_ad(i, k)
1509  gz_ad(i, k) = gz_ad(i, k) + tmp_ad
1510  CALL poprealarray(gzh(i))
1511  delz_ad(i, j, k) = delz_ad(i, j, k) - g2*gz_ad(i, k) - &
1512 & grav*gzh_ad(i)
1513  CALL poprealarray(te(i, k))
1514  cvm_ad(i) = cvm_ad(i) + t0(i, k)*te_ad(i, k)
1515  t0_ad(i, k) = t0_ad(i, k) + cpm(i)*hd_ad(i, k) + cvm(i)*&
1516 & te_ad(i, k)
1517  te_ad(i, k) = 0.0
1518  CALL poprealarray(hd(i, k))
1519  cpm_ad(i) = cpm_ad(i) + t0(i, k)*hd_ad(i, k)
1520  hd_ad(i, k) = 0.0
1521  temp_ad5 = 0.5*tmp_ad
1522  u0_ad(i, k) = u0_ad(i, k) + 2*u0(i, k)*temp_ad5
1523  v0_ad(i, k) = v0_ad(i, k) + 2*v0(i, k)*temp_ad5
1524  w0_ad(i, k) = w0_ad(i, k) + 2*w0(i, k)*temp_ad5
1525  CALL poprealarray(gz(i, k))
1526  gzh_ad(i) = gzh_ad(i) + gz_ad(i, k)
1527  gz_ad(i, k) = 0.0
1528  CALL poprealarray(w0(i, k))
1529  w_ad(i, j, k) = w_ad(i, j, k) + w0_ad(i, k)
1530  w0_ad(i, k) = 0.0
1531  END DO
1532  CALL popcontrol(3,branch)
1533  IF (branch .LT. 3) THEN
1534  IF (branch .EQ. 0) THEN
1535  DO i=ie,is,-1
1536  temp_ad4 = cp_air*cpm_ad(i)
1537  CALL poprealarray(cvm(i))
1538  temp_ad3 = cv_air*cvm_ad(i)
1539  q_liq_ad = c_liq*cpm_ad(i) - temp_ad4 + c_liq*cvm_ad(i&
1540 & ) - temp_ad3
1541  q_sol_ad = c_ice*cpm_ad(i) - temp_ad4 + c_ice*cvm_ad(i&
1542 & ) - temp_ad3
1543  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + cp_vapor*&
1544 & cpm_ad(i) - temp_ad4 + cv_vap*cvm_ad(i) - temp_ad3
1545  cvm_ad(i) = 0.0
1546  CALL poprealarray(cpm(i))
1547  cpm_ad(i) = 0.0
1548  q0_ad(i, k, ice_wat) = q0_ad(i, k, ice_wat) + q_sol_ad
1549  q0_ad(i, k, snowwat) = q0_ad(i, k, snowwat) + q_sol_ad
1550  q0_ad(i, k, graupel) = q0_ad(i, k, graupel) + q_sol_ad
1551  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + q_liq_ad
1552  q0_ad(i, k, rainwat) = q0_ad(i, k, rainwat) + q_liq_ad
1553  END DO
1554  ELSE IF (branch .EQ. 1) THEN
1555  DO i=ie,is,-1
1556  CALL poprealarray(cvm(i))
1557  q_liq_ad = (c_liq-cp_air)*cpm_ad(i) + (c_liq-cv_air)*&
1558 & cvm_ad(i)
1559  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + (cp_vapor-&
1560 & cp_air)*cpm_ad(i) + (cv_vap-cv_air)*cvm_ad(i)
1561  cvm_ad(i) = 0.0
1562  CALL poprealarray(cpm(i))
1563  cpm_ad(i) = 0.0
1564  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + q_liq_ad
1565  q0_ad(i, k, rainwat) = q0_ad(i, k, rainwat) + q_liq_ad
1566  END DO
1567  ELSE
1568  DO i=ie,is,-1
1569  temp_ad2 = cp_air*cpm_ad(i)
1570  CALL poprealarray(cvm(i))
1571  temp_ad1 = cv_air*cvm_ad(i)
1572  q_liq_ad = c_liq*cpm_ad(i) - temp_ad2 + c_liq*cvm_ad(i&
1573 & ) - temp_ad1
1574  q_sol_ad = c_ice*cpm_ad(i) - temp_ad2 + c_ice*cvm_ad(i&
1575 & ) - temp_ad1
1576  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + cp_vapor*&
1577 & cpm_ad(i) - temp_ad2 + cv_vap*cvm_ad(i) - temp_ad1
1578  cvm_ad(i) = 0.0
1579  CALL poprealarray(cpm(i))
1580  cpm_ad(i) = 0.0
1581  q0_ad(i, k, ice_wat) = q0_ad(i, k, ice_wat) + q_sol_ad
1582  q0_ad(i, k, liq_wat) = q0_ad(i, k, liq_wat) + q_liq_ad
1583  END DO
1584  END IF
1585  ELSE IF (branch .EQ. 3) THEN
1586  DO i=ie,is,-1
1587  CALL poprealarray(cvm(i))
1588  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + (cp_vapor-&
1589 & cp_air)*cpm_ad(i) + (cv_vap-cv_air)*cvm_ad(i)
1590  cvm_ad(i) = 0.0
1591  CALL poprealarray(cpm(i))
1592  cpm_ad(i) = 0.0
1593  END DO
1594  ELSE IF (branch .EQ. 4) THEN
1595  DO i=ie,is,-1
1596  CALL poprealarray(cvm(i))
1597  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + (cp_vapor-&
1598 & cp_air)*cpm_ad(i) + (cv_vap-cv_air)*cvm_ad(i)
1599  cvm_ad(i) = 0.0
1600  CALL poprealarray(cpm(i))
1601  cpm_ad(i) = 0.0
1602  END DO
1603  ELSE
1604  DO i=ie,is,-1
1605  CALL poprealarray(cvm(i))
1606  cvm_ad(i) = 0.0
1607  CALL poprealarray(cpm(i))
1608  cpm_ad(i) = 0.0
1609  END DO
1610  END IF
1611  END DO
1612  ELSE
1613  DO k=1,kbot,1
1614  DO i=ie,is,-1
1615  temp1 = pm(i, k)
1616  temp0 = pe(i, k, j)/temp1
1617  gz_ad(i, k) = gz_ad(i, k) + hd_ad(i, k)
1618  CALL poprealarray(gzh(i))
1619  tv_ad = (1.-temp0)*gz_ad(i, k) + (peln(i, k+1, j)-peln(i, &
1620 & k, j))*gzh_ad(i)
1621  peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + tv*gzh_ad(i)
1622  peln_ad(i, k, j) = peln_ad(i, k, j) - tv*gzh_ad(i)
1623  CALL poprealarray(hd(i, k))
1624  tvm_ad(i, k) = tvm_ad(i, k) + rdgas*tv_ad + cp_air*hd_ad(i&
1625 & , k)
1626  u0_ad(i, k) = u0_ad(i, k) + 0.5*2*u0(i, k)*hd_ad(i, k)
1627  v0_ad(i, k) = v0_ad(i, k) + 0.5*2*v0(i, k)*hd_ad(i, k)
1628  hd_ad(i, k) = 0.0
1629  CALL poprealarray(gz(i, k))
1630  temp_ad0 = -(tv*gz_ad(i, k)/temp1)
1631  gzh_ad(i) = gzh_ad(i) + gz_ad(i, k)
1632  pe_ad(i, k, j) = pe_ad(i, k, j) + temp_ad0
1633  pm_ad(i, k) = pm_ad(i, k) - temp0*temp_ad0
1634  gz_ad(i, k) = 0.0
1635  CALL poprealarray(tv)
1636  END DO
1637  END DO
1638  END IF
1639  DO i=ie,is,-1
1640  CALL poprealarray(gzh(i))
1641  gzh_ad(i) = 0.0
1642  END DO
1643  DO k=kbot,1,-1
1644  DO i=ie,is,-1
1645  CALL poprealarray(pm(i, k))
1646  temp = peln(i, k+1, j) - peln(i, k, j)
1647  temp_ad = -(delp(i, j, k)*pm_ad(i, k)/temp**2)
1648  delp_ad(i, j, k) = delp_ad(i, j, k) + pm_ad(i, k)/temp
1649  peln_ad(i, k+1, j) = peln_ad(i, k+1, j) + temp_ad
1650  peln_ad(i, k, j) = peln_ad(i, k, j) - temp_ad
1651  pm_ad(i, k) = 0.0
1652  CALL poprealarray(v0(i, k))
1653  va_ad(i, j, k) = va_ad(i, j, k) + v0_ad(i, k)
1654  v0_ad(i, k) = 0.0
1655  CALL poprealarray(u0(i, k))
1656  ua_ad(i, j, k) = ua_ad(i, j, k) + u0_ad(i, k)
1657  u0_ad(i, k) = 0.0
1658  t0_ad(i, k) = t0_ad(i, k) + (xvir*q0(i, k, sphum)+1.)*tvm_ad&
1659 & (i, k)
1660  q0_ad(i, k, sphum) = q0_ad(i, k, sphum) + t0(i, k)*xvir*&
1661 & tvm_ad(i, k)
1662  tvm_ad(i, k) = 0.0
1663  CALL poprealarray(t0(i, k))
1664  ta_ad(i, j, k) = ta_ad(i, j, k) + t0_ad(i, k)
1665  t0_ad(i, k) = 0.0
1666  END DO
1667  END DO
1668  DO iq=nq,1,-1
1669  DO k=kbot,1,-1
1670  DO i=ie,is,-1
1671  CALL poprealarray(q0(i, k, iq))
1672  qa_ad(i, j, k, iq) = qa_ad(i, j, k, iq) + q0_ad(i, k, iq)
1673  q0_ad(i, k, iq) = 0.0
1674  END DO
1675  END DO
1676  END DO
1677  END DO
1678  CALL popcontrol(1,branch)
1679  CALL popcontrol(1,branch)
1680  END IF
1681  END SUBROUTINE fv_subgrid_z_bwd
1682  SUBROUTINE fv_subgrid_z(isd, ied, jsd, jed, is, ie, js, je, km, nq, dt&
1683 & , tau, nwat, delp, pe, peln, pkz, ta, qa, ua, va, hydrostatic, w, &
1684 & delz, u_dt, v_dt, t_dt, k_bot)
1685  IMPLICIT NONE
1686 ! Dry convective adjustment-mixing
1687 !-------------------------------------------
1688  INTEGER, INTENT(IN) :: is, ie, js, je, km, nq, nwat
1689  INTEGER, INTENT(IN) :: isd, ied, jsd, jed
1690 ! Relaxation time scale
1691  INTEGER, INTENT(IN) :: tau
1692 ! model time step
1693  REAL, INTENT(IN) :: dt
1694  REAL, INTENT(IN) :: pe(is-1:ie+1, km+1, js-1:je+1)
1695  REAL, INTENT(IN) :: peln(is:ie, km+1, js:je)
1696 ! Delta p at each model level
1697  REAL, INTENT(IN) :: delp(isd:ied, jsd:jed, km)
1698 ! Delta z at each model level
1699  REAL, INTENT(IN) :: delz(isd:, jsd:, :)
1700  REAL, INTENT(IN) :: pkz(is:ie, js:je, km)
1701  LOGICAL, INTENT(IN) :: hydrostatic
1702  INTEGER, INTENT(IN), OPTIONAL :: k_bot
1703 !
1704  REAL, INTENT(INOUT) :: ua(isd:ied, jsd:jed, km)
1705  REAL, INTENT(INOUT) :: va(isd:ied, jsd:jed, km)
1706  REAL, INTENT(INOUT) :: w(isd:, jsd:, :)
1707 ! Temperature
1708  REAL, INTENT(INOUT) :: ta(isd:ied, jsd:jed, km)
1709 ! Specific humidity & tracers
1710  REAL, INTENT(INOUT) :: qa(isd:ied, jsd:jed, km, nq)
1711  REAL, INTENT(INOUT) :: u_dt(isd:ied, jsd:jed, km)
1712  REAL, INTENT(INOUT) :: v_dt(isd:ied, jsd:jed, km)
1713  REAL, INTENT(INOUT) :: t_dt(is:ie, js:je, km)
1714 !---------------------------Local variables-----------------------------
1715  REAL, DIMENSION(is:ie, km) :: u0, v0, w0, t0, hd, te, gz, tvm, pm, &
1716 & den
1717  REAL :: q0(is:ie, km, nq), qcon(is:ie, km)
1718  REAL, DIMENSION(is:ie) :: gzh, lcp2, icp2, cvm, cpm, qs
1719  REAL :: ri_ref, ri, pt1, pt2, ratio, tv, cv, tmp, q_liq, q_sol
1720  REAL :: tv1, tv2, g2, h0, mc, fra, rk, rz, rdt, tvd, tv_surf
1721  REAL :: dh, dq, qsw, dqsdt, tcp3, t_max, t_min
1722  INTEGER :: i, j, k, kk, n, m, iq, km1, im, kbot
1723  REAL, PARAMETER :: ustar2=1.e-4
1724  REAL :: cv_air, xvir
1725  INTEGER :: sphum, liq_wat, rainwat, snowwat, graupel, ice_wat, &
1726 & cld_amt
1727  INTRINSIC present
1728  INTRINSIC min
1729  INTRINSIC real
1730  INTRINSIC max
1731  INTEGER :: min1
1732  REAL :: max1
1733  REAL :: max2
1734  REAL :: y1
1735 ! = rdgas * (7/2-1) = 2.5*rdgas=717.68
1736  cv_air = cp_air - rdgas
1737  rk = cp_air/rdgas + 1.
1738  cv = cp_air - rdgas
1739  g2 = 0.5*grav
1740  rdt = 1./dt
1741  im = ie - is + 1
1742  IF (PRESENT(k_bot)) THEN
1743  IF (k_bot .LT. 3) THEN
1744  RETURN
1745  ELSE
1746  kbot = k_bot
1747  END IF
1748  ELSE
1749  kbot = km
1750  END IF
1751  IF (pe(is, 1, js) .LT. 2.) THEN
1752  t_min = t1_min
1753  ELSE
1754  t_min = t2_min
1755  END IF
1756  IF (km .GT. 24) THEN
1757  min1 = 24
1758  ELSE
1759  min1 = km
1760  END IF
1761  IF (k_bot .LT. min1) THEN
1762  t_max = t2_max
1763  ELSE
1764  t_max = t3_max
1765  END IF
1766  sphum = 1
1767  rainwat = -1
1768  snowwat = -1
1769  graupel = -1
1770  IF (nwat .EQ. 0) THEN
1771  xvir = 0.
1772  rz = 0.
1773  ELSE
1774  xvir = zvir
1775 ! rz = zvir * rdgas
1776  rz = rvgas - rdgas
1777  IF (nwat .EQ. 3) THEN
1778  liq_wat = 2
1779  ice_wat = 3
1780  END IF
1781  END IF
1782 !------------------------------------------------------------------------
1783 ! The nonhydrostatic pressure changes if there is heating (under constant
1784 ! volume and mass is locally conserved).
1785 !------------------------------------------------------------------------
1786  m = 3
1787  fra = dt/REAL(tau)
1788 !$OMP parallel do default(none) shared(im,is,ie,js,je,nq,kbot,qa,ta,sphum,ua,va,delp,peln, &
1789 !$OMP hydrostatic,pe,delz,g2,w,liq_wat,rainwat,ice_wat, &
1790 !$OMP snowwat,cv_air,m,graupel,pkz,rk,rz,fra, t_max, t_min, &
1791 !$OMP u_dt,rdt,v_dt,xvir,nwat) &
1792 !$OMP private(kk,lcp2,icp2,tcp3,dh,dq,den,qs,qsw,dqsdt,qcon,q0, &
1793 !$OMP t0,u0,v0,w0,h0,pm,gzh,tvm,tmp,cpm,cvm,q_liq,q_sol, &
1794 !$OMP tv,gz,hd,te,ratio,pt1,pt2,tv1,tv2,ri_ref, ri,mc,km1)
1795  DO j=js,je
1796  DO iq=1,nq
1797  DO k=1,kbot
1798  DO i=is,ie
1799  q0(i, k, iq) = qa(i, j, k, iq)
1800  END DO
1801  END DO
1802  END DO
1803  DO k=1,kbot
1804  DO i=is,ie
1805  t0(i, k) = ta(i, j, k)
1806  tvm(i, k) = t0(i, k)*(1.+xvir*q0(i, k, sphum))
1807  u0(i, k) = ua(i, j, k)
1808  v0(i, k) = va(i, j, k)
1809  pm(i, k) = delp(i, j, k)/(peln(i, k+1, j)-peln(i, k, j))
1810  END DO
1811  END DO
1812  DO i=is,ie
1813  gzh(i) = 0.
1814  END DO
1815  IF (hydrostatic) THEN
1816  DO k=kbot,1,-1
1817  DO i=is,ie
1818  tv = rdgas*tvm(i, k)
1819  den(i, k) = pm(i, k)/tv
1820  gz(i, k) = gzh(i) + tv*(1.-pe(i, k, j)/pm(i, k))
1821  hd(i, k) = cp_air*tvm(i, k) + gz(i, k) + 0.5*(u0(i, k)**2+v0&
1822 & (i, k)**2)
1823  gzh(i) = gzh(i) + tv*(peln(i, k+1, j)-peln(i, k, j))
1824  END DO
1825  END DO
1826  ELSE
1827  DO k=kbot,1,-1
1828  IF (nwat .EQ. 0) THEN
1829  DO i=is,ie
1830  cpm(i) = cp_air
1831  cvm(i) = cv_air
1832  END DO
1833  ELSE IF (nwat .EQ. 1) THEN
1834  DO i=is,ie
1835  cpm(i) = (1.-q0(i, k, sphum))*cp_air + q0(i, k, sphum)*&
1836 & cp_vapor
1837  cvm(i) = (1.-q0(i, k, sphum))*cv_air + q0(i, k, sphum)*&
1838 & cv_vap
1839  END DO
1840  ELSE IF (nwat .EQ. 2) THEN
1841 ! GFS
1842  DO i=is,ie
1843  cpm(i) = (1.-q0(i, k, sphum))*cp_air + q0(i, k, sphum)*&
1844 & cp_vapor
1845  cvm(i) = (1.-q0(i, k, sphum))*cv_air + q0(i, k, sphum)*&
1846 & cv_vap
1847  END DO
1848  ELSE IF (nwat .EQ. 3) THEN
1849  DO i=is,ie
1850  q_liq = q0(i, k, liq_wat)
1851  q_sol = q0(i, k, ice_wat)
1852  cpm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cp_air + q0(i&
1853 & , k, sphum)*cp_vapor + q_liq*c_liq + q_sol*c_ice
1854  cvm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cv_air + q0(i&
1855 & , k, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
1856  END DO
1857  ELSE IF (nwat .EQ. 4) THEN
1858  DO i=is,ie
1859  q_liq = q0(i, k, liq_wat) + q0(i, k, rainwat)
1860  cpm(i) = (1.-(q0(i, k, sphum)+q_liq))*cp_air + q0(i, k, &
1861 & sphum)*cp_vapor + q_liq*c_liq
1862  cvm(i) = (1.-(q0(i, k, sphum)+q_liq))*cv_air + q0(i, k, &
1863 & sphum)*cv_vap + q_liq*c_liq
1864  END DO
1865  ELSE
1866  DO i=is,ie
1867  q_liq = q0(i, k, liq_wat) + q0(i, k, rainwat)
1868  q_sol = q0(i, k, ice_wat) + q0(i, k, snowwat) + q0(i, k, &
1869 & graupel)
1870  cpm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cp_air + q0(i&
1871 & , k, sphum)*cp_vapor + q_liq*c_liq + q_sol*c_ice
1872  cvm(i) = (1.-(q0(i, k, sphum)+q_liq+q_sol))*cv_air + q0(i&
1873 & , k, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
1874  END DO
1875  END IF
1876  DO i=is,ie
1877  den(i, k) = -(delp(i, j, k)/(grav*delz(i, j, k)))
1878  w0(i, k) = w(i, j, k)
1879  gz(i, k) = gzh(i) - g2*delz(i, j, k)
1880  tmp = gz(i, k) + 0.5*(u0(i, k)**2+v0(i, k)**2+w0(i, k)**2)
1881  hd(i, k) = cpm(i)*t0(i, k) + tmp
1882  te(i, k) = cvm(i)*t0(i, k) + tmp
1883  gzh(i) = gzh(i) - grav*delz(i, j, k)
1884  END DO
1885  END DO
1886  END IF
1887  DO n=1,m
1888  IF (m .EQ. 3) THEN
1889  IF (n .EQ. 1) ratio = 0.25
1890  IF (n .EQ. 2) ratio = 0.5
1891  IF (n .EQ. 3) ratio = 0.999
1892  ELSE
1893  ratio = REAL(n)/REAL(m)
1894  END IF
1895  DO i=is,ie
1896  gzh(i) = 0.
1897  END DO
1898 ! Compute total condensate
1899  IF (nwat .LT. 2) THEN
1900  DO k=1,kbot
1901  DO i=is,ie
1902  qcon(i, k) = 0.
1903  END DO
1904  END DO
1905  ELSE IF (nwat .EQ. 2) THEN
1906 ! GFS_2015
1907  DO k=1,kbot
1908  DO i=is,ie
1909  qcon(i, k) = q0(i, k, liq_wat)
1910  END DO
1911  END DO
1912  ELSE IF (nwat .EQ. 3) THEN
1913  DO k=1,kbot
1914  DO i=is,ie
1915  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, ice_wat)
1916  END DO
1917  END DO
1918  ELSE IF (nwat .EQ. 4) THEN
1919  DO k=1,kbot
1920  DO i=is,ie
1921  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, rainwat)
1922  END DO
1923  END DO
1924  ELSE
1925  DO k=1,kbot
1926  DO i=is,ie
1927  qcon(i, k) = q0(i, k, liq_wat) + q0(i, k, ice_wat) + q0(i&
1928 & , k, snowwat) + q0(i, k, rainwat) + q0(i, k, graupel)
1929  END DO
1930  END DO
1931  END IF
1932  DO k=kbot,2,-1
1933  km1 = k - 1
1934  DO i=is,ie
1935 ! Richardson number = g*delz * del_theta/theta / (del_u**2 + del_v**2)
1936 ! Use exact form for "density temperature"
1937  tv1 = t0(i, km1)*(1.+xvir*q0(i, km1, sphum)-qcon(i, km1))
1938  tv2 = t0(i, k)*(1.+xvir*q0(i, k, sphum)-qcon(i, k))
1939  pt1 = tv1/pkz(i, j, km1)
1940  pt2 = tv2/pkz(i, j, k)
1941 !
1942  ri = (gz(i, km1)-gz(i, k))*(pt1-pt2)/(0.5*(pt1+pt2)*((u0(i, &
1943 & km1)-u0(i, k))**2+(v0(i, km1)-v0(i, k))**2+ustar2))
1944  IF (tv1 .GT. t_max .AND. tv1 .GT. tv2) THEN
1945 ! top layer unphysically warm
1946  ri = 0.
1947  ELSE IF (tv2 .LT. t_min) THEN
1948  IF (ri .GT. 0.2) THEN
1949  ri = 0.2
1950  ELSE
1951  ri = ri
1952  END IF
1953  END IF
1954  IF (400.e2 - pm(i, k) .LT. 0.) THEN
1955  max2 = 0.
1956  ELSE
1957  max2 = 400.e2 - pm(i, k)
1958  END IF
1959  y1 = ri_min + (ri_max-ri_min)*max2/200.e2
1960  IF (ri_max .GT. y1) THEN
1961  ri_ref = y1
1962  ELSE
1963  ri_ref = ri_max
1964  END IF
1965 ! Enhancing mixing at the model top
1966  IF (k .EQ. 2) THEN
1967  ri_ref = 4.*ri_ref
1968  ELSE IF (k .EQ. 3) THEN
1969  ri_ref = 2.*ri_ref
1970  ELSE IF (k .EQ. 4) THEN
1971  ri_ref = 1.5*ri_ref
1972  END IF
1973  IF (ri .LT. ri_ref) THEN
1974  IF (0.0 .LT. ri/ri_ref) THEN
1975  max1 = ri/ri_ref
1976  ELSE
1977  max1 = 0.0
1978  END IF
1979  mc = ratio*delp(i, j, km1)*delp(i, j, k)/(delp(i, j, km1)+&
1980 & delp(i, j, k))*(1.-max1)**2
1981  DO iq=1,nq
1982  h0 = mc*(q0(i, k, iq)-q0(i, km1, iq))
1983  q0(i, km1, iq) = q0(i, km1, iq) + h0/delp(i, j, km1)
1984  q0(i, k, iq) = q0(i, k, iq) - h0/delp(i, j, k)
1985  END DO
1986 ! Recompute qcon
1987  IF (nwat .LT. 2) THEN
1988  qcon(i, km1) = 0.
1989  ELSE IF (nwat .EQ. 2) THEN
1990 ! GFS_2015
1991  qcon(i, km1) = q0(i, km1, liq_wat)
1992  ELSE IF (nwat .EQ. 3) THEN
1993 ! AM3/AM4
1994  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, ice_wat)
1995  ELSE IF (nwat .EQ. 4) THEN
1996 ! K_warm_rain scheme with fake ice
1997  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, rainwat)
1998  ELSE
1999  qcon(i, km1) = q0(i, km1, liq_wat) + q0(i, km1, ice_wat)&
2000 & + q0(i, km1, snowwat) + q0(i, km1, rainwat) + q0(i, &
2001 & km1, graupel)
2002  END IF
2003 ! u:
2004  h0 = mc*(u0(i, k)-u0(i, k-1))
2005  u0(i, k-1) = u0(i, k-1) + h0/delp(i, j, k-1)
2006  u0(i, k) = u0(i, k) - h0/delp(i, j, k)
2007 ! v:
2008  h0 = mc*(v0(i, k)-v0(i, k-1))
2009  v0(i, k-1) = v0(i, k-1) + h0/delp(i, j, k-1)
2010  v0(i, k) = v0(i, k) - h0/delp(i, j, k)
2011  IF (hydrostatic) THEN
2012 ! Static energy
2013  h0 = mc*(hd(i, k)-hd(i, k-1))
2014  hd(i, k-1) = hd(i, k-1) + h0/delp(i, j, k-1)
2015  hd(i, k) = hd(i, k) - h0/delp(i, j, k)
2016  ELSE
2017 ! Total energy
2018  h0 = mc*(hd(i, k)-hd(i, k-1))
2019  te(i, k-1) = te(i, k-1) + h0/delp(i, j, k-1)
2020  te(i, k) = te(i, k) - h0/delp(i, j, k)
2021 ! w:
2022  h0 = mc*(w0(i, k)-w0(i, k-1))
2023  w0(i, k-1) = w0(i, k-1) + h0/delp(i, j, k-1)
2024  w0(i, k) = w0(i, k) - h0/delp(i, j, k)
2025  END IF
2026  END IF
2027  END DO
2028 !--------------
2029 ! Retrive Temp:
2030 !--------------
2031  IF (hydrostatic) THEN
2032  kk = k
2033  DO i=is,ie
2034  t0(i, kk) = (hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, kk)&
2035 & **2))/(rk-pe(i, kk, j)/pm(i, kk))
2036  gzh(i) = gzh(i) + t0(i, kk)*(peln(i, kk+1, j)-peln(i, kk, &
2037 & j))
2038  t0(i, kk) = t0(i, kk)/(rdgas+rz*q0(i, kk, sphum))
2039  END DO
2040  kk = k - 1
2041  DO i=is,ie
2042  t0(i, kk) = (hd(i, kk)-gzh(i)-0.5*(u0(i, kk)**2+v0(i, kk)&
2043 & **2))/((rk-pe(i, kk, j)/pm(i, kk))*(rdgas+rz*q0(i, kk, &
2044 & sphum)))
2045  END DO
2046  ELSE
2047 ! Non-hydrostatic under constant volume heating/cooling
2048  DO kk=k-1,k
2049  IF (nwat .EQ. 0) THEN
2050  DO i=is,ie
2051  cpm(i) = cp_air
2052  cvm(i) = cv_air
2053  END DO
2054  ELSE IF (nwat .EQ. 1) THEN
2055  DO i=is,ie
2056  cpm(i) = (1.-q0(i, kk, sphum))*cp_air + q0(i, kk, &
2057 & sphum)*cp_vapor
2058  cvm(i) = (1.-q0(i, kk, sphum))*cv_air + q0(i, kk, &
2059 & sphum)*cv_vap
2060  END DO
2061  ELSE IF (nwat .EQ. 2) THEN
2062  DO i=is,ie
2063  cpm(i) = (1.-q0(i, kk, sphum))*cp_air + q0(i, kk, &
2064 & sphum)*cp_vapor
2065  cvm(i) = (1.-q0(i, kk, sphum))*cv_air + q0(i, kk, &
2066 & sphum)*cv_vap
2067  END DO
2068  ELSE IF (nwat .EQ. 3) THEN
2069  DO i=is,ie
2070  q_liq = q0(i, kk, liq_wat)
2071  q_sol = q0(i, kk, ice_wat)
2072  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cp_air + &
2073 & q0(i, kk, sphum)*cp_vapor + q_liq*c_liq + q_sol*&
2074 & c_ice
2075  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cv_air + &
2076 & q0(i, kk, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
2077  END DO
2078  ELSE IF (nwat .EQ. 4) THEN
2079  DO i=is,ie
2080  q_liq = q0(i, kk, liq_wat) + q0(i, kk, rainwat)
2081  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq))*cp_air + q0(i, &
2082 & kk, sphum)*cp_vapor + q_liq*c_liq
2083  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq))*cv_air + q0(i, &
2084 & kk, sphum)*cv_vap + q_liq*c_liq
2085  END DO
2086  ELSE
2087  DO i=is,ie
2088  q_liq = q0(i, kk, liq_wat) + q0(i, kk, rainwat)
2089  q_sol = q0(i, kk, ice_wat) + q0(i, kk, snowwat) + q0(i&
2090 & , kk, graupel)
2091  cpm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cp_air + &
2092 & q0(i, kk, sphum)*cp_vapor + q_liq*c_liq + q_sol*&
2093 & c_ice
2094  cvm(i) = (1.-(q0(i, kk, sphum)+q_liq+q_sol))*cv_air + &
2095 & q0(i, kk, sphum)*cv_vap + q_liq*c_liq + q_sol*c_ice
2096  END DO
2097  END IF
2098  DO i=is,ie
2099  tv = gz(i, kk) + 0.5*(u0(i, kk)**2+v0(i, kk)**2+w0(i, kk&
2100 & )**2)
2101  t0(i, kk) = (te(i, kk)-tv)/cvm(i)
2102  hd(i, kk) = cpm(i)*t0(i, kk) + tv
2103  END DO
2104  END DO
2105  END IF
2106  END DO
2107  END DO
2108 ! k-loop
2109 ! n-loop
2110 !--------------------
2111  IF (fra .LT. 1.) THEN
2112  DO k=1,kbot
2113  DO i=is,ie
2114  t0(i, k) = ta(i, j, k) + (t0(i, k)-ta(i, j, k))*fra
2115  u0(i, k) = ua(i, j, k) + (u0(i, k)-ua(i, j, k))*fra
2116  v0(i, k) = va(i, j, k) + (v0(i, k)-va(i, j, k))*fra
2117  END DO
2118  END DO
2119  IF (.NOT.hydrostatic) THEN
2120  DO k=1,kbot
2121  DO i=is,ie
2122  w0(i, k) = w(i, j, k) + (w0(i, k)-w(i, j, k))*fra
2123  END DO
2124  END DO
2125  END IF
2126  DO iq=1,nq
2127  DO k=1,kbot
2128  DO i=is,ie
2129  q0(i, k, iq) = qa(i, j, k, iq) + (q0(i, k, iq)-qa(i, j, k&
2130 & , iq))*fra
2131  END DO
2132  END DO
2133  END DO
2134  END IF
2135  DO k=1,kbot
2136  DO i=is,ie
2137  u_dt(i, j, k) = rdt*(u0(i, k)-ua(i, j, k))
2138  v_dt(i, j, k) = rdt*(v0(i, k)-va(i, j, k))
2139 ! *** temperature updated ***
2140  ta(i, j, k) = t0(i, k)
2141  ua(i, j, k) = u0(i, k)
2142  va(i, j, k) = v0(i, k)
2143  END DO
2144  DO iq=1,nq
2145  DO i=is,ie
2146  qa(i, j, k, iq) = q0(i, k, iq)
2147  END DO
2148  END DO
2149  END DO
2150  IF (.NOT.hydrostatic) THEN
2151  DO k=1,kbot
2152  DO i=is,ie
2153 ! w updated
2154  w(i, j, k) = w0(i, k)
2155  END DO
2156  END DO
2157  END IF
2158  END DO
2159  END SUBROUTINE fv_subgrid_z
2160  SUBROUTINE neg_adj3(is, ie, js, je, ng, kbot, hydrostatic, peln, delz&
2161 & , pt, dp, qv, ql, qr, qi, qs, qg, qa, check_negative)
2162  IMPLICIT NONE
2163 !Stubbed, would require careful checking of nonlinearity
2164 ! This is designed for 6-class micro-physics schemes
2165  INTEGER, INTENT(IN) :: is, ie, js, je, ng, kbot
2166  LOGICAL, INTENT(IN) :: hydrostatic
2167 ! total delp-p
2168  REAL, INTENT(IN) :: dp(is-ng:ie+ng, js-ng:je+ng, kbot)
2169  REAL, INTENT(IN) :: delz(is-ng:, js-ng:, :)
2170 ! ln(pe)
2171  REAL, INTENT(IN) :: peln(is:ie, kbot+1, js:je)
2172  LOGICAL, INTENT(IN), OPTIONAL :: check_negative
2173  REAL, DIMENSION(is-ng:ie+ng, js-ng:je+ng, kbot), INTENT(INOUT) :: pt&
2174 & , qv, ql, qr, qi, qs, qg
2175  REAL, DIMENSION(is-ng:ie+ng, js-ng:je+ng, kbot), INTENT(INOUT), &
2176 & OPTIONAL :: qa
2177  end subroutine neg_adj3
2178 
2179 end module fv_sg_adm_mod
fv_sg_adm_mod::neg_adj3
subroutine, public neg_adj3(is, ie, js, je, ng, kbot, hydrostatic, peln, delz, pt, dp, qv, ql, qr, qi, qs, qg, qa, check_negative)
Definition: fv_sg_adm.F90:2162
tracer_manager_mod::get_tracer_index
Definition: tracer_manager.F90:124
fv_mp_nlm_mod
Definition: fv_mp_nlm_mod.F90:24
constants_mod
Definition: constants.F90:24
lin_cld_microphys_mod::wqs2
real function, public wqs2(ta, den, dqdt)
Definition: lin_cloud_microphys_nlm.F90:3569
tapenade_iter::poprealarray_adm
Definition: tapenade_iter.F90:186
field_manager_mod::model_atmos
integer, parameter, public model_atmos
Definition: field_manager.F90:299
fv_sg_adm_mod::c_liq
real, parameter c_liq
Definition: fv_sg_adm.F90:43
fv_sg_adm_mod::dc_ice
real, parameter dc_ice
Definition: fv_sg_adm.F90:50
tapenade_iter::pushrealarray
Definition: tapenade_iter.F90:146
fv_sg_adm_mod::li0
real, parameter li0
Definition: fv_sg_adm.F90:64
fv_sg_adm_mod
Definition: fv_sg_adm.F90:20
field_manager_mod
Definition: field_manager.F90:20
fv_sg_adm_mod::lv0
real, parameter lv0
Definition: fv_sg_adm.F90:63
tapenade_iter::pushinteger
Definition: tapenade_iter.F90:136
constants_mod::hlv
real, parameter, public hlv
Latent heat of evaporation [J/kg].
Definition: constants.F90:80
fv_sg_adm_mod::tice
real, parameter tice
Definition: fv_sg_adm.F90:40
fv_sg_adm_mod::ri_max
real, parameter ri_max
Definition: fv_sg_adm.F90:57
tapenade_iter::pushcontrol
subroutine, public pushcontrol(ctype, field)
Definition: tapenade_iter.F90:861
constants_mod::rdgas
real, parameter, public rdgas
Gas constant for dry air [J/kg/deg].
Definition: constants.F90:77
fv_sg_adm_mod::fv_subgrid_z
subroutine, public fv_subgrid_z(isd, ied, jsd, jed, is, ie, js, je, km, nq, dt, tau, nwat, delp, pe, peln, pkz, ta, qa, ua, va, hydrostatic, w, delz, u_dt, v_dt, t_dt, k_bot)
Definition: fv_sg_adm.F90:1685
constants_mod::cp_vapor
real, parameter, public cp_vapor
Specific heat capacity of water vapor at constant pressure [J/kg/deg].
Definition: constants.F90:89
fv_sg_adm_mod::fv_subgrid_z_fwd
subroutine, public fv_subgrid_z_fwd(isd, ied, jsd, jed, is, ie, js, je, km, nq, dt, tau, nwat, delp, pe, peln, pkz, ta, qa, ua, va, hydrostatic, w, delz, u_dt, v_dt, t_dt, k_bot)
Definition: fv_sg_adm.F90:95
lin_cld_microphys_mod
Definition: lin_cloud_microphys_nlm.F90:7
fv_sg_adm_mod::c_con
real, parameter c_con
Definition: fv_sg_adm.F90:46
fv_sg_adm_mod::d0_vap
real d0_vap
Definition: fv_sg_adm.F90:68
constants_mod::rvgas
real, parameter, public rvgas
Gas constant for water vapor [J/kg/deg].
Definition: constants.F90:78
tapenade_iter::pushrealarray_adm
Definition: tapenade_iter.F90:173
fv_sg_adm_mod::dc_vap
real, parameter dc_vap
Definition: fv_sg_adm.F90:49
tapenade_iter
Definition: tapenade_iter.F90:1
constants_mod::cp_air
real, parameter, public cp_air
Specific heat capacity of dry air at constant pressure [J/kg/deg].
Definition: constants.F90:83
fv_sg_adm_mod::t1_min
real, parameter t1_min
Definition: fv_sg_adm.F90:59
fv_sg_adm_mod::hlv0
real, parameter hlv0
Definition: fv_sg_adm.F90:52
tapenade_iter::poprealarray
Definition: tapenade_iter.F90:159
fv_sg_adm_mod::zvir
real, parameter zvir
Definition: fv_sg_adm.F90:66
fv_sg_adm_mod::des
real, dimension(:), allocatable des
Definition: fv_sg_adm.F90:67
fv_sg_adm_mod::t2_min
real, parameter t2_min
Definition: fv_sg_adm.F90:60
constants_mod::hlf
real, parameter, public hlf
Latent heat of fusion [J/kg].
Definition: constants.F90:81
constants_mod::grav
real, parameter, public grav
Acceleration due to gravity [m/s^2].
Definition: constants.F90:76
tracer_manager_mod
Definition: tracer_manager.F90:20
fv_sg_adm_mod::lv00
real lv00
Definition: fv_sg_adm.F90:68
fv_sg_adm_mod::t_ice
real, parameter t_ice
Definition: fv_sg_adm.F90:56
fv_sg_adm_mod::t2_max
real, parameter t2_max
Definition: fv_sg_adm.F90:61
fv_sg_adm_mod::hlf0
real, parameter hlf0
Definition: fv_sg_adm.F90:53
max
#define max(a, b)
Definition: mosaic_util.h:33
fv_sg_adm_mod::t3_max
real, parameter t3_max
Definition: fv_sg_adm.F90:62
fv_sg_adm_mod::cv_vap
real, parameter cv_vap
Definition: fv_sg_adm.F90:45
fv_sg_adm_mod::fv_subgrid_z_bwd
subroutine, public fv_subgrid_z_bwd(isd, ied, jsd, jed, is, ie, js, je, km, nq, dt, tau, nwat, delp, delp_ad, pe, pe_ad, peln, peln_ad, pkz, pkz_ad, ta, ta_ad, qa, qa_ad, ua, ua_ad, va, va_ad, hydrostatic, w, w_ad, delz, delz_ad, u_dt, v_dt, t_dt, k_bot)
Definition: fv_sg_adm.F90:786
fv_sg_adm_mod::table
real, dimension(:), allocatable table
Definition: fv_sg_adm.F90:67
fv_sg_adm_mod::ri_min
real, parameter ri_min
Definition: fv_sg_adm.F90:58
fv_sg_adm_mod::esl
real, parameter esl
Definition: fv_sg_adm.F90:39
min
#define min(a, b)
Definition: mosaic_util.h:32
lin_cld_microphys_mod::wqsat2_moist
real function, public wqsat2_moist(ta, qv, pa, dqdt)
Definition: lin_cloud_microphys_nlm.F90:3666
constants_mod::kappa
real, parameter, public kappa
RDGAS / CP_AIR [dimensionless].
Definition: constants.F90:82
tapenade_iter::popcontrol
subroutine, public popcontrol(ctype, field)
Definition: tapenade_iter.F90:922
fv_sg_adm_mod::c_ice
real, parameter c_ice
Definition: fv_sg_adm.F90:42
tapenade_iter::popinteger
Definition: tapenade_iter.F90:141