module_chem_utilities.F
References to this file elsewhere.
1 MODULE module_chem_utilities
2 USE module_domain
3 USE module_model_constants
4 USE module_state_description
5 USE module_configure
6
7 CONTAINS
8 SUBROUTINE chem_prep ( config_flags, & ! input
9 u, v, p, pb, alt, ph, & ! input
10 phb, t, moist, n_moist, & ! input
11 rho, p_phy , & ! output
12 u_phy, v_phy, p8w, t_phy, t8w, & ! output
13 z, z_at_w, dz8w, & ! output
14 fzm, fzp, & ! params
15 ids, ide, jds, jde, kds, kde, &
16 ims, ime, jms, jme, kms, kme, &
17 its, ite, jts, jte, kts, kte )
18 !----------------------------------------------------------------------
19 IMPLICIT NONE
20 !----------------------------------------------------------------------
21
22 TYPE(grid_config_rec_type) , INTENT(IN ) :: config_flags
23 INTEGER , INTENT(IN ) :: ids, ide, jds, jde, kds, kde, &
24 ims, ime, jms, jme, kms, kme, &
25 its, ite, jts, jte, kts, kte
26 INTEGER , INTENT(IN ) :: n_moist
27
28 REAL, DIMENSION( ims:ime, kms:kme , jms:jme , n_moist ), INTENT(IN) :: moist
29
30
31
32 REAL , DIMENSION( ims:ime , kms:kme , jms:jme ) , &
33 INTENT( OUT) :: u_phy, &
34 v_phy, &
35 p_phy, &
36 p8w, &
37 t_phy, &
38 t8w, &
39 rho, &
40 z, &
41 dz8w, &
42 z_at_w
43
44 REAL , DIMENSION( ims:ime , kms:kme , jms:jme ) , &
45 INTENT(IN ) :: pb, &
46 p, &
47 u, &
48 v, &
49 alt, &
50 ph, &
51 phb, &
52 t
53
54
55 REAL , DIMENSION( kms:kme ) , INTENT(IN ) :: fzm, &
56 fzp
57
58 INTEGER :: i_start, i_end, j_start, j_end, k_start, k_end
59 INTEGER :: i, j, k
60 REAL :: w1, w2, z0, z1, z2
61
62 !-----------------------------------------------------------------------
63
64 ! set up loop bounds for this grid's boundary conditions
65
66 i_start = its
67 i_end = min( ite,ide-1 )
68 j_start = jts
69 j_end = min( jte,jde-1 )
70
71 k_start = kts
72 k_end = min( kte, kde-1 )
73
74 ! compute thermodynamics and velocities at pressure points
75 do j = j_start,j_end
76 do k = k_start, k_end
77 do i = i_start, i_end
78
79 p_phy(i,k,j) = p(i,k,j) + pb(i,k,j)
80 t_phy(i,k,j) = (t(i,k,j)+t0)*(p_phy(i,k,j)/p1000mb)**rcp
81 rho(i,k,j) = 1./alt(i,k,j)*(1.+moist(i,k,j,P_QV))
82 u_phy(i,k,j) = 0.5*(u(i,k,j)+u(i+1,k,j))
83 v_phy(i,k,j) = 0.5*(v(i,k,j)+v(i,k,j+1))
84
85 enddo
86 enddo
87 enddo
88
89 ! wig: added to make sure there is no junk in the top level even
90 ! though it should not be used
91 do j = j_start,j_end
92 do i = i_start, i_end
93 p_phy(i,kte,j) = p_phy(i,k_end,j)
94 t_phy(i,kte,j) = t_phy(i,k_end,j)
95 rho(i,kte,j) = rho(i,k_end,j)
96 u_phy(i,kte,j) = u_phy(i,k_end,j)
97 v_phy(i,kte,j) = v_phy(i,k_end,j)
98 enddo
99 enddo
100
101 ! compute z at w points
102
103 do j = j_start,j_end
104 do k = k_start, kte
105 do i = i_start, i_end
106 z_at_w(i,k,j) = (phb(i,k,j)+ph(i,k,j))/g
107 enddo
108 enddo
109 enddo
110
111 do j = j_start,j_end
112 do k = k_start, kte-1
113 do i = i_start, i_end
114 dz8w(i,k,j) = z_at_w(i,k+1,j)-z_at_w(i,k,j)
115 enddo
116 enddo
117 enddo
118
119 do j = j_start,j_end
120 do i = i_start, i_end
121 dz8w(i,kte,j) = 0.
122 enddo
123 enddo
124
125 ! compute z at p points (average of z at w points)
126 do j = j_start,j_end
127 do k = k_start, k_end
128 do i = i_start, i_end
129 z(i,k,j) = 0.5*(z_at_w(i,k,j) +z_at_w(i,k+1,j) )
130 enddo
131 enddo
132 enddo
133
134 ! interp t and p at w points
135
136 do j = j_start,j_end
137 do k = 2, k_end
138 do i = i_start, i_end
139 p8w(i,k,j) = fzm(k)*p_phy(i,k,j)+fzp(k)*p_phy(i,k-1,j)
140 t8w(i,k,j) = fzm(k)*t_phy(i,k,j)+fzp(k)*t_phy(i,k-1,j)
141 enddo
142 enddo
143 enddo
144
145 ! extrapolate p and t to surface and top.
146 ! we'll use an extrapolation in z for now
147
148 do j = j_start,j_end
149 do i = i_start, i_end
150
151 ! bottom
152
153 z0 = z_at_w(i,1,j)
154 z1 = z(i,1,j)
155 z2 = z(i,2,j)
156 w1 = (z0 - z2)/(z1 - z2)
157 w2 = 1. - w1
158 p8w(i,1,j) = w1*p_phy(i,1,j)+w2*p_phy(i,2,j)
159 t8w(i,1,j) = w1*t_phy(i,1,j)+w2*t_phy(i,2,j)
160
161 ! top
162
163 z0 = z_at_w(i,kte,j)
164 z1 = z(i,k_end,j)
165 z2 = z(i,k_end-1,j)
166 w1 = (z0 - z2)/(z1 - z2)
167 w2 = 1. - w1
168
169 ! p8w(i,kde,j) = w1*p_phy(i,kde-1,j)+w2*p_phy(i,kde-2,j)
170 !!! bug fix extrapolate ln(p) so p is positive definite
171 p8w(i,kde,j) = exp(w1*log(p_phy(i,kde-1,j))+w2*log(p_phy(i,kde-2,j)))
172 t8w(i,kde,j) = w1*t_phy(i,kde-1,j)+w2*t_phy(i,kde-2,j)
173
174 enddo
175 enddo
176 END SUBROUTINE chem_prep
177 END MODULE module_chem_utilities