module_cbmz_initmixrats.F
References to this file elsewhere.
1 !**********************************************************************************
2 ! This computer software was prepared by Battelle Memorial Institute, hereinafter
3 ! the Contractor, under Contract No. DE-AC05-76RL0 1830 with the Department of
4 ! Energy (DOE). NEITHER THE GOVERNMENT NOR THE CONTRACTOR MAKES ANY WARRANTY,
5 ! EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF THIS SOFTWARE.
6 !
7 ! CBMZ module: see module_cbmz.F for information and terms of use
8 !**********************************************************************************
9 module module_cbmz_initmixrats
10
11
12 use module_peg_util
13
14
15 integer, parameter :: cbmz_init_wrf_mixrats_flagaa = 1
16 ! turns subr cbmz_init_wrf_mixrats on/off
17
18
19 contains
20
21
22 !-----------------------------------------------------------------------
23 subroutine cbmz_init_wrf_mixrats( &
24 config_flags, &
25 z_at_w, g, &
26 chem, numgas, &
27 ids,ide, jds,jde, kds,kde, &
28 ims,ime, jms,jme, kms,kme, &
29 its,ite, jts,jte, kts,kte )
30 !
31 ! provides initial values for cbmz gas species
32 ! for gas species that are common to both cbmz and radm2, the initial
33 ! values are provided via the run initialization file
34 ! (The radm2 gas species are initialized from this file
35 ! when chem_in_opt==anything. This ought to be changed!)
36 ! for gas species that are in cbmz but not in radm2, the initial values
37 ! are provided here
38 ! currently only hcl and "par" have non-zero initial values,
39 ! and other species are near-zero
40 !
41 ! when (gas_ic_opt == gas_ic_pnnl) AND (chem_in_opt == 0),
42 ! ozone is set to "Texas August 2000" values
43 !
44 ! setting cbmz_init_wrf_mixrats_flagaa = 1/0 turns this subr on/off.
45 !
46
47 USE module_configure, only: grid_config_rec_type, num_chem, &
48 p_o3, p_ald, p_hc3, p_hc5, p_hc8, p_ket, p_oli, p_olt, p_ora2, &
49 p_hcl, p_par
50 USE module_state_description, only: param_first_scalar, &
51 gas_ic_pnnl
52 USE module_input_chem_data, only: bdy_chem_value
53
54 IMPLICIT NONE
55
56
57 !-----------------------------------------------------------------------
58 ! subr arguments
59
60 INTEGER, INTENT(IN) :: numgas, &
61 ids,ide, jds,jde, kds,kde, &
62 ims,ime, jms,jme, kms,kme, &
63 its,ite, jts,jte, kts,kte
64
65 REAL, INTENT(IN) :: g
66
67 ! perturbation and base geopotential at layer boundaries
68 REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), &
69 INTENT(IN) :: z_at_w
70
71 ! advected chemical tracers
72 REAL, DIMENSION( ims:ime, kms:kme, jms:jme, num_chem ), &
73 INTENT(INOUT) :: chem
74
75 TYPE(grid_config_rec_type), INTENT(IN) :: config_flags
76 !-----------------------------------------------------------------------
77
78 ! local variables
79 integer i, j, k, kp1
80 real, dimension( its:ite, kts:kte, jts:jte ) :: z
81
82 if (cbmz_init_wrf_mixrats_flagaa <= 0) return
83
84 !
85 ! calculate the mid-level heights
86 !
87 do j = jts, min(jte,jde-1)
88 do k = kts, kte
89 kp1 = min(k+1, kte)
90 do i = its, min(ite,ide-1)
91 z(i,k,j) = (z_at_w(i,k,j)+z_at_w(i,kp1,j))*0.5
92 end do
93 end do
94 end do
95
96 !
97 ! when (gas_ic_opt == gas_ic_pnnl) AND (chem_in_opt == 0),
98 ! set ozone (and other gas?) to "Texas August 2000" values
99 !
100 if ( (config_flags%chem_in_opt == 0) .and. &
101 (config_flags%gas_ic_opt == gas_ic_pnnl) ) then
102 do j = jts, min(jte,jde-1)
103 do k = kts, kte
104 do i = its, min(ite,ide-1)
105 call bdy_chem_value( chem(i,k,j,p_o3),z(i,k,j), p_o3, numgas )
106 end do
107 end do
108 end do
109 end if
110
111 !
112 ! compute hcl initial mixing ratio based on the article:
113 ! Graedel TE and WC Keene, 1995: Troposhperic budget of reactive chlorine.
114 ! Global Biogeochemical Cycles. 9, (1), 47-77.
115 ! This calculation should mimic the hcl profile in bdy_chem_value_cbmz,
116 ! below.
117 !
118 ! Height(m) HCl concentration
119 ! --------- -----------------
120 ! <=1000 0.4 ppbv
121 ! 1000<z<2500 linear transision zone
122 ! >=2500 0.1 ppbv
123 !
124 do j = jts, min(jte,jde-1)
125 do k = kts, kte
126 do i = its, min(ite,ide-1)
127 if( z(i,k,j) <= 1000. ) then
128 chem(i,k,j,p_hcl) = 0.4*1e-3
129 elseif( z(i,k,j) > 1000. &
130 .and. z(i,k,j) <= 2500. ) then
131 chem(i,k,j,p_hcl) = (0.4*1e-3) + (z(i,k,j)-1000.)* &
132 ((0.1*1e-3)-(0.4*1e-3)) / (2500.-1000.)
133 else
134 chem(i,k,j,p_hcl) = 0.1*1e-3
135 end if
136 end do
137 end do
138 end do
139
140 !wig, 2-May-2007: Moved this logic to make_chem_profile so this species
141 ! will now come in via the real.exe generated wrfinput.
142 !!$!
143 !!$! compute par initial mixing ratio from radm2 hydrocarbon species
144 !!$! using same formula as for par emissions
145 !!$!
146 !!$ do j = jts, min(jte,jde-1)
147 !!$ do k = kts, kte
148 !!$ do i = its, min(ite,ide-1)
149 !!$ chem(i,k,j,p_par) = &
150 !!$ 0.4*chem(i,k,j,p_ald) + 2.9*chem(i,k,j,p_hc3) &
151 !!$ + 4.8*chem(i,k,j,p_hc5) + 7.9*chem(i,k,j,p_hc8) &
152 !!$ + 0.9*chem(i,k,j,p_ket) + 2.8*chem(i,k,j,p_oli) &
153 !!$ + 1.8*chem(i,k,j,p_olt) + 1.0*chem(i,k,j,p_ora2)
154 !!$ end do
155 !!$ end do
156 !!$ end do
157
158 return
159 end subroutine cbmz_init_wrf_mixrats
160
161
162
163 !-----------------------------------------------------------------------
164 end module module_cbmz_initmixrats
165
166
167
168 !-----------------------------------------------------------------------
169 subroutine bdy_chem_value_cbmz ( chem_bv, z, nch, numgas )
170 !
171 ! provides boundary values for cbmz gas species
172 ! for gas species that are common to both cbmz and radm2, the boundary
173 ! values are provided by subr bdy_chem_value
174 ! for gas species that are in cbmz but not in radm2, the boundary values
175 ! are provided here, except for par
176 ! currently only "hcl" has a non-zero boundary value,
177 ! and other species are near-zero
178 !
179 ! this is outside of the module declaration because of potential
180 ! module1 --> module2 --> module1 use conflicts
181 !
182 use module_configure, only: grid_config_rec_type, &
183 p_o3, p_ald, p_hc3, p_hc5, p_hc8, p_ket, p_oli, &
184 p_olt, p_ora2, p_hcl, p_par
185 use module_input_chem_data, only: bdy_chem_value
186
187 implicit none
188
189 ! arguments
190 REAL, INTENT(OUT) :: chem_bv ! boundary value for chem(-,-,-,nch)
191 REAL, INTENT(IN) :: z ! height
192 INTEGER, INTENT(IN) :: nch ! index number of chemical species
193 INTEGER, INTENT(IN) :: numgas ! index number of last gas species
194 ! local variables
195 real chem_bv_ald, chem_bv_hc3, chem_bv_hc5, &
196 chem_bv_hc8, chem_bv_ket, chem_bv_oli, &
197 chem_bv_olt, chem_bv_ora2
198 real, parameter :: chem_bv_def = 1.0e-20
199
200
201 if( nch == p_hcl ) then
202 !This calculation should mimic the hcl profile in
203 !cbmz_init_wrf_mixrats, above.
204 if( z <= 1000. ) then
205 chem_bv = 0.4*1e-3
206 elseif( z > 1000. &
207 .and. z <= 2500. ) then
208 chem_bv = (0.4*1e-3) + (z-1000.)* &
209 ((0.1*1e-3)-(0.4*1e-3)) / (2500.-1000.)
210 else
211 chem_bv = 0.1*1e-3
212 end if
213
214 else if( nch == p_par ) then
215 call bdy_chem_value( chem_bv_ald, z, p_ald, numgas )
216 !wig, 2-May-2007: begin
217 ! The extra +1 offsets the blank first index for p_XXX
218 call bdy_chem_value( chem_bv_hc3, z, numgas+1+1, numgas )
219 call bdy_chem_value( chem_bv_hc5, z, numgas+1+2, numgas )
220 call bdy_chem_value( chem_bv_hc8, z, numgas+1+3, numgas )
221 !wig, 2-May-2007: end
222 call bdy_chem_value( chem_bv_ket, z, p_ket, numgas )
223 call bdy_chem_value( chem_bv_oli, z, p_oli, numgas )
224 call bdy_chem_value( chem_bv_olt, z, p_olt, numgas )
225 call bdy_chem_value( chem_bv_ora2, z, p_ora2, numgas )
226
227 chem_bv = 0.4*chem_bv_ald + 2.9*chem_bv_hc3 &
228 + 4.8*chem_bv_hc5 + 7.9*chem_bv_hc8 &
229 + 0.9*chem_bv_ket + 2.8*chem_bv_oli &
230 + 1.8*chem_bv_olt + 1.0*chem_bv_ora2
231
232 else
233 ! chem_bv=0 for all other species
234 chem_bv = chem_bv_def
235
236 end if
237
238 return
239 end subroutine bdy_chem_value_cbmz