Line data Source code
1 : !================================================================================
2 : ! manages mapping of CLM generated wild fire emissions to chemical constituents
3 : !================================================================================
4 : module fire_emissions
5 :
6 : use shr_kind_mod, only : r8 => shr_kind_r8
7 : use shr_fire_emis_mod, only : shr_fire_emis_mechcomps, shr_fire_emis_mechcomps_n, shr_fire_emis_elevated
8 : use srf_field_check, only : active_Fall_flxfire
9 : use shr_const_mod, only : pi => SHR_CONST_PI
10 : use shr_const_mod, only : avogad => SHR_CONST_AVOGAD ! Avogadro's number ~ molecules/kmole
11 : use cam_abortutils, only : endrun
12 : use cam_history, only : addfld, horiz_only, outfld, fieldname_len
13 : use cam_logfile, only : iulog
14 : use ppgrid, only : pver, pverp
15 : use constituents, only : cnst_get_ind
16 : use rad_constituents, only : rad_cnst_get_aer_props, rad_cnst_num_name
17 : use mo_chem_utls, only : get_spc_ndx, get_extfrc_ndx
18 : use chem_mods, only : adv_mass ! g/mole
19 : use infnan, only : nan, assignment(=)
20 :
21 : implicit none
22 : private
23 : save
24 :
25 : public :: fire_emissions_init
26 : public :: fire_emissions_srf
27 : public :: fire_emissions_vrt
28 :
29 : ! for surface emissions
30 : integer, allocatable :: fire_emis_indices_map(:)
31 :
32 : ! for vertically distributed forcings
33 : integer, allocatable :: frc_spc_map(:)
34 : integer, allocatable :: chm_spc_map(:)
35 : integer, allocatable :: frc_num_map(:)
36 : real(r8), allocatable :: spc_mass_factor(:)
37 : real(r8), allocatable :: num_mass_factor(:)
38 : character(len=fieldname_len), allocatable :: fire_frc_name(:)
39 : character(len=fieldname_len), allocatable :: fire_numfrc_name(:)
40 : character(len=fieldname_len), allocatable :: fire_sflx_name(:)
41 : character(len=fieldname_len), allocatable :: fire_vflx_name(:)
42 :
43 : !================================================================================
44 : contains
45 : !================================================================================
46 :
47 : !------------------------------------------------------------------------------
48 : !------------------------------------------------------------------------------
49 0 : subroutine fire_emissions_init()
50 :
51 : ! local vars
52 : integer :: n, ii
53 :
54 : integer :: frc_ndx, spc_ndx
55 : integer :: mode, spec
56 : character(len=16) :: name
57 : character(len=32) :: num_name
58 :
59 : real(r8), parameter :: demis_acc = 0.134e-6_r8 ! meters
60 : ! volume-mean emissions diameter of primary BC/OM aerosols, see :
61 : ! Liu et al, Toward a minimal representation of aerosols in climate models:
62 : ! Description and evaluation in the Community Atmosphere Model CAM5.
63 : ! Geosci. Model Dev., 5, 709–739, doi:10.5194/gmd-5-709-2012
64 : ! and Table S1 in Supplement: http://www.geosci-model-dev.net/5/709/2012/gmd-5-709-2012-supplement.pdf
65 :
66 : real(r8), parameter :: x_numfact = 1.e-6_r8 * avogad * 6.0_r8 / (pi*(demis_acc**3)) ! 1.e-6 converts m-3 to cm-3.
67 : real(r8) :: specdens ! kg/m3
68 : logical :: found
69 : character(len=12) :: units
70 :
71 0 : if (shr_fire_emis_mechcomps_n<1) return
72 :
73 0 : if (shr_fire_emis_elevated) then ! initialize elevated forcings
74 :
75 0 : allocate( frc_spc_map(shr_fire_emis_mechcomps_n) )
76 0 : allocate( chm_spc_map(shr_fire_emis_mechcomps_n) )
77 0 : allocate( frc_num_map(shr_fire_emis_mechcomps_n) )
78 0 : allocate( spc_mass_factor(shr_fire_emis_mechcomps_n) )
79 0 : allocate( num_mass_factor(shr_fire_emis_mechcomps_n) )
80 0 : allocate( fire_frc_name(shr_fire_emis_mechcomps_n) )
81 0 : allocate( fire_numfrc_name(shr_fire_emis_mechcomps_n) )
82 0 : allocate( fire_sflx_name(shr_fire_emis_mechcomps_n) )
83 0 : allocate( fire_vflx_name(shr_fire_emis_mechcomps_n) )
84 :
85 0 : frc_spc_map(:) = -1
86 0 : frc_num_map(:) = -1
87 0 : spc_mass_factor(:) = nan
88 0 : num_mass_factor(:) = nan
89 :
90 0 : call addfld ('Fire_ZTOP', horiz_only, 'A', 'm', 'top of vertical fire emissions' )
91 :
92 0 : do n=1,shr_fire_emis_mechcomps_n
93 :
94 0 : name = shr_fire_emis_mechcomps(n)%name
95 0 : fire_frc_name(n) = 'FireFrc_'//trim(name)
96 0 : fire_sflx_name(n) = 'FireSFLX_'//trim(name)
97 0 : fire_vflx_name(n) = 'FireVFLX_'//trim(name)
98 :
99 0 : call addfld (fire_frc_name(n), (/'lev'/), 'A','molecules/cm^3/s', 'vertical fire emissions for '//trim(name))
100 0 : call addfld (fire_sflx_name(n),horiz_only,'A','kg/m^2/s', 'surface fire emissions for '//trim(name))
101 0 : call addfld (fire_vflx_name(n),horiz_only,'A','kg/m^2/s', 'vertically integrated fire emissions for '//trim(name))
102 :
103 0 : frc_ndx = get_extfrc_ndx( name )
104 0 : spc_ndx = get_spc_ndx( name )
105 :
106 0 : if (frc_ndx>0 .and. spc_ndx>0) then
107 0 : frc_spc_map(n) = frc_ndx
108 0 : chm_spc_map(n) = spc_ndx
109 : else
110 0 : write(iulog,*) 'fire_emissions_init: not able to map '//trim(name)//' to chem species/forcing '
111 0 : write(iulog,*) 'fire_emissions_init: ... frc_ndx = ',frc_ndx
112 0 : write(iulog,*) 'fire_emissions_init: ... spc_ndx = ',spc_ndx
113 0 : call endrun('fire_emissions_init: not able to map '//trim(name)//' to chem species/forcing ')
114 : endif
115 :
116 0 : spc_mass_factor(n) = 1.e-6_r8 * avogad / adv_mass(spc_ndx) ! 1.e-6 converts m-3 to cm-3.
117 : ! (molecules/kmole) / (g/mole) --> molecules/kg
118 :
119 : ! for MAM need to include cooresponding forcings of number densities
120 :
121 : found = rad_cnst_num_name(0, name, num_name, mode_out=mode, spec_out=spec )
122 :
123 0 : if ( found ) then
124 :
125 : frc_ndx = get_extfrc_ndx( num_name )
126 :
127 : call rad_cnst_get_aer_props(0, mode, spec, density_aer=specdens)
128 : frc_num_map(n) = frc_ndx
129 : num_mass_factor(n) = x_numfact / specdens
130 :
131 : fire_numfrc_name(n) = 'FireFrc_'//trim(name)//'_'//trim(num_name)
132 : call addfld (fire_numfrc_name(n),(/'lev'/), 'A', 'molecules/cm^3/s', &
133 : 'vertical fire emissions for '//trim(num_name)//' due to component '//trim(name))
134 :
135 : endif
136 :
137 : enddo
138 :
139 : else ! initialize surface emissions
140 :
141 0 : allocate( fire_emis_indices_map(shr_fire_emis_mechcomps_n) )
142 :
143 0 : do n=1,shr_fire_emis_mechcomps_n
144 0 : call cnst_get_ind (shr_fire_emis_mechcomps(n)%name, fire_emis_indices_map(n), abort=.false. )
145 :
146 0 : ii = get_spc_ndx(shr_fire_emis_mechcomps(n)%name)
147 0 : if (ii<1) then
148 : call endrun('gas_phase_chemdr_inti: Fire emissions compound not in chemistry mechanism : '&
149 0 : //trim(shr_fire_emis_mechcomps(n)%name))
150 : endif
151 :
152 0 : if (index(shr_fire_emis_mechcomps(n)%name,'num_')>0) then
153 0 : units = '1/m2/sec'
154 : else
155 0 : units = 'kg/m2/sec'
156 : endif
157 :
158 : ! Fire emis history fields
159 : call addfld( 'FireSF_'//trim(shr_fire_emis_mechcomps(n)%name),horiz_only,'A',units,&
160 0 : trim(shr_fire_emis_mechcomps(n)%name)//' Fire emissions flux')
161 :
162 : enddo
163 :
164 : endif
165 :
166 : end subroutine fire_emissions_init
167 :
168 : !------------------------------------------------------------------------------
169 : ! sets surface emissions
170 : !------------------------------------------------------------------------------
171 0 : subroutine fire_emissions_srf( lchnk, ncol, fireflx, sflx )
172 :
173 : ! dummy args
174 : integer, intent(in) :: lchnk, ncol
175 : real(r8), pointer, intent(in) :: fireflx(:,:)
176 : real(r8), intent(inout) :: sflx(:,:)
177 :
178 : ! local vars
179 : integer :: i, n
180 :
181 : ! fire surface emissions if not elevated forcing
182 0 : if ((.not.shr_fire_emis_elevated) .and. active_Fall_flxfire .and. shr_fire_emis_mechcomps_n>0 ) then
183 :
184 : ! set Fire Emis fluxes ( add to other emis )
185 0 : do i =1,ncol
186 0 : do n = 1,shr_fire_emis_mechcomps_n
187 0 : sflx(i,fire_emis_indices_map(n)) &
188 0 : = sflx(i,fire_emis_indices_map(n)) + fireflx(i,n)
189 : enddo
190 : end do
191 :
192 : ! output fire emis fluxes to history
193 0 : do n = 1,shr_fire_emis_mechcomps_n
194 0 : call outfld('FireSF_'//trim(shr_fire_emis_mechcomps(n)%name), fireflx(:ncol,n), ncol, lchnk)
195 : enddo
196 :
197 : endif
198 :
199 0 : end subroutine fire_emissions_srf
200 :
201 : !------------------------------------------------------------------------------
202 : ! sets vertical emissions (forcings)
203 : ! vertically distributes wild fire emissions
204 : !------------------------------------------------------------------------------
205 0 : subroutine fire_emissions_vrt( ncol, lchnk, zint, fire_sflx, fire_ztop, frcing )
206 :
207 : ! args
208 : integer, intent(in) :: ncol,lchnk
209 : real(r8), intent(in) :: zint(:,:) ! interface geopot above surface (km)
210 : real(r8),pointer, intent(in) :: fire_sflx(:,:) ! fire surface emissions (kg/m2/sec)
211 : real(r8),pointer, intent(in) :: fire_ztop(:) ! top of vert distribution of fire surface emissions (m)
212 : real(r8), intent(inout) :: frcing(:,:,:) ! insitu forcings (molecules/cm3/sec)
213 :
214 : ! local vars
215 0 : real(r8) :: vertical_fire(ncol,pver), ztop
216 : integer :: n, i,k
217 0 : real(r8) :: fire_frc(ncol,pver)
218 0 : real(r8) :: sflx(ncol)
219 :
220 0 : if (.not.shr_fire_emis_elevated) return
221 0 : if (shr_fire_emis_mechcomps_n<1) return
222 :
223 : ! define vertical_fire from Dentener units /m
224 0 : do k=1,pver
225 0 : do i=1,ncol
226 0 : ztop = fire_ztop(i)*1.e-3_r8 ! convert m to km
227 0 : if(zint(i,k)<ztop)then
228 0 : vertical_fire(i,k)=1.e-3_r8/(ztop-zint(i,pverp))
229 0 : elseif(zint(i,k)>ztop.and.zint(i,k+1)<ztop)then
230 0 : vertical_fire(i,k)=1.e-3_r8*(ztop-zint(i,k+1))/(zint(i,k)-zint(i,k+1))/(ztop-zint(i,pverp))
231 : else
232 0 : vertical_fire(i,k)=0._r8
233 : endif
234 : end do
235 : end do
236 :
237 : ! extend the surface emissions in the vertical ....
238 :
239 : do n=1,shr_fire_emis_mechcomps_n
240 :
241 0 : fire_frc(:,:) = 0._r8
242 :
243 0 : do k = 1,pver ! (kg/m2/sec) (/m) (molecules/kg)*1.e-6 --> molecules/cm3/s
244 0 : fire_frc(:ncol,k) = fire_sflx(:ncol,n) * vertical_fire(:ncol,k) * spc_mass_factor(n) ! molecules/cm3/s
245 : enddo
246 0 : call outfld( fire_frc_name(n), fire_frc, ncol, lchnk )
247 0 : frcing(:ncol,:,frc_spc_map(n)) = frcing(:ncol,:,frc_spc_map(n)) + fire_frc(:ncol,:)
248 :
249 : ! for debugging ...
250 : ! vertical intergration of the forcing should get back the surface flux
251 0 : sflx(:) = 0._r8
252 0 : do k = 1,pver
253 0 : sflx(:ncol) = sflx(:ncol) + 1.e5_r8*(zint(:ncol,k)-zint(:ncol,k+1))*fire_frc(:ncol,k) ! molecules/cm3/s --> molecules/cm2/sec
254 : enddo
255 0 : sflx(:ncol) = sflx(:ncol) * 1.e4_r8 * adv_mass(chm_spc_map(n))/avogad ! molecules/cm2/sec --> kg/m2/sec
256 : ! / avogad --> kmoles/cm2/sec
257 : ! * adv_mass --> kg/cm2/sec
258 : ! * 1.e4 --> kg/ m2/sec
259 : call outfld( fire_vflx_name(n), sflx(:ncol ), ncol, lchnk )
260 : call outfld( fire_sflx_name(n), fire_sflx(:ncol,n), ncol, lchnk )
261 :
262 : ! for MAM need to include corresponding forcings of number densities
263 0 : if (frc_num_map(n)>0) then
264 : do k = 1,pver
265 : fire_frc(:ncol,k) = fire_sflx(:ncol,n) * vertical_fire(:ncol,k) * num_mass_factor(n) ! molecules/cm3/s
266 : enddo
267 : call outfld( fire_numfrc_name(n), fire_frc, ncol, lchnk )
268 : frcing(:ncol,:,frc_num_map(n)) = frcing(:ncol,:,frc_num_map(n)) + fire_frc(:ncol,:)
269 : endif
270 :
271 : enddo
272 :
273 : call outfld( 'Fire_ZTOP', fire_ztop(:ncol), ncol, lchnk )
274 :
275 : end subroutine fire_emissions_vrt
276 :
277 : end module fire_emissions
|
