1 | !> @file chem_emissions_mod.f90 |
---|
2 | !--------------------------------------------------------------------------------------------------! |
---|
3 | ! This file is part of PALM model system. |
---|
4 | ! |
---|
5 | ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General |
---|
6 | ! Public License as published by the Free Software Foundation, either version 3 of the License, or |
---|
7 | ! (at your option) any later version. |
---|
8 | ! |
---|
9 | ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the |
---|
10 | ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General |
---|
11 | ! Public License for more details. |
---|
12 | ! |
---|
13 | ! You should have received a copy of the GNU General Public License along with PALM. If not, see |
---|
14 | ! <http://www.gnu.org/licenses/>. |
---|
15 | ! |
---|
16 | ! Copyright 2018-2020 Leibniz Universitaet Hannover |
---|
17 | ! Copyright 2018-2020 Freie Universitaet Berlin |
---|
18 | ! Copyright 2018-2020 Karlsruhe Institute of Technology |
---|
19 | !--------------------------------------------------------------------------------------------------! |
---|
20 | ! |
---|
21 | ! Current revisions: |
---|
22 | ! ------------------ |
---|
23 | ! |
---|
24 | ! |
---|
25 | ! Former revisions: |
---|
26 | ! ----------------- |
---|
27 | ! $Id: chem_emissions_mod.f90 4559 2020-06-11 08:51:48Z suehring $ |
---|
28 | ! file re-formatted to follow the PALM coding standard |
---|
29 | ! |
---|
30 | ! 4481 2020-03-31 18:55:54Z maronga |
---|
31 | ! Implemented on-demand read mode for LOD 2 (ECC) |
---|
32 | ! - added following module global variables |
---|
33 | ! - input_file_chem (namesake in netcdf_data_input_mod is local) |
---|
34 | ! - timestamps |
---|
35 | ! - added following public subroutines / interfaces |
---|
36 | ! - chem_emisisons_header_init |
---|
37 | ! - chem_emisisons_update_on_demand |
---|
38 | ! - added following local subroutines |
---|
39 | ! - chem_emisisons_header_init_lod2 |
---|
40 | ! - chem_emisisons_update_on_demand_lod2 |
---|
41 | ! - added following local auxiliary subroutines |
---|
42 | ! - chem_emissions_init_species ( ) |
---|
43 | ! - chem_emissions_init_timestamps ( ) |
---|
44 | ! - chem_emissions_assign_surface_flux ( ) |
---|
45 | ! - added following local functions |
---|
46 | ! - chem_emisisons_convert_base_units ( ) |
---|
47 | ! - chem_emissions_mass_2_molar_flux ( ) |
---|
48 | ! - chem_emissions_locate_species ( ) |
---|
49 | ! - chem_emissions_locate_timestep ( ) |
---|
50 | ! - added following error messages |
---|
51 | ! - CM0468 - LOD mismatch (namelist / chemistry file) |
---|
52 | ! - CM0469 - Timestamps no in choronological order |
---|
53 | ! - depreciated unused module variable filename_emis |
---|
54 | ! |
---|
55 | ! 4356 2019-12-20 17:09:33Z suehring |
---|
56 | ! Minor formatting adjustment |
---|
57 | ! |
---|
58 | ! 4242 2019-09-27 12:59:10Z suehring |
---|
59 | ! Adjust index_hh access to new definition accompanied with new |
---|
60 | ! palm_date_time_mod. Note, this is just a preleminary fix. (E Chan) |
---|
61 | ! |
---|
62 | ! 4230 2019-09-11 13:58:14Z suehring |
---|
63 | ! Bugfix, consider that time_since_reference_point can be also negative when |
---|
64 | ! time indices are determined. |
---|
65 | ! |
---|
66 | ! 4227 2019-09-10 18:04:34Z gronemeier |
---|
67 | ! implement new palm_date_time_mod |
---|
68 | ! |
---|
69 | ! 4223 2019-09-10 09:20:47Z gronemeier |
---|
70 | ! Unused routine chem_emissions_check_parameters commented out due to uninitialized content |
---|
71 | ! |
---|
72 | ! 4182 2019-08-22 15:20:23Z scharf |
---|
73 | ! Corrected "Former revisions" section |
---|
74 | ! |
---|
75 | ! 4157 2019-08-14 09:19:12Z suehring |
---|
76 | ! Replace global arrays also in mode_emis branch |
---|
77 | ! |
---|
78 | ! 4154 2019-08-13 13:35:59Z suehring |
---|
79 | ! Replace global arrays for emissions by local ones. |
---|
80 | ! |
---|
81 | ! 4144 2019-08-06 09:11:47Z raasch |
---|
82 | ! relational operators .EQ., .NE., etc. replaced by ==, /=, etc. |
---|
83 | ! |
---|
84 | ! 4055 2019-06-27 09:47:29Z suehring |
---|
85 | ! - replaced local thermo. constants w/ module definitions in basic_constants_and_equations_mod |
---|
86 | ! (rgas_univ, p_0, r_d_cp) |
---|
87 | ! - initialize array emis_distribution immediately following allocation |
---|
88 | ! - lots of minor formatting changes based on review sesson in 20190325 (E.C. Chan) |
---|
89 | ! |
---|
90 | ! 3968 2019-05-13 11:04:01Z suehring |
---|
91 | ! - in subroutine chem_emissions_match replace all decision structures relating to mode_emis to |
---|
92 | ! emiss_lod |
---|
93 | ! - in subroutine chem_check_parameters replace emt%nspec with emt%n_emiss_species |
---|
94 | ! - spring cleaning (E.C. Chan) |
---|
95 | ! |
---|
96 | ! 3885 2019-04-11 11:29:34Z kanani |
---|
97 | ! Changes related to global restructuring of location messages and introduction of additional debug |
---|
98 | ! messages |
---|
99 | ! |
---|
100 | ! 3831 2019-03-28 09:11:22Z forkel |
---|
101 | ! added nvar to USE chem_gasphase_mod (chem_modules will not include nvar anymore) |
---|
102 | ! |
---|
103 | ! 3788 2019-03-07 11:40:09Z banzhafs |
---|
104 | ! Removed unused variables from chem_emissions_mod |
---|
105 | ! |
---|
106 | ! 3772 2019-02-28 15:51:57Z suehring |
---|
107 | ! - In case of parametrized emissions, assure that emissions are only on natural surfaces |
---|
108 | ! (i.e. streets) and not on urban surfaces. |
---|
109 | ! - some unnecessary if clauses removed |
---|
110 | ! |
---|
111 | ! 3685 2019 -01-21 01:02:11Z knoop |
---|
112 | ! Some interface calls moved to module_interface + cleanup |
---|
113 | ! 3286 2018-09-28 07:41:39Z forkel |
---|
114 | ! |
---|
115 | ! Authors: |
---|
116 | ! -------- |
---|
117 | ! @author Emmanuele Russo (FU-Berlin) |
---|
118 | ! @author Sabine Banzhaf (FU-Berlin) |
---|
119 | ! @author Martijn Schaap (FU-Berlin, TNO Utrecht) |
---|
120 | ! |
---|
121 | ! Description: |
---|
122 | ! ------------ |
---|
123 | !> MODULE for reading-in Chemistry Emissions data |
---|
124 | !> |
---|
125 | !> @todo Rename nspec to n_emis to avoid inteferece with nspec from chem_gasphase_mod |
---|
126 | !> @todo Check_parameters routine should be developed: for now it includes just one condition |
---|
127 | !> @todo Use of Restart files not contempled at the moment |
---|
128 | !> @todo revise indices of files read from the netcdf: preproc_emission_data and expert_emission_data |
---|
129 | !> @todo for now emission data may be passed on a singular vertical level: need to be more flexible |
---|
130 | !> @todo fill/activate restart option in chem_emissions_init |
---|
131 | !> @todo discuss dt_emis |
---|
132 | !> @note <Enter notes on the module> |
---|
133 | !> @bug <Enter known bugs here> |
---|
134 | !--------------------------------------------------------------------------------------------------! |
---|
135 | |
---|
136 | MODULE chem_emissions_mod |
---|
137 | |
---|
138 | USE arrays_3d, & |
---|
139 | ONLY: rho_air |
---|
140 | |
---|
141 | USE basic_constants_and_equations_mod, & |
---|
142 | ONLY: p_0, rd_d_cp, rgas_univ |
---|
143 | |
---|
144 | USE control_parameters, & |
---|
145 | ONLY: debug_output, end_time, initializing_actions, intermediate_timestep_count, & |
---|
146 | message_string, dt_3d |
---|
147 | |
---|
148 | USE indices |
---|
149 | |
---|
150 | USE kinds |
---|
151 | |
---|
152 | #if defined( __netcdf ) |
---|
153 | USE netcdf |
---|
154 | #endif |
---|
155 | |
---|
156 | USE netcdf_data_input_mod, & |
---|
157 | ONLY: chem_emis_att_type, chem_emis_val_type |
---|
158 | |
---|
159 | USE chem_gasphase_mod, & |
---|
160 | ONLY: nvar, spc_names |
---|
161 | |
---|
162 | USE chem_modules |
---|
163 | |
---|
164 | USE statistics, & |
---|
165 | ONLY: weight_pres |
---|
166 | |
---|
167 | ! |
---|
168 | !-- 20200203 (ECC) |
---|
169 | !-- Added new palm_date_time_mod for on-demand emission reading |
---|
170 | |
---|
171 | USE palm_date_time_mod, & |
---|
172 | ONLY: get_date_time |
---|
173 | |
---|
174 | IMPLICIT NONE |
---|
175 | |
---|
176 | ! |
---|
177 | !-- Declare all global variables within the module |
---|
178 | |
---|
179 | ! |
---|
180 | !-- 20200203 (ECC) - variable unused |
---|
181 | ! CHARACTER (LEN=80) :: filename_emis !< Variable for the name of the netcdf input file |
---|
182 | |
---|
183 | ! |
---|
184 | !-- 20200203 (ECC) new variables for on-demand read mode |
---|
185 | |
---|
186 | CHARACTER(LEN=*), PARAMETER :: input_file_chem = 'PIDS_CHEM' !< chemistry file |
---|
187 | |
---|
188 | CHARACTER(LEN=512), ALLOCATABLE, DIMENSION(:) :: timestamps !< timestamps in chemistry file |
---|
189 | |
---|
190 | |
---|
191 | INTEGER(iwp) :: dt_emis !< Time Step Emissions |
---|
192 | INTEGER(iwp) :: i !< index 1st selected dimension (some dims are not spatial) |
---|
193 | INTEGER(iwp) :: j !< index 2nd selected dimension |
---|
194 | INTEGER(iwp) :: i_end !< Index to end read variable from netcdf in one dims |
---|
195 | INTEGER(iwp) :: i_start !< Index to start read variable from netcdf along one dims |
---|
196 | INTEGER(iwp) :: j_end !< Index to end read variable from netcdf in additional dims |
---|
197 | INTEGER(iwp) :: j_start !< Index to start read variable from netcdf in additional dims |
---|
198 | INTEGER(iwp) :: len_index !< length of index (used for several indices) |
---|
199 | INTEGER(iwp) :: len_index_pm !< length of PMs index |
---|
200 | INTEGER(iwp) :: len_index_voc !< length of voc index |
---|
201 | INTEGER(iwp) :: previous_timestamp_index !< index for current timestamp (20200203 ECC) |
---|
202 | INTEGER(iwp) :: z_end !< Index to end read variable from netcdf in additional dims |
---|
203 | INTEGER(iwp) :: z_start !< Index to start read variable from netcdf in additional dims |
---|
204 | |
---|
205 | REAL(wp) :: conversion_factor !< Units Conversion Factor |
---|
206 | |
---|
207 | SAVE |
---|
208 | |
---|
209 | ! ! |
---|
210 | ! !-- Checks Input parameters |
---|
211 | ! INTERFACE chem_emissions_check_parameters |
---|
212 | ! MODULE PROCEDURE chem_emissions_check_parameters |
---|
213 | ! END INTERFACE chem_emissions_check_parameters |
---|
214 | ! |
---|
215 | !-- Matching Emissions actions |
---|
216 | INTERFACE chem_emissions_match |
---|
217 | MODULE PROCEDURE chem_emissions_match |
---|
218 | END INTERFACE chem_emissions_match |
---|
219 | ! |
---|
220 | !-- Initialization actions |
---|
221 | INTERFACE chem_emissions_init |
---|
222 | MODULE PROCEDURE chem_emissions_init |
---|
223 | END INTERFACE chem_emissions_init |
---|
224 | ! |
---|
225 | !-- Setup of Emissions |
---|
226 | INTERFACE chem_emissions_setup |
---|
227 | MODULE PROCEDURE chem_emissions_setup |
---|
228 | END INTERFACE chem_emissions_setup |
---|
229 | |
---|
230 | ! |
---|
231 | !-- 20200203 (ECC) new interfaces for on-demand mode |
---|
232 | |
---|
233 | ! |
---|
234 | !-- initialization actions for on-demand mode |
---|
235 | INTERFACE chem_emissions_header_init |
---|
236 | MODULE PROCEDURE chem_emissions_header_init |
---|
237 | END INTERFACE chem_emissions_header_init |
---|
238 | ! |
---|
239 | !-- load emission data for on-demand mode |
---|
240 | INTERFACE chem_emissions_update_on_demand |
---|
241 | MODULE PROCEDURE chem_emissions_update_on_demand |
---|
242 | END INTERFACE chem_emissions_update_on_demand |
---|
243 | |
---|
244 | ! |
---|
245 | !-- 20200203 (ECC) update public routines |
---|
246 | |
---|
247 | ! PUBLIC chem_emissions_init, chem_emissions_match, chem_emissions_setup |
---|
248 | |
---|
249 | PUBLIC chem_emissions_init, chem_emissions_match, chem_emissions_setup, & |
---|
250 | chem_emissions_header_init, chem_emissions_update_on_demand |
---|
251 | ! |
---|
252 | !-- Public Variables |
---|
253 | PUBLIC conversion_factor, len_index, len_index_pm, len_index_voc |
---|
254 | |
---|
255 | CONTAINS |
---|
256 | |
---|
257 | ! !------------------------------------------------------------------------------------------------! |
---|
258 | ! ! Description: |
---|
259 | ! ! ------------ |
---|
260 | ! !> Routine for checking input parameters |
---|
261 | ! !------------------------------------------------------------------------------------------------! |
---|
262 | ! SUBROUTINE chem_emissions_check_parameters |
---|
263 | ! |
---|
264 | ! IMPLICIT NONE |
---|
265 | ! |
---|
266 | ! TYPE(chem_emis_att_type) :: emt |
---|
267 | ! |
---|
268 | ! ! |
---|
269 | ! !-- Check if species count matches the number of names |
---|
270 | ! !-- passed for the chemiscal species |
---|
271 | ! |
---|
272 | ! IF ( SIZE(emt%species_name) /= emt%n_emiss_species ) THEN |
---|
273 | ! ! IF ( SIZE(emt%species_name) /= emt%nspec ) THEN |
---|
274 | ! |
---|
275 | ! message_string = 'Numbers of input emission species names and number of species' // & |
---|
276 | ! 'for which emission values are given do not match' |
---|
277 | ! CALL message( 'chem_emissions_check_parameters', 'CM0437', 2, 2, 0, 6, 0 ) |
---|
278 | ! |
---|
279 | ! ENDIF |
---|
280 | ! |
---|
281 | ! END SUBROUTINE chem_emissions_check_parameters |
---|
282 | |
---|
283 | |
---|
284 | !--------------------------------------------------------------------------------------------------! |
---|
285 | ! Description: |
---|
286 | ! ------------ |
---|
287 | !> Matching the chemical species indices. The routine checks what are the indices of the emission |
---|
288 | !> input species and the corresponding ones of the model species. The routine gives as output a |
---|
289 | !> vector containing the number of common species: it is important to note that while the model |
---|
290 | !> species are distinct, their values could be given to a single species in input. |
---|
291 | !> For example, in the case of NO2 and NO, values may be passed in input as NOX values. |
---|
292 | !--------------------------------------------------------------------------------------------------! |
---|
293 | |
---|
294 | SUBROUTINE chem_emissions_match( emt_att,len_index ) |
---|
295 | |
---|
296 | INTEGER(iwp) :: ind_inp !< Parameters for cycling through chemical input species |
---|
297 | INTEGER(iwp) :: ind_mod !< Parameters for cycling through chemical model species |
---|
298 | INTEGER(iwp) :: ind_voc !< Indices to check whether a split for voc should be done |
---|
299 | INTEGER(iwp) :: ispec !< index for cycle over effective number of emission species |
---|
300 | INTEGER(iwp) :: nspec_emis_inp !< Variable where to store # of emission species in input |
---|
301 | |
---|
302 | INTEGER(iwp), INTENT(INOUT) :: len_index !< number of common species between input dataset & model species |
---|
303 | |
---|
304 | TYPE(chem_emis_att_type), INTENT(INOUT) :: emt_att !< Chemistry Emission Array (decl. netcdf_data_input.f90) |
---|
305 | |
---|
306 | |
---|
307 | IF ( debug_output ) CALL debug_message( 'chem_emissions_match', 'start' ) |
---|
308 | |
---|
309 | ! |
---|
310 | !-- Number of input emission species |
---|
311 | |
---|
312 | nspec_emis_inp = emt_att%n_emiss_species |
---|
313 | ! nspec_emis_inp=emt_att%nspec |
---|
314 | |
---|
315 | ! |
---|
316 | !-- Check the emission LOD: 0 (PARAMETERIZED), 1 (DEFAULT), 2 (PRE-PROCESSED) |
---|
317 | SELECT CASE (emiss_lod) |
---|
318 | |
---|
319 | ! |
---|
320 | !-- LOD 0 (PARAMETERIZED mode) |
---|
321 | CASE (0) |
---|
322 | |
---|
323 | len_index = 0 |
---|
324 | ! |
---|
325 | !-- Number of species and number of matched species can be different but call is only made if |
---|
326 | !-- both are greater than zero. |
---|
327 | IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN |
---|
328 | |
---|
329 | ! |
---|
330 | !-- Cycle over model species |
---|
331 | DO ind_mod = 1, nvar |
---|
332 | ind_inp = 1 |
---|
333 | DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) !< 'novalue' is the default |
---|
334 | IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN |
---|
335 | len_index = len_index + 1 |
---|
336 | ENDIF |
---|
337 | ind_inp = ind_inp + 1 |
---|
338 | ENDDO |
---|
339 | ENDDO |
---|
340 | |
---|
341 | IF ( len_index > 0 ) THEN |
---|
342 | |
---|
343 | ! |
---|
344 | !-- Allocation of Arrays of the matched species |
---|
345 | ALLOCATE ( match_spec_input(len_index) ) |
---|
346 | ALLOCATE ( match_spec_model(len_index) ) |
---|
347 | |
---|
348 | ! |
---|
349 | !-- Pass species indices to declared arrays |
---|
350 | len_index = 0 |
---|
351 | |
---|
352 | DO ind_mod = 1, nvar |
---|
353 | ind_inp = 1 |
---|
354 | DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) |
---|
355 | IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM(spc_names(ind_mod) ) ) & |
---|
356 | THEN |
---|
357 | len_index = len_index + 1 |
---|
358 | match_spec_input(len_index) = ind_inp |
---|
359 | match_spec_model(len_index) = ind_mod |
---|
360 | ENDIF |
---|
361 | ind_inp = ind_inp + 1 |
---|
362 | END DO |
---|
363 | END DO |
---|
364 | |
---|
365 | ! |
---|
366 | !-- Check |
---|
367 | DO ispec = 1, len_index |
---|
368 | |
---|
369 | IF ( emiss_factor_main(match_spec_input(ispec) ) < 0 .AND. & |
---|
370 | emiss_factor_side(match_spec_input(ispec) ) < 0 ) THEN |
---|
371 | |
---|
372 | message_string = 'PARAMETERIZED emissions mode selected:' // & |
---|
373 | ' EMISSIONS POSSIBLE ONLY ON STREET SURFACES' // & |
---|
374 | ' but values of scaling factors for street types' // & |
---|
375 | ' emiss_factor_main AND emiss_factor_side' // & |
---|
376 | ' not provided for each of the emissions species' // & |
---|
377 | ' or not provided at all: PLEASE set a finite value' // & |
---|
378 | ' for these parameters in the chemistry namelist' |
---|
379 | CALL message( 'chem_emissions_matching', 'CM0442', 2, 2, 0, 6, 0 ) |
---|
380 | |
---|
381 | ENDIF |
---|
382 | |
---|
383 | END DO |
---|
384 | |
---|
385 | |
---|
386 | ELSE |
---|
387 | |
---|
388 | message_string = 'Non of given Emission Species' // & |
---|
389 | ' matches' // & |
---|
390 | ' model chemical species' // & |
---|
391 | ' Emission routine is not called' |
---|
392 | CALL message( 'chem_emissions_matching', 'CM0443', 0, 0, 0, 6, 0 ) |
---|
393 | |
---|
394 | ENDIF |
---|
395 | |
---|
396 | ELSE |
---|
397 | |
---|
398 | message_string = 'Array of Emission species not allocated: ' // & |
---|
399 | ' Either no emission species are provided as input or' // & |
---|
400 | ' no chemical species are used by PALM.' // & |
---|
401 | ' Emission routine is not called' |
---|
402 | CALL message( 'chem_emissions_matching', 'CM0444', 0, 2, 0, 6, 0 ) |
---|
403 | |
---|
404 | ENDIF |
---|
405 | |
---|
406 | ! |
---|
407 | !-- LOD 1 (DEFAULT mode) |
---|
408 | CASE (1) |
---|
409 | |
---|
410 | len_index = 0 ! total number of species (to be accumulated) |
---|
411 | len_index_voc = 0 ! total number of VOCs (to be accumulated) |
---|
412 | len_index_pm = 3 ! total number of PMs: PM1, PM2.5, PM10. |
---|
413 | |
---|
414 | ! |
---|
415 | !-- Number of model species and input species could be different but process this only when both are |
---|
416 | !-- non-zero |
---|
417 | IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN |
---|
418 | |
---|
419 | ! |
---|
420 | !-- Cycle over model species |
---|
421 | DO ind_mod = 1, nvar |
---|
422 | |
---|
423 | ! |
---|
424 | !-- Cycle over input species |
---|
425 | DO ind_inp = 1, nspec_emis_inp |
---|
426 | |
---|
427 | ! |
---|
428 | !-- Check for VOC Species |
---|
429 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN |
---|
430 | DO ind_voc= 1, emt_att%nvoc |
---|
431 | |
---|
432 | IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) & |
---|
433 | THEN |
---|
434 | len_index = len_index + 1 |
---|
435 | len_index_voc = len_index_voc + 1 |
---|
436 | ENDIF |
---|
437 | |
---|
438 | END DO |
---|
439 | ENDIF |
---|
440 | |
---|
441 | ! |
---|
442 | !-- PMs: There is one input species name for all PM. This variable has 3 dimensions, one for PM1, |
---|
443 | !-- PM2.5 and PM10 |
---|
444 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" ) THEN |
---|
445 | |
---|
446 | IF ( TRIM( spc_names(ind_mod) ) == "PM1" ) THEN |
---|
447 | len_index = len_index + 1 |
---|
448 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" ) THEN |
---|
449 | len_index = len_index + 1 |
---|
450 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" ) THEN |
---|
451 | len_index = len_index + 1 |
---|
452 | ENDIF |
---|
453 | |
---|
454 | ENDIF |
---|
455 | |
---|
456 | ! |
---|
457 | !-- NOX: NO2 and NO |
---|
458 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOX" ) THEN |
---|
459 | |
---|
460 | IF ( TRIM( spc_names(ind_mod) ) == "NO" ) THEN |
---|
461 | len_index = len_index + 1 |
---|
462 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" ) THEN |
---|
463 | len_index = len_index + 1 |
---|
464 | ENDIF |
---|
465 | |
---|
466 | ENDIF |
---|
467 | |
---|
468 | ! |
---|
469 | !-- SOX: SO2 and SO4 |
---|
470 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOX" ) THEN |
---|
471 | |
---|
472 | IF ( TRIM( spc_names(ind_mod) ) == "SO2" ) THEN |
---|
473 | len_index = len_index + 1 |
---|
474 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" ) THEN |
---|
475 | len_index = len_index + 1 |
---|
476 | ENDIF |
---|
477 | |
---|
478 | ENDIF |
---|
479 | |
---|
480 | ! |
---|
481 | !-- Other Species |
---|
482 | IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN |
---|
483 | len_index = len_index + 1 |
---|
484 | ENDIF |
---|
485 | |
---|
486 | END DO ! ind_inp ... |
---|
487 | |
---|
488 | END DO ! ind_mod ... |
---|
489 | |
---|
490 | |
---|
491 | ! |
---|
492 | !-- Allocate arrays |
---|
493 | IF ( len_index > 0 ) THEN |
---|
494 | |
---|
495 | ALLOCATE ( match_spec_input(len_index) ) |
---|
496 | ALLOCATE ( match_spec_model(len_index) ) |
---|
497 | |
---|
498 | IF ( len_index_voc > 0 ) THEN |
---|
499 | ! |
---|
500 | !-- Contains indices of the VOC model species |
---|
501 | ALLOCATE( match_spec_voc_model(len_index_voc) ) |
---|
502 | ! |
---|
503 | !-- Contains the indices of different values of VOC composition of input variable |
---|
504 | !-- VOC_composition |
---|
505 | ALLOCATE( match_spec_voc_input(len_index_voc) ) |
---|
506 | |
---|
507 | ENDIF |
---|
508 | |
---|
509 | ! |
---|
510 | !-- Pass the species indices to declared arrays |
---|
511 | len_index = 0 |
---|
512 | len_index_voc = 0 |
---|
513 | |
---|
514 | DO ind_mod = 1, nvar |
---|
515 | DO ind_inp = 1, nspec_emis_inp |
---|
516 | |
---|
517 | ! |
---|
518 | !-- VOCs |
---|
519 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND. & |
---|
520 | ALLOCATED( match_spec_voc_input ) ) THEN |
---|
521 | |
---|
522 | DO ind_voc = 1, emt_att%nvoc |
---|
523 | |
---|
524 | IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) )& |
---|
525 | THEN |
---|
526 | |
---|
527 | len_index = len_index + 1 |
---|
528 | len_index_voc = len_index_voc + 1 |
---|
529 | |
---|
530 | match_spec_input(len_index) = ind_inp |
---|
531 | match_spec_model(len_index) = ind_mod |
---|
532 | |
---|
533 | match_spec_voc_input(len_index_voc) = ind_voc |
---|
534 | match_spec_voc_model(len_index_voc) = ind_mod |
---|
535 | |
---|
536 | ENDIF |
---|
537 | |
---|
538 | END DO |
---|
539 | |
---|
540 | ENDIF |
---|
541 | |
---|
542 | ! |
---|
543 | !-- PMs |
---|
544 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" ) THEN |
---|
545 | |
---|
546 | IF ( TRIM( spc_names(ind_mod) ) == "PM1" ) THEN |
---|
547 | len_index = len_index + 1 |
---|
548 | match_spec_input(len_index) = ind_inp |
---|
549 | match_spec_model(len_index) = ind_mod |
---|
550 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" ) THEN |
---|
551 | len_index = len_index + 1 |
---|
552 | match_spec_input(len_index) = ind_inp |
---|
553 | match_spec_model(len_index) = ind_mod |
---|
554 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" ) THEN |
---|
555 | len_index = len_index + 1 |
---|
556 | match_spec_input(len_index) = ind_inp |
---|
557 | match_spec_model(len_index) = ind_mod |
---|
558 | ENDIF |
---|
559 | |
---|
560 | ENDIF |
---|
561 | |
---|
562 | ! |
---|
563 | !-- NOX |
---|
564 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOX" ) THEN |
---|
565 | |
---|
566 | IF ( TRIM( spc_names(ind_mod) ) == "NO" ) THEN |
---|
567 | len_index = len_index + 1 |
---|
568 | |
---|
569 | match_spec_input(len_index) = ind_inp |
---|
570 | match_spec_model(len_index) = ind_mod |
---|
571 | |
---|
572 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" ) THEN |
---|
573 | len_index = len_index + 1 |
---|
574 | |
---|
575 | match_spec_input(len_index) = ind_inp |
---|
576 | match_spec_model(len_index) = ind_mod |
---|
577 | |
---|
578 | ENDIF |
---|
579 | |
---|
580 | ENDIF |
---|
581 | |
---|
582 | |
---|
583 | ! |
---|
584 | !-- SOX |
---|
585 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOX" ) THEN |
---|
586 | |
---|
587 | IF ( TRIM( spc_names(ind_mod) ) == "SO2" ) THEN |
---|
588 | len_index = len_index + 1 |
---|
589 | match_spec_input(len_index) = ind_inp |
---|
590 | match_spec_model(len_index) = ind_mod |
---|
591 | ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" ) THEN |
---|
592 | len_index = len_index + 1 |
---|
593 | match_spec_input(len_index) = ind_inp |
---|
594 | match_spec_model(len_index) = ind_mod |
---|
595 | ENDIF |
---|
596 | |
---|
597 | ENDIF |
---|
598 | |
---|
599 | ! |
---|
600 | !-- Other Species |
---|
601 | IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) & |
---|
602 | THEN |
---|
603 | len_index = len_index + 1 |
---|
604 | match_spec_input(len_index) = ind_inp |
---|
605 | match_spec_model(len_index) = ind_mod |
---|
606 | ENDIF |
---|
607 | |
---|
608 | END DO ! inp_ind |
---|
609 | |
---|
610 | END DO ! inp_mod |
---|
611 | |
---|
612 | ! |
---|
613 | !-- Error reporting (no matching) |
---|
614 | ELSE |
---|
615 | |
---|
616 | message_string = 'None of given Emission Species matches' // & |
---|
617 | ' model chemical species' // & |
---|
618 | ' Emission routine is not called' |
---|
619 | CALL message( 'chem_emissions_matching', 'CM0440', 0, 0, 0, 6, 0 ) |
---|
620 | |
---|
621 | ENDIF |
---|
622 | |
---|
623 | ! |
---|
624 | !-- Error reporting (no species) |
---|
625 | ELSE |
---|
626 | |
---|
627 | message_string = 'Array of Emission species not allocated: ' // & |
---|
628 | ' Either no emission species are provided as input or' // & |
---|
629 | ' no chemical species are used by PALM:' // & |
---|
630 | ' Emission routine is not called' |
---|
631 | CALL message( 'chem_emissions_matching', 'CM0441', 0, 2, 0, 6, 0 ) |
---|
632 | |
---|
633 | ENDIF |
---|
634 | |
---|
635 | ! |
---|
636 | !-- LOD 2 (PRE-PROCESSED mode) |
---|
637 | CASE (2) |
---|
638 | |
---|
639 | len_index = 0 |
---|
640 | len_index_voc = 0 |
---|
641 | |
---|
642 | IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN |
---|
643 | ! |
---|
644 | !-- Cycle over model species |
---|
645 | DO ind_mod = 1, nvar |
---|
646 | |
---|
647 | ! |
---|
648 | !-- Cycle over input species |
---|
649 | DO ind_inp = 1, nspec_emis_inp |
---|
650 | |
---|
651 | ! |
---|
652 | !-- Check for VOC Species |
---|
653 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN |
---|
654 | DO ind_voc = 1, emt_att%nvoc |
---|
655 | IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) & |
---|
656 | THEN |
---|
657 | len_index = len_index + 1 |
---|
658 | len_index_voc = len_index_voc + 1 |
---|
659 | ENDIF |
---|
660 | END DO |
---|
661 | ENDIF |
---|
662 | |
---|
663 | ! |
---|
664 | !-- Other Species |
---|
665 | IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN |
---|
666 | len_index = len_index + 1 |
---|
667 | ENDIF |
---|
668 | ENDDO |
---|
669 | ENDDO |
---|
670 | |
---|
671 | ! |
---|
672 | !-- Allocate array for storing the indices of the matched species |
---|
673 | IF ( len_index > 0 ) THEN |
---|
674 | |
---|
675 | ALLOCATE ( match_spec_input(len_index) ) |
---|
676 | |
---|
677 | ALLOCATE ( match_spec_model(len_index) ) |
---|
678 | |
---|
679 | IF ( len_index_voc > 0 ) THEN |
---|
680 | ! |
---|
681 | !-- Contains indices of the VOC model species |
---|
682 | ALLOCATE( match_spec_voc_model(len_index_voc) ) |
---|
683 | ! |
---|
684 | !-- Contains the indices of different values of VOC composition of input variable |
---|
685 | !-- VOC_composition |
---|
686 | ALLOCATE( match_spec_voc_input(len_index_voc) ) |
---|
687 | |
---|
688 | ENDIF |
---|
689 | |
---|
690 | ! |
---|
691 | !-- Pass the species indices to declared arrays |
---|
692 | len_index = 0 |
---|
693 | |
---|
694 | ! |
---|
695 | !-- Cycle over model species |
---|
696 | DO ind_mod = 1, nvar |
---|
697 | |
---|
698 | ! |
---|
699 | !-- Cycle over Input species |
---|
700 | DO ind_inp = 1, nspec_emis_inp |
---|
701 | |
---|
702 | ! |
---|
703 | !-- VOCs |
---|
704 | IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND. & |
---|
705 | ALLOCATED( match_spec_voc_input ) ) THEN |
---|
706 | |
---|
707 | DO ind_voc= 1, emt_att%nvoc |
---|
708 | IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) )& |
---|
709 | THEN |
---|
710 | len_index = len_index + 1 |
---|
711 | len_index_voc = len_index_voc + 1 |
---|
712 | |
---|
713 | match_spec_input(len_index) = ind_inp |
---|
714 | match_spec_model(len_index) = ind_mod |
---|
715 | |
---|
716 | match_spec_voc_input(len_index_voc) = ind_voc |
---|
717 | match_spec_voc_model(len_index_voc) = ind_mod |
---|
718 | ENDIF |
---|
719 | END DO |
---|
720 | ENDIF |
---|
721 | |
---|
722 | ! |
---|
723 | !-- Other Species |
---|
724 | IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) & |
---|
725 | THEN |
---|
726 | len_index = len_index + 1 |
---|
727 | match_spec_input(len_index) = ind_inp |
---|
728 | match_spec_model(len_index) = ind_mod |
---|
729 | ENDIF |
---|
730 | |
---|
731 | END DO ! ind_inp |
---|
732 | END DO ! ind_mod |
---|
733 | |
---|
734 | ELSE ! if len_index_voc <= 0 |
---|
735 | |
---|
736 | ! |
---|
737 | !-- In case there are no species matching (just informational message) |
---|
738 | message_string = 'Non of given emission species' // & |
---|
739 | ' matches' // & |
---|
740 | ' model chemical species:' // & |
---|
741 | ' Emission routine is not called' |
---|
742 | CALL message( 'chem_emissions_matching', 'CM0438', 0, 0, 0, 6, 0 ) |
---|
743 | ENDIF |
---|
744 | |
---|
745 | ! |
---|
746 | !-- Error check (no matching) |
---|
747 | ELSE |
---|
748 | |
---|
749 | ! |
---|
750 | !-- Either spc_names is zero or nspec_emis_inp is not allocated |
---|
751 | message_string = 'Array of Emission species not allocated:' // & |
---|
752 | ' Either no emission species are provided as input or' // & |
---|
753 | ' no chemical species are used by PALM:' // & |
---|
754 | ' Emission routine is not called' |
---|
755 | CALL message( 'chem_emissions_matching', 'CM0439', 0, 2, 0, 6, 0 ) |
---|
756 | |
---|
757 | ENDIF |
---|
758 | |
---|
759 | ! |
---|
760 | !-- If emission module is switched on but mode_emis is not specified or it is given the wrong name |
---|
761 | |
---|
762 | ! |
---|
763 | !-- Error check (no species) |
---|
764 | CASE DEFAULT |
---|
765 | |
---|
766 | message_string = 'Emission Module switched ON, but' // & |
---|
767 | ' either no emission mode specified or incorrectly given :' // & |
---|
768 | ' please, pass the correct value to the namelist parameter "mode_emis"' |
---|
769 | CALL message( 'chem_emissions_matching', 'CM0445', 2, 2, 0, 6, 0 ) |
---|
770 | |
---|
771 | END SELECT |
---|
772 | |
---|
773 | IF ( debug_output ) CALL debug_message( 'chem_emissions_match', 'end' ) |
---|
774 | |
---|
775 | END SUBROUTINE chem_emissions_match |
---|
776 | |
---|
777 | |
---|
778 | !--------------------------------------------------------------------------------------------------! |
---|
779 | ! Description: |
---|
780 | ! ------------ |
---|
781 | !> Initialization: |
---|
782 | !> Netcdf reading, arrays allocation and first calculation of cssws fluxes at timestep 0 |
---|
783 | !--------------------------------------------------------------------------------------------------! |
---|
784 | |
---|
785 | SUBROUTINE chem_emissions_init |
---|
786 | |
---|
787 | USE netcdf_data_input_mod, & |
---|
788 | ONLY: chem_emis, chem_emis_att |
---|
789 | |
---|
790 | IMPLICIT NONE |
---|
791 | |
---|
792 | INTEGER(iwp) :: ispec !< running index |
---|
793 | |
---|
794 | ! |
---|
795 | !-- Actions for initial runs |
---|
796 | ! IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN |
---|
797 | !-- ... |
---|
798 | ! |
---|
799 | ! |
---|
800 | !-- Actions for restart runs |
---|
801 | ! ELSE |
---|
802 | !-- ... |
---|
803 | ! |
---|
804 | ! ENDIF |
---|
805 | |
---|
806 | |
---|
807 | IF ( debug_output ) CALL debug_message( 'chem_emissions_init', 'start' ) |
---|
808 | |
---|
809 | ! |
---|
810 | !-- Matching |
---|
811 | CALL chem_emissions_match( chem_emis_att, n_matched_vars ) |
---|
812 | |
---|
813 | IF ( n_matched_vars == 0 ) THEN |
---|
814 | |
---|
815 | emission_output_required = .FALSE. |
---|
816 | |
---|
817 | ELSE |
---|
818 | |
---|
819 | emission_output_required = .TRUE. |
---|
820 | |
---|
821 | |
---|
822 | ! |
---|
823 | !-- Set molecule masses (in kg/mol) |
---|
824 | ALLOCATE( chem_emis_att%xm(n_matched_vars) ) |
---|
825 | |
---|
826 | DO ispec = 1, n_matched_vars |
---|
827 | SELECT CASE ( TRIM( spc_names(match_spec_model(ispec)) ) ) |
---|
828 | CASE ( 'SO2' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 2 |
---|
829 | CASE ( 'SO4' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 4 |
---|
830 | CASE ( 'NO' ); chem_emis_att%xm(ispec) = xm_N + xm_O |
---|
831 | CASE ( 'NO2' ); chem_emis_att%xm(ispec) = xm_N + xm_O * 2 |
---|
832 | CASE ( 'NH3' ); chem_emis_att%xm(ispec) = xm_N + xm_H * 3 |
---|
833 | CASE ( 'CO' ); chem_emis_att%xm(ispec) = xm_C + xm_O |
---|
834 | CASE ( 'CO2' ); chem_emis_att%xm(ispec) = xm_C + xm_O * 2 |
---|
835 | CASE ( 'CH4' ); chem_emis_att%xm(ispec) = xm_C + xm_H * 4 |
---|
836 | CASE ( 'HNO3' ); chem_emis_att%xm(ispec) = xm_H + xm_N + xm_O*3 |
---|
837 | CASE DEFAULT |
---|
838 | chem_emis_att%xm(ispec) = 1.0_wp |
---|
839 | END SELECT |
---|
840 | ENDDO |
---|
841 | |
---|
842 | |
---|
843 | ! |
---|
844 | !-- Get emissions for the first time step base on LOD (if defined) or emission mode |
---|
845 | !-- (if no LOD defined) |
---|
846 | |
---|
847 | ! |
---|
848 | !-- NOTE - I could use a combined if ( lod = xxx .or. mode = 'XXX' ) type of decision structure but |
---|
849 | ! I think it is much better to implement it this way (i.e., conditional on lod if it is |
---|
850 | ! defined, and mode if not) as we can easily take out the case structure for mode_emis |
---|
851 | ! later on. |
---|
852 | |
---|
853 | IF ( emiss_lod < 0 ) THEN !-- no LOD defined (not likely) |
---|
854 | |
---|
855 | SELECT CASE ( TRIM( mode_emis ) ) |
---|
856 | |
---|
857 | CASE ( 'PARAMETERIZED' ) ! LOD 0 |
---|
858 | |
---|
859 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
860 | ALLOCATE( emis_distribution(1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
861 | ENDIF |
---|
862 | |
---|
863 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars) |
---|
864 | |
---|
865 | CASE ( 'DEFAULT' ) ! LOD 1 |
---|
866 | |
---|
867 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
868 | ALLOCATE( emis_distribution(1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
869 | ENDIF |
---|
870 | |
---|
871 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars ) |
---|
872 | |
---|
873 | CASE ( 'PRE-PROCESSED' ) ! LOD 2 |
---|
874 | |
---|
875 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
876 | ! |
---|
877 | !-- Note, at the moment emissions are considered only by surface fluxes rather than |
---|
878 | !-- by volume sources. Therefore, no vertical dimension is required and is thus |
---|
879 | !-- allocated with 1. Later when volume sources are considered, the vertical |
---|
880 | !-- dimension will increase. |
---|
881 | !ALLOCATE( emis_distribution(nzb:nzt+1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
882 | ALLOCATE( emis_distribution(1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
883 | ENDIF |
---|
884 | |
---|
885 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars ) |
---|
886 | |
---|
887 | END SELECT |
---|
888 | |
---|
889 | ELSE ! if LOD is defined |
---|
890 | |
---|
891 | SELECT CASE ( emiss_lod ) |
---|
892 | |
---|
893 | CASE ( 0 ) ! parameterized mode |
---|
894 | |
---|
895 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
896 | ALLOCATE( emis_distribution(1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
897 | ENDIF |
---|
898 | |
---|
899 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars) |
---|
900 | |
---|
901 | CASE ( 1 ) ! default mode |
---|
902 | |
---|
903 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
904 | ALLOCATE( emis_distribution(1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
905 | ENDIF |
---|
906 | |
---|
907 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars ) |
---|
908 | |
---|
909 | CASE ( 2 ) ! pre-processed mode |
---|
910 | |
---|
911 | IF ( .NOT. ALLOCATED( emis_distribution) ) THEN |
---|
912 | ALLOCATE( emis_distribution(nzb:nzt+1,nys:nyn,nxl:nxr,n_matched_vars) ) |
---|
913 | ENDIF |
---|
914 | |
---|
915 | CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars ) |
---|
916 | |
---|
917 | END SELECT |
---|
918 | |
---|
919 | ENDIF |
---|
920 | |
---|
921 | ! |
---|
922 | ! -- Initialize |
---|
923 | emis_distribution = 0.0_wp |
---|
924 | |
---|
925 | ENDIF |
---|
926 | |
---|
927 | IF ( debug_output ) CALL debug_message( 'chem_emissions_init', 'end' ) |
---|
928 | |
---|
929 | END SUBROUTINE chem_emissions_init |
---|
930 | |
---|
931 | |
---|
932 | |
---|
933 | !--------------------------------------------------------------------------------------------------! |
---|
934 | ! Description: |
---|
935 | ! ------------ |
---|
936 | !> Routine for Update of Emission values at each timestep. |
---|
937 | !> |
---|
938 | !> @todo Clarify the correct usage of index_dd, index_hh and index_mm. Consider renaming of these |
---|
939 | !> variables. |
---|
940 | !> @todo Clarify time used in emis_lod=2 mode. ATM, the used time seems strange. |
---|
941 | !--------------------------------------------------------------------------------------------------! |
---|
942 | |
---|
943 | SUBROUTINE chem_emissions_setup( emt_att, emt, n_matched_vars ) |
---|
944 | |
---|
945 | USE surface_mod, & |
---|
946 | ONLY: surf_def_h, surf_lsm_h, surf_usm_h |
---|
947 | |
---|
948 | USE netcdf_data_input_mod, & |
---|
949 | ONLY: street_type_f |
---|
950 | |
---|
951 | USE arrays_3d, & |
---|
952 | ONLY: hyp, pt |
---|
953 | |
---|
954 | USE control_parameters, & |
---|
955 | ONLY: time_since_reference_point |
---|
956 | |
---|
957 | USE palm_date_time_mod, & |
---|
958 | ONLY: days_per_week, get_date_time, hours_per_day, months_per_year, seconds_per_day |
---|
959 | |
---|
960 | IMPLICIT NONE |
---|
961 | |
---|
962 | INTEGER(iwp) :: day_of_month !< day of the month |
---|
963 | INTEGER(iwp) :: day_of_week !< day of the week |
---|
964 | INTEGER(iwp) :: day_of_year !< day of the year |
---|
965 | INTEGER(iwp) :: days_since_reference_point !< days since reference point |
---|
966 | INTEGER(iwp) :: i !< running index for grid in x-direction |
---|
967 | INTEGER(iwp) :: i_pm_comp !< index for number of PM components |
---|
968 | INTEGER(iwp) :: icat !< Index for number of categories |
---|
969 | INTEGER(iwp) :: index_dd !< index day |
---|
970 | INTEGER(iwp) :: index_hh !< index hour |
---|
971 | INTEGER(iwp) :: index_mm !< index month |
---|
972 | INTEGER(iwp) :: ispec !< index for number of species |
---|
973 | INTEGER(iwp) :: ivoc !< Index for number of VOCs |
---|
974 | INTEGER(iwp) :: hour_of_day !< hour of the day |
---|
975 | INTEGER(iwp) :: j !< running index for grid in y-direction |
---|
976 | INTEGER(iwp) :: k !< running index for grid in z-direction |
---|
977 | INTEGER(iwp) :: m !< running index for horizontal surfaces |
---|
978 | INTEGER(iwp) :: month_of_year !< month of the year |
---|
979 | |
---|
980 | INTEGER,INTENT(IN) :: n_matched_vars !< Output of matching routine with number |
---|
981 | !< of matched species |
---|
982 | |
---|
983 | REAL(wp) :: time_utc_init !< second of day of initial date |
---|
984 | |
---|
985 | TYPE(chem_emis_att_type), INTENT(INOUT) :: emt_att !< variable to store emission information |
---|
986 | |
---|
987 | TYPE(chem_emis_val_type), INTENT(INOUT), ALLOCATABLE, DIMENSION(:) :: emt !< variable to store emission input values, |
---|
988 | !< depending on the considered species |
---|
989 | ! |
---|
990 | !-- CONVERSION FACTORS: TIME |
---|
991 | REAL(wp), PARAMETER :: hour_per_year = 8760.0_wp !< number of hours in a year of 365 days |
---|
992 | REAL(wp), PARAMETER :: s_per_hour = 3600.0_wp !< number of sec per hour (s)/(hour) |
---|
993 | REAL(wp), PARAMETER :: s_per_day = 86400.0_wp !< number of sec per day (s)/(day) |
---|
994 | |
---|
995 | REAL(wp), PARAMETER :: day_to_s = 1.0_wp/s_per_day !< conversion day -> sec |
---|
996 | REAL(wp), PARAMETER :: hour_to_s = 1.0_wp/s_per_hour !< conversion hours -> sec |
---|
997 | REAL(wp), PARAMETER :: year_to_s = 1.0_wp/(s_per_hour*hour_per_year) !< conversion year -> sec |
---|
998 | |
---|
999 | |
---|
1000 | ! |
---|
1001 | !-- CONVERSION FACTORS: MASS |
---|
1002 | REAL(wp), PARAMETER :: g_to_kg = 1.0E-03_wp !< Conversion from g to kg (kg/g) |
---|
1003 | REAL(wp), PARAMETER :: miug_to_kg = 1.0E-09_wp !< Conversion from g to kg (kg/g) |
---|
1004 | REAL(wp), PARAMETER :: tons_to_kg = 100.0_wp !< Conversion from tons to kg (kg/tons) |
---|
1005 | ! |
---|
1006 | !-- CONVERSION FACTORS: PPM |
---|
1007 | REAL(wp), PARAMETER :: ratio2ppm = 1.0E+06_wp |
---|
1008 | |
---|
1009 | REAL(wp), DIMENSION(24) :: par_emis_time_factor !< time factors for the parameterized mode: |
---|
1010 | !< fixed houlry profile for example day |
---|
1011 | REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: conv_to_ratio !< factor used for converting input |
---|
1012 | !< to concentration ratio |
---|
1013 | REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: tmp_temp !< temporary variable for abs. temperature |
---|
1014 | |
---|
1015 | REAL(wp), DIMENSION(:), ALLOCATABLE :: time_factor !< factor for time scaling of emissions |
---|
1016 | |
---|
1017 | REAL(wp), DIMENSION(:,:), ALLOCATABLE :: delta_emis !< incremental emission factor |
---|
1018 | REAL(wp), DIMENSION(:,:), ALLOCATABLE :: emis !< emission factor |
---|
1019 | |
---|
1020 | |
---|
1021 | IF ( emission_output_required ) THEN |
---|
1022 | |
---|
1023 | ! |
---|
1024 | !-- Set emis_dt to be used - since chemistry ODEs can be stiff, the option to solve them at every |
---|
1025 | !-- RK substep is present to help improve stability should the need arises |
---|
1026 | |
---|
1027 | IF ( call_chem_at_all_substeps ) THEN |
---|
1028 | |
---|
1029 | dt_emis = dt_3d * weight_pres(intermediate_timestep_count) |
---|
1030 | |
---|
1031 | ELSE |
---|
1032 | |
---|
1033 | dt_emis = dt_3d |
---|
1034 | |
---|
1035 | ENDIF |
---|
1036 | |
---|
1037 | ! |
---|
1038 | !-- Conversion of units to the ones employed in PALM |
---|
1039 | !-- In PARAMETERIZED mode no conversion is performed: in this case input units are fixed |
---|
1040 | IF ( TRIM( mode_emis ) == "DEFAULT" .OR. TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN |
---|
1041 | |
---|
1042 | SELECT CASE ( TRIM( emt_att%units ) ) |
---|
1043 | |
---|
1044 | CASE ( 'kg/m2/s', 'KG/M2/S' ); conversion_factor = 1.0_wp ! kg |
---|
1045 | CASE ( 'kg/m2/hour', 'KG/M2/HOUR' ); conversion_factor = hour_to_s |
---|
1046 | CASE ( 'kg/m2/day', 'KG/M2/DAY' ); conversion_factor = day_to_s |
---|
1047 | CASE ( 'kg/m2/year', 'KG/M2/YEAR' ); conversion_factor = year_to_s |
---|
1048 | |
---|
1049 | CASE ( 'ton/m2/s', 'TON/M2/S' ); conversion_factor = tons_to_kg ! tonnes |
---|
1050 | CASE ( 'ton/m2/hour', 'TON/M2/HOUR' ); conversion_factor = tons_to_kg*hour_to_s |
---|
1051 | CASE ( 'ton/m2/year', 'TON/M2/YEAR' ); conversion_factor = tons_to_kg*year_to_s |
---|
1052 | |
---|
1053 | CASE ( 'g/m2/s', 'G/M2/S' ); conversion_factor = g_to_kg ! grams |
---|
1054 | CASE ( 'g/m2/hour', 'G/M2/HOUR' ); conversion_factor = g_to_kg*hour_to_s |
---|
1055 | CASE ( 'g/m2/year', 'G/M2/YEAR' ); conversion_factor = g_to_kg*year_to_s |
---|
1056 | |
---|
1057 | CASE ( 'micrograms/m2/s', 'MICROGRAMS/M2/S' ); conversion_factor = miug_to_kg ! ug |
---|
1058 | CASE ( 'micrograms/m2/hour', 'MICROGRAMS/M2/HOUR' ); conversion_factor = miug_to_kg*hour_to_s |
---|
1059 | CASE ( 'micrograms/m2/year', 'MICROGRAMS/M2/YEAR' ); conversion_factor = miug_to_kg*year_to_s |
---|
1060 | |
---|
1061 | ! |
---|
1062 | !-- Error check (need units) |
---|
1063 | CASE DEFAULT |
---|
1064 | message_string = 'The Units of the provided emission input' // & |
---|
1065 | ' are not the ones required by PALM-4U: please check ' // & |
---|
1066 | ' emission module documentation.' |
---|
1067 | CALL message( 'chem_emissions_setup', 'CM0446', 2, 2, 0, 6, 0 ) |
---|
1068 | |
---|
1069 | END SELECT |
---|
1070 | |
---|
1071 | ENDIF |
---|
1072 | |
---|
1073 | ! |
---|
1074 | !-- Conversion factor to convert kg/m**2/s to ppm/s |
---|
1075 | DO i = nxl, nxr |
---|
1076 | DO j = nys, nyn |
---|
1077 | |
---|
1078 | ! |
---|
1079 | !-- Derive Temperature from Potential Temperature |
---|
1080 | tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * ( hyp(nzb:nzt+1) / p_0 )**rd_d_cp |
---|
1081 | |
---|
1082 | ! |
---|
1083 | !-- We need to pass to cssws <- (ppm/s) * dz |
---|
1084 | !-- Input is Nmole/(m^2*s) |
---|
1085 | !-- To go to ppm*dz multiply the input by (m**2/N)*dz |
---|
1086 | !-- (m**2/N)*dz == V/N |
---|
1087 | !-- V/N = RT/P |
---|
1088 | conv_to_ratio(nzb:nzt+1,j,i) = rgas_univ * & ! J K-1 mol-1 |
---|
1089 | tmp_temp(nzb:nzt+1,j,i) / & ! K |
---|
1090 | hyp(nzb:nzt+1) ! Pa |
---|
1091 | |
---|
1092 | ! (ecc) for reference |
---|
1093 | ! m**3/Nmole (J/mol)*K^-1 K Pa |
---|
1094 | ! conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) ) |
---|
1095 | |
---|
1096 | ENDDO |
---|
1097 | ENDDO |
---|
1098 | |
---|
1099 | |
---|
1100 | ! (ecc) moved initialization immediately after allocation |
---|
1101 | ! |
---|
1102 | !-- Initialize |
---|
1103 | |
---|
1104 | ! emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp |
---|
1105 | |
---|
1106 | |
---|
1107 | ! |
---|
1108 | !-- LOD 2 (PRE-PROCESSED MODE) |
---|
1109 | |
---|
1110 | IF ( emiss_lod == 2 ) THEN |
---|
1111 | |
---|
1112 | ! for reference (ecc) |
---|
1113 | ! IF ( TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN |
---|
1114 | |
---|
1115 | ! |
---|
1116 | !-- Update time indices |
---|
1117 | CALL get_date_time( 0.0_wp, second_of_day=time_utc_init ) |
---|
1118 | CALL get_date_time( MAX( 0.0_wp, time_since_reference_point ), hour=hour_of_day ) |
---|
1119 | |
---|
1120 | days_since_reference_point = INT( FLOOR( ( time_utc_init + & |
---|
1121 | MAX( 0.0_wp, time_since_reference_point ) ) & |
---|
1122 | / seconds_per_day ) ) |
---|
1123 | |
---|
1124 | index_hh = days_since_reference_point * hours_per_day + hour_of_day |
---|
1125 | |
---|
1126 | ! |
---|
1127 | !-- LOD 1 (DEFAULT MODE) |
---|
1128 | ELSEIF ( emiss_lod == 1 ) THEN |
---|
1129 | |
---|
1130 | ! for reference (ecc) |
---|
1131 | ! ELSEIF ( TRIM( mode_emis ) == "DEFAULT" ) THEN |
---|
1132 | |
---|
1133 | ! |
---|
1134 | !-- Allocate array where to store temporary emission values |
---|
1135 | IF ( .NOT. ALLOCATED(emis) ) ALLOCATE( emis(nys:nyn,nxl:nxr) ) |
---|
1136 | |
---|
1137 | ! |
---|
1138 | !-- Allocate time factor per category |
---|
1139 | ALLOCATE( time_factor(emt_att%ncat) ) |
---|
1140 | |
---|
1141 | ! |
---|
1142 | !-- Read-in hourly emission time factor |
---|
1143 | IF ( TRIM( time_fac_type ) == "HOUR" ) THEN |
---|
1144 | |
---|
1145 | ! |
---|
1146 | !-- Update time indices |
---|
1147 | CALL get_date_time( MAX( time_since_reference_point, 0.0_wp ), & |
---|
1148 | day_of_year=day_of_year, hour=hour_of_day ) |
---|
1149 | index_hh = ( day_of_year - 1_iwp ) * hour_of_day |
---|
1150 | |
---|
1151 | ! |
---|
1152 | !-- Check if the index is less or equal to the temporal dimension of HOURLY emission files |
---|
1153 | IF ( index_hh <= SIZE( emt_att%hourly_emis_time_factor(1,:) ) ) THEN |
---|
1154 | |
---|
1155 | ! |
---|
1156 | !-- Read-in the correspondant time factor |
---|
1157 | time_factor(:) = emt_att%hourly_emis_time_factor(:,index_hh+1) |
---|
1158 | |
---|
1159 | ! |
---|
1160 | !-- Error check (time out of range) |
---|
1161 | ELSE |
---|
1162 | |
---|
1163 | message_string = 'The "HOUR" time-factors in the DEFAULT mode ' // & |
---|
1164 | ' are not provided for each hour of the total simulation time' |
---|
1165 | CALL message( 'chem_emissions_setup', 'CM0448', 2, 2, 0, 6, 0 ) |
---|
1166 | |
---|
1167 | ENDIF |
---|
1168 | |
---|
1169 | ! |
---|
1170 | !-- Read-in MDH emissions time factors |
---|
1171 | ELSEIF ( TRIM( time_fac_type ) == "MDH" ) THEN |
---|
1172 | |
---|
1173 | ! |
---|
1174 | !-- Update time indices |
---|
1175 | CALL get_date_time( MAX( time_since_reference_point, 0.0_wp ), & |
---|
1176 | month = month_of_year, & |
---|
1177 | day = day_of_month, & |
---|
1178 | hour = hour_of_day, & |
---|
1179 | day_of_week = day_of_week & |
---|
1180 | ) |
---|
1181 | index_mm = month_of_year |
---|
1182 | index_dd = months_per_year + day_of_week |
---|
1183 | SELECT CASE( TRIM( daytype_mdh ) ) |
---|
1184 | |
---|
1185 | CASE ("workday") |
---|
1186 | index_hh = months_per_year + days_per_week + hour_of_day |
---|
1187 | |
---|
1188 | CASE ("weekend") |
---|
1189 | index_hh = months_per_year + days_per_week + hours_per_day + hour_of_day |
---|
1190 | |
---|
1191 | CASE ("holiday") |
---|
1192 | index_hh = months_per_year + days_per_week + 2*hours_per_day + hour_of_day |
---|
1193 | |
---|
1194 | END SELECT |
---|
1195 | |
---|
1196 | ! |
---|
1197 | !-- Check if the index is less or equal to the temporal dimension of MDH emission files |
---|
1198 | IF ( ( index_hh + index_dd + index_mm) <= SIZE( emt_att%mdh_emis_time_factor(1,:) ) )& |
---|
1199 | THEN |
---|
1200 | ! |
---|
1201 | !-- Read in corresponding time factor |
---|
1202 | time_factor(:) = emt_att%mdh_emis_time_factor(:,index_mm) * & |
---|
1203 | emt_att%mdh_emis_time_factor(:,index_dd) * & |
---|
1204 | emt_att%mdh_emis_time_factor(:,index_hh+1) |
---|
1205 | |
---|
1206 | ! |
---|
1207 | !-- Error check (MDH time factor not provided) |
---|
1208 | ELSE |
---|
1209 | |
---|
1210 | message_string = 'The "MDH" time-factors in the DEFAULT mode ' // & |
---|
1211 | ' are not provided for each hour/day/month of the total simulation time' |
---|
1212 | CALL message( 'chem_emissions_setup', 'CM0449', 2, 2, 0, 6, 0 ) |
---|
1213 | |
---|
1214 | ENDIF |
---|
1215 | |
---|
1216 | ! |
---|
1217 | !-- Error check (no time factor defined) |
---|
1218 | ELSE |
---|
1219 | |
---|
1220 | message_string = 'In the DEFAULT mode the time factor' // & |
---|
1221 | ' has to be defined in the NAMELIST' |
---|
1222 | CALL message( 'chem_emissions_setup', 'CM0450', 2, 2, 0, 6, 0 ) |
---|
1223 | |
---|
1224 | ENDIF |
---|
1225 | |
---|
1226 | ! |
---|
1227 | !-- PARAMETERIZED MODE |
---|
1228 | ELSEIF ( emiss_lod == 0 ) THEN |
---|
1229 | |
---|
1230 | |
---|
1231 | ! for reference (ecc) |
---|
1232 | ! ELSEIF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN |
---|
1233 | |
---|
1234 | ! |
---|
1235 | !-- Assign constant values of time factors, diurnal profile for traffic sector |
---|
1236 | par_emis_time_factor(:) = (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, & |
---|
1237 | 0.056, 0.053, 0.051, 0.051, 0.052, 0.055, & |
---|
1238 | 0.059, 0.061, 0.064, 0.067, 0.069, 0.069, & |
---|
1239 | 0.049, 0.039, 0.039, 0.029, 0.024, 0.019 /) |
---|
1240 | |
---|
1241 | IF ( .NOT. ALLOCATED( time_factor ) ) ALLOCATE( time_factor(1) ) |
---|
1242 | |
---|
1243 | ! |
---|
1244 | !-- Get time-factor for specific hour |
---|
1245 | CALL get_date_time( MAX( time_since_reference_point, 0.0_wp ), hour = hour_of_day ) |
---|
1246 | |
---|
1247 | index_hh = hour_of_day |
---|
1248 | time_factor(1) = par_emis_time_factor(index_hh+1) |
---|
1249 | |
---|
1250 | ENDIF ! emiss_lod |
---|
1251 | |
---|
1252 | |
---|
1253 | ! |
---|
1254 | !-- Emission distribution calculation |
---|
1255 | |
---|
1256 | ! |
---|
1257 | !-- LOD 0 (PARAMETERIZED mode) |
---|
1258 | IF ( emiss_lod == 0 ) THEN |
---|
1259 | |
---|
1260 | ! for reference (ecc) |
---|
1261 | ! IF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN |
---|
1262 | |
---|
1263 | DO ispec = 1, n_matched_vars |
---|
1264 | |
---|
1265 | ! |
---|
1266 | !-- Units are micromoles/m**2*day (or kilograms/m**2*day for PMs) |
---|
1267 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = surface_csflux( match_spec_input(ispec) )& |
---|
1268 | * time_factor(1) * hour_to_s |
---|
1269 | |
---|
1270 | ENDDO |
---|
1271 | |
---|
1272 | |
---|
1273 | ! |
---|
1274 | !-- LOD 1 (DEFAULT mode) |
---|
1275 | ELSEIF ( emiss_lod == 1 ) THEN |
---|
1276 | |
---|
1277 | ! for referene (ecc) |
---|
1278 | ! ELSEIF ( TRIM( mode_emis ) == "DEFAULT" ) THEN |
---|
1279 | |
---|
1280 | ! |
---|
1281 | !-- Allocate array for the emission value corresponding to a specific category and time factor |
---|
1282 | ALLOCATE (delta_emis(nys:nyn,nxl:nxr)) |
---|
1283 | |
---|
1284 | ! |
---|
1285 | !-- Cycle over categories |
---|
1286 | DO icat = 1, emt_att%ncat |
---|
1287 | |
---|
1288 | ! |
---|
1289 | !-- Cycle over Species: n_matched_vars represents the number of species in common between |
---|
1290 | !-- the emission input data and the chemistry mechanism used |
---|
1291 | DO ispec = 1, n_matched_vars |
---|
1292 | |
---|
1293 | emis(nys:nyn,nxl:nxr) = emt( match_spec_input(ispec) )% & |
---|
1294 | default_emission_data(icat,nys+1:nyn+1,nxl+1:nxr+1) |
---|
1295 | |
---|
1296 | ! |
---|
1297 | !-- NO |
---|
1298 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "NO" ) THEN |
---|
1299 | |
---|
1300 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1301 | time_factor(icat) * & |
---|
1302 | emt_att%nox_comp(icat,1) * & |
---|
1303 | conversion_factor * hours_per_day |
---|
1304 | |
---|
1305 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1306 | emis_distribution(1,nys:nyn,nxl:nxr,ispec)& |
---|
1307 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1308 | ! |
---|
1309 | !-- NO2 |
---|
1310 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "NO2" ) THEN |
---|
1311 | |
---|
1312 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1313 | time_factor(icat) * & |
---|
1314 | emt_att%nox_comp(icat,2) * & |
---|
1315 | conversion_factor * hours_per_day |
---|
1316 | |
---|
1317 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1318 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1319 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1320 | |
---|
1321 | ! |
---|
1322 | !-- SO2 |
---|
1323 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "SO2" ) THEN |
---|
1324 | |
---|
1325 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1326 | time_factor(icat) * & |
---|
1327 | emt_att%sox_comp(icat,1) * & |
---|
1328 | conversion_factor * hours_per_day |
---|
1329 | |
---|
1330 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1331 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1332 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1333 | |
---|
1334 | ! |
---|
1335 | !-- SO4 |
---|
1336 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "SO4" ) THEN |
---|
1337 | |
---|
1338 | |
---|
1339 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1340 | time_factor(icat) * & |
---|
1341 | emt_att%sox_comp(icat,2) * & |
---|
1342 | conversion_factor * hours_per_day |
---|
1343 | |
---|
1344 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1345 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1346 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1347 | |
---|
1348 | |
---|
1349 | ! |
---|
1350 | !-- PM1 |
---|
1351 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" ) THEN |
---|
1352 | |
---|
1353 | DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,1) ) ! cycle through components |
---|
1354 | |
---|
1355 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1356 | time_factor(icat) * & |
---|
1357 | emt_att%pm_comp(icat,i_pm_comp,1) * & |
---|
1358 | conversion_factor * hours_per_day |
---|
1359 | |
---|
1360 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1361 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1362 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1363 | |
---|
1364 | ENDDO |
---|
1365 | |
---|
1366 | ! |
---|
1367 | !-- PM2.5 |
---|
1368 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" ) THEN |
---|
1369 | |
---|
1370 | DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,2) ) ! cycle through components |
---|
1371 | |
---|
1372 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1373 | time_factor(icat) * & |
---|
1374 | emt_att%pm_comp(icat,i_pm_comp,2) * & |
---|
1375 | conversion_factor * hours_per_day |
---|
1376 | |
---|
1377 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1378 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1379 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1380 | |
---|
1381 | ENDDO |
---|
1382 | |
---|
1383 | ! |
---|
1384 | !-- PM10 |
---|
1385 | ELSEIF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1386 | |
---|
1387 | DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,3) ) ! cycle through components |
---|
1388 | |
---|
1389 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s |
---|
1390 | time_factor(icat) * & |
---|
1391 | emt_att%pm_comp(icat,i_pm_comp,3) * & |
---|
1392 | conversion_factor * hours_per_day |
---|
1393 | |
---|
1394 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1395 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1396 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1397 | |
---|
1398 | ENDDO |
---|
1399 | |
---|
1400 | ! |
---|
1401 | !-- VOCs |
---|
1402 | ELSEIF ( SIZE( match_spec_voc_input ) > 0 ) THEN |
---|
1403 | |
---|
1404 | DO ivoc = 1, SIZE( match_spec_voc_input ) ! cycle through components |
---|
1405 | |
---|
1406 | IF ( TRIM( spc_names(match_spec_model(ispec) ) ) == & |
---|
1407 | TRIM( emt_att%voc_name(ivoc) ) ) THEN |
---|
1408 | |
---|
1409 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & |
---|
1410 | time_factor(icat) * & |
---|
1411 | emt_att%voc_comp(icat,match_spec_voc_input(ivoc)) * & |
---|
1412 | conversion_factor * hours_per_day |
---|
1413 | |
---|
1414 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1415 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1416 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1417 | |
---|
1418 | ENDIF |
---|
1419 | |
---|
1420 | ENDDO |
---|
1421 | |
---|
1422 | ! |
---|
1423 | !-- Any other species |
---|
1424 | ELSE |
---|
1425 | |
---|
1426 | delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * time_factor(icat) * & |
---|
1427 | conversion_factor * hours_per_day |
---|
1428 | |
---|
1429 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1430 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) & |
---|
1431 | + delta_emis(nys:nyn,nxl:nxr) |
---|
1432 | |
---|
1433 | ENDIF ! TRIM spc_names |
---|
1434 | |
---|
1435 | emis = 0 |
---|
1436 | |
---|
1437 | ENDDO |
---|
1438 | |
---|
1439 | delta_emis = 0 |
---|
1440 | |
---|
1441 | ENDDO |
---|
1442 | |
---|
1443 | ! |
---|
1444 | !-- LOD 2 (PRE-PROCESSED mode) |
---|
1445 | ELSEIF ( emiss_lod == 2 ) THEN |
---|
1446 | |
---|
1447 | ! for reference (ecc) |
---|
1448 | ! ELSEIF ( TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN |
---|
1449 | |
---|
1450 | ! |
---|
1451 | !-- Cycle over species: n_matched_vars represents the number of species in common between the |
---|
1452 | !-- emission input data and the chemistry mechanism used |
---|
1453 | DO ispec = 1, n_matched_vars |
---|
1454 | |
---|
1455 | ! (ecc) |
---|
1456 | emis_distribution(1,nys:nyn,nxl:nxr,ispec) = emt(match_spec_input(ispec))% & |
---|
1457 | preproc_emission_data(index_hh+1,1,nys+1:nyn+1,nxl+1:nxr+1) & |
---|
1458 | * conversion_factor |
---|
1459 | |
---|
1460 | |
---|
1461 | ! emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & |
---|
1462 | ! emt(match_spec_input(ispec))% & |
---|
1463 | ! preproc_emission_data(index_hh,1,:,:) * & |
---|
1464 | ! conversion_factor |
---|
1465 | ENDDO |
---|
1466 | |
---|
1467 | ENDIF ! emiss_lod |
---|
1468 | |
---|
1469 | |
---|
1470 | ! |
---|
1471 | !-- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to |
---|
1472 | !-- cssws |
---|
1473 | |
---|
1474 | ! |
---|
1475 | !-- LOD 0 (PARAMETERIZED mode) |
---|
1476 | !-- Units of inputs are micromoles/m2/s |
---|
1477 | IF ( emiss_lod == 0 ) THEN |
---|
1478 | ! for reference (ecc) |
---|
1479 | ! IF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN |
---|
1480 | |
---|
1481 | IF (street_type_f%from_file) THEN |
---|
1482 | |
---|
1483 | ! |
---|
1484 | !-- Streets are lsm surfaces, hence, no usm surface treatment required. |
---|
1485 | !-- However, urban surface may be initialized via default initialization in surface_mod, e.g. at |
---|
1486 | !-- horizontal urban walls that are at k == 0 (building is lower than the first grid point). Hence, |
---|
1487 | !-- in order to have only emissions at streets, set the surfaces emissions to zero at urban walls. |
---|
1488 | IF ( surf_usm_h%ns >=1 ) surf_usm_h%cssws = 0.0_wp |
---|
1489 | |
---|
1490 | ! |
---|
1491 | !-- Treat land-surfaces. |
---|
1492 | DO m = 1, surf_lsm_h%ns |
---|
1493 | |
---|
1494 | i = surf_lsm_h%i(m) |
---|
1495 | j = surf_lsm_h%j(m) |
---|
1496 | k = surf_lsm_h%k(m) |
---|
1497 | |
---|
1498 | ! |
---|
1499 | !-- Set everything to zero then reassign according to street type |
---|
1500 | surf_lsm_h%cssws(:,m) = 0.0_wp |
---|
1501 | |
---|
1502 | IF ( street_type_f%var(j,i) >= main_street_id .AND. & |
---|
1503 | street_type_f%var(j,i) < max_street_id ) THEN |
---|
1504 | |
---|
1505 | ! |
---|
1506 | !-- Cycle over matched species |
---|
1507 | DO ispec = 1, n_matched_vars |
---|
1508 | |
---|
1509 | ! |
---|
1510 | !-- PMs are already in kilograms |
---|
1511 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" .OR. & |
---|
1512 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" .OR. & |
---|
1513 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1514 | |
---|
1515 | ! |
---|
1516 | !-- kg/(m^2*s) * kg/m^3 |
---|
1517 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & |
---|
1518 | emiss_factor_main(match_spec_input(ispec)) * & |
---|
1519 | emis_distribution(1,j,i,ispec) * & ! kg/(m^2*s) |
---|
1520 | rho_air(k) ! kg/m^3 |
---|
1521 | |
---|
1522 | ! |
---|
1523 | !-- Other Species |
---|
1524 | !-- Inputs are micromoles |
---|
1525 | ELSE |
---|
1526 | |
---|
1527 | ! |
---|
1528 | !-- ppm/s *m *kg/m^3 |
---|
1529 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & |
---|
1530 | emiss_factor_main( match_spec_input(ispec) ) * & |
---|
1531 | emis_distribution(1,j,i,ispec) * & ! micromoles/(m^2*s) |
---|
1532 | conv_to_ratio(k,j,i) * & ! m^3/Nmole |
---|
1533 | rho_air(k) ! kg/m^3 |
---|
1534 | |
---|
1535 | ENDIF |
---|
1536 | |
---|
1537 | ENDDO ! ispec |
---|
1538 | |
---|
1539 | |
---|
1540 | ELSEIF ( street_type_f%var(j,i) >= side_street_id .AND. & |
---|
1541 | street_type_f%var(j,i) < main_street_id ) THEN |
---|
1542 | |
---|
1543 | ! |
---|
1544 | !-- Cycle over matched species |
---|
1545 | DO ispec = 1, n_matched_vars |
---|
1546 | |
---|
1547 | ! |
---|
1548 | !-- PMs are already in kilograms |
---|
1549 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" .OR. & |
---|
1550 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" .OR. & |
---|
1551 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1552 | |
---|
1553 | ! |
---|
1554 | !-- kg/(m^2*s) * kg/m^3 |
---|
1555 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & |
---|
1556 | emiss_factor_side( match_spec_input(ispec) ) * & |
---|
1557 | emis_distribution(1,j,i,ispec) * & ! kg/(m^2*s) |
---|
1558 | rho_air(k) ! kg/m^3 |
---|
1559 | ! |
---|
1560 | !-- Other species |
---|
1561 | !-- Inputs are micromoles |
---|
1562 | ELSE |
---|
1563 | |
---|
1564 | ! |
---|
1565 | !-- ppm/s *m *kg/m^3 |
---|
1566 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & |
---|
1567 | emiss_factor_side( match_spec_input(ispec) ) * & |
---|
1568 | emis_distribution(1,j,i,ispec) * & ! micromoles/(m^2*s) |
---|
1569 | conv_to_ratio(k,j,i) * & ! m^3/Nmole |
---|
1570 | rho_air(k) ! kg/m^3 |
---|
1571 | |
---|
1572 | ENDIF |
---|
1573 | |
---|
1574 | ENDDO ! ispec |
---|
1575 | |
---|
1576 | ! |
---|
1577 | !-- If no street type is defined, then assign zero emission to all the species |
---|
1578 | ! (ecc) moved to front (for reference) |
---|
1579 | ! ELSE |
---|
1580 | ! |
---|
1581 | ! surf_lsm_h%cssws(:,m) = 0.0_wp |
---|
1582 | |
---|
1583 | ENDIF ! street type |
---|
1584 | |
---|
1585 | ENDDO ! m |
---|
1586 | |
---|
1587 | ENDIF ! street_type_f%from_file |
---|
1588 | |
---|
1589 | |
---|
1590 | ! |
---|
1591 | !-- LOD 1 (DEFAULT) and LOD 2 (PRE-PROCESSED) |
---|
1592 | ELSE |
---|
1593 | |
---|
1594 | |
---|
1595 | DO ispec = 1, n_matched_vars |
---|
1596 | |
---|
1597 | ! |
---|
1598 | !-- Default surfaces |
---|
1599 | DO m = 1, surf_def_h(0)%ns |
---|
1600 | |
---|
1601 | i = surf_def_h(0)%i(m) |
---|
1602 | j = surf_def_h(0)%j(m) |
---|
1603 | |
---|
1604 | IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN |
---|
1605 | |
---|
1606 | ! |
---|
1607 | !-- PMs |
---|
1608 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" .OR. & |
---|
1609 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" .OR. & |
---|
1610 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1611 | |
---|
1612 | surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3 |
---|
1613 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1614 | rho_air(nzb) ! kg/m^3 |
---|
1615 | |
---|
1616 | ELSE |
---|
1617 | |
---|
1618 | ! |
---|
1619 | !-- VOCs |
---|
1620 | IF ( len_index_voc > 0 .AND. & |
---|
1621 | emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN |
---|
1622 | |
---|
1623 | surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1624 | emis_distribution(1,j,i,ispec) * & ! mole/m2/s |
---|
1625 | conv_to_ratio(nzb,j,i) * & ! m^3/mole |
---|
1626 | ratio2ppm * & ! ppm |
---|
1627 | rho_air(nzb) ! kg/m^3 |
---|
1628 | |
---|
1629 | |
---|
1630 | ! |
---|
1631 | !-- Other species |
---|
1632 | ELSE |
---|
1633 | |
---|
1634 | surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1635 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1636 | ( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg |
---|
1637 | conv_to_ratio(nzb,j,i) * & ! m^3/mole |
---|
1638 | ratio2ppm * & ! ppm |
---|
1639 | rho_air(nzb) ! kg/m^3 |
---|
1640 | |
---|
1641 | ENDIF ! VOC |
---|
1642 | |
---|
1643 | ENDIF ! PM |
---|
1644 | |
---|
1645 | ENDIF ! emis_distribution > 0 |
---|
1646 | |
---|
1647 | ENDDO ! m |
---|
1648 | |
---|
1649 | ! |
---|
1650 | !-- LSM surfaces |
---|
1651 | DO m = 1, surf_lsm_h%ns |
---|
1652 | |
---|
1653 | i = surf_lsm_h%i(m) |
---|
1654 | j = surf_lsm_h%j(m) |
---|
1655 | k = surf_lsm_h%k(m) |
---|
1656 | |
---|
1657 | IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN |
---|
1658 | |
---|
1659 | ! |
---|
1660 | !-- PMs |
---|
1661 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" .OR. & |
---|
1662 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" .OR. & |
---|
1663 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1664 | |
---|
1665 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3 |
---|
1666 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1667 | rho_air(k) ! kg/m^3 |
---|
1668 | |
---|
1669 | ELSE |
---|
1670 | |
---|
1671 | ! |
---|
1672 | !-- VOCs |
---|
1673 | IF ( len_index_voc > 0 .AND. & |
---|
1674 | emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN |
---|
1675 | |
---|
1676 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1677 | emis_distribution(1,j,i,ispec) * & ! mole/m2/s |
---|
1678 | conv_to_ratio(k,j,i) * & ! m^3/mole |
---|
1679 | ratio2ppm * & ! ppm |
---|
1680 | rho_air(k) ! kg/m^3 |
---|
1681 | |
---|
1682 | ! |
---|
1683 | !-- Other species |
---|
1684 | ELSE |
---|
1685 | |
---|
1686 | surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1687 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1688 | ( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg |
---|
1689 | conv_to_ratio(k,j,i) * & ! m^3/mole |
---|
1690 | ratio2ppm * & ! ppm |
---|
1691 | rho_air(k) ! kg/m^3 |
---|
1692 | |
---|
1693 | ENDIF ! VOC |
---|
1694 | |
---|
1695 | ENDIF ! PM |
---|
1696 | |
---|
1697 | ENDIF ! emis_distribution |
---|
1698 | |
---|
1699 | ENDDO ! m |
---|
1700 | |
---|
1701 | |
---|
1702 | |
---|
1703 | ! |
---|
1704 | !-- USM surfaces |
---|
1705 | DO m = 1, surf_usm_h%ns |
---|
1706 | |
---|
1707 | i = surf_usm_h%i(m) |
---|
1708 | j = surf_usm_h%j(m) |
---|
1709 | k = surf_usm_h%k(m) |
---|
1710 | |
---|
1711 | IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN |
---|
1712 | |
---|
1713 | ! |
---|
1714 | !-- PMs |
---|
1715 | IF ( TRIM( spc_names( match_spec_model(ispec) ) ) == "PM1" .OR. & |
---|
1716 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM25" .OR. & |
---|
1717 | TRIM( spc_names( match_spec_model(ispec) ) ) == "PM10" ) THEN |
---|
1718 | |
---|
1719 | surf_usm_h%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3 |
---|
1720 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1721 | rho_air(k) ! kg/m^3 |
---|
1722 | |
---|
1723 | ELSE |
---|
1724 | |
---|
1725 | ! |
---|
1726 | !-- VOCs |
---|
1727 | IF ( len_index_voc > 0 .AND. & |
---|
1728 | emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN |
---|
1729 | |
---|
1730 | surf_usm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1731 | emis_distribution(1,j,i,ispec) * & ! m2/s |
---|
1732 | conv_to_ratio(k,j,i) * & ! m^3/mole |
---|
1733 | ratio2ppm * & ! ppm |
---|
1734 | rho_air(k) ! kg/m^3 |
---|
1735 | |
---|
1736 | ! |
---|
1737 | !-- Other species |
---|
1738 | ELSE |
---|
1739 | |
---|
1740 | surf_usm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3 |
---|
1741 | emis_distribution(1,j,i,ispec) * & ! kg/m2/s |
---|
1742 | ( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg |
---|
1743 | conv_to_ratio(k,j,i) * & ! m^3/mole |
---|
1744 | ratio2ppm * & ! ppm |
---|
1745 | rho_air(k) ! kg/m^3 |
---|
1746 | |
---|
1747 | |
---|
1748 | ENDIF ! VOC |
---|
1749 | |
---|
1750 | ENDIF ! PM |
---|
1751 | |
---|
1752 | ENDIF ! emis_distribution |
---|
1753 | |
---|
1754 | ENDDO ! m |
---|
1755 | |
---|
1756 | ENDDO |
---|
1757 | |
---|
1758 | ENDIF |
---|
1759 | |
---|
1760 | ! |
---|
1761 | !-- Deallocate time_factor in case of DEFAULT mode) |
---|
1762 | IF ( ALLOCATED( time_factor ) ) DEALLOCATE( time_factor ) |
---|
1763 | |
---|
1764 | ENDIF |
---|
1765 | |
---|
1766 | END SUBROUTINE chem_emissions_setup |
---|
1767 | |
---|
1768 | |
---|
1769 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1770 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1771 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1772 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1773 | !! |
---|
1774 | !! 20200203 (ECC) - ON DEMAND EMISSION UPDATE MODE |
---|
1775 | !! |
---|
1776 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1777 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1778 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1779 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1780 | |
---|
1781 | |
---|
1782 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1783 | !! |
---|
1784 | !! WRAPPER / INTERFACE FUNCTIONS |
---|
1785 | !! |
---|
1786 | !! NOTE - I find using an explicity wrapper provides much better flow control |
---|
1787 | !! over an interface block |
---|
1788 | !! |
---|
1789 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1790 | |
---|
1791 | ! |
---|
1792 | !-- 20200203 (ECC) |
---|
1793 | ! |
---|
1794 | !--------------------------------------------------------------------------------------------------! |
---|
1795 | ! Description: |
---|
1796 | ! ------------ |
---|
1797 | !> interface for initiation of emission arrays based on emission LOD |
---|
1798 | ! |
---|
1799 | !--------------------------------------------------------------------------------------------------! |
---|
1800 | |
---|
1801 | SUBROUTINE chem_emissions_header_init |
---|
1802 | |
---|
1803 | IMPLICIT NONE |
---|
1804 | |
---|
1805 | SELECT CASE ( emiss_lod ) |
---|
1806 | CASE ( 0 ) |
---|
1807 | ! do nothing at the moment |
---|
1808 | CASE ( 1 ) |
---|
1809 | ! do nothing at the moment |
---|
1810 | CASE ( 2 ) |
---|
1811 | CALL chem_emissions_header_init_lod2 |
---|
1812 | END SELECT |
---|
1813 | |
---|
1814 | END SUBROUTINE chem_emissions_header_init |
---|
1815 | |
---|
1816 | |
---|
1817 | ! |
---|
1818 | !-- 20200203 (ECC) |
---|
1819 | ! |
---|
1820 | !--------------------------------------------------------------------------------------------------! |
---|
1821 | ! Description: |
---|
1822 | ! ------------ |
---|
1823 | !> interface for initiation of emission arrays based on emission LOD |
---|
1824 | !--------------------------------------------------------------------------------------------------! |
---|
1825 | |
---|
1826 | SUBROUTINE chem_emissions_update_on_demand |
---|
1827 | |
---|
1828 | IMPLICIT NONE |
---|
1829 | |
---|
1830 | SELECT CASE ( emiss_lod ) |
---|
1831 | CASE ( 0 ) |
---|
1832 | ! do nothing at the moment |
---|
1833 | CASE ( 1 ) |
---|
1834 | ! do nothing at the moment |
---|
1835 | CASE ( 2 ) |
---|
1836 | CALL chem_emissions_update_on_demand_lod2 |
---|
1837 | END SELECT |
---|
1838 | |
---|
1839 | END SUBROUTINE ! chem_emisisons_update_on_demand |
---|
1840 | |
---|
1841 | |
---|
1842 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1843 | !! |
---|
1844 | !! SUBROUTINES SPECIFIC FOR LOD 2 |
---|
1845 | !! |
---|
1846 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
1847 | |
---|
1848 | ! |
---|
1849 | !-- 20200203 (ECC) |
---|
1850 | ! |
---|
1851 | !--------------------------------------------------------------------------------------------------! |
---|
1852 | ! Description: |
---|
1853 | ! ------------ |
---|
1854 | !> Initiates header for emissions data attributes for LOD 2 |
---|
1855 | !--------------------------------------------------------------------------------------------------! |
---|
1856 | |
---|
1857 | SUBROUTINE chem_emissions_header_init_lod2 |
---|
1858 | |
---|
1859 | USE control_parameters, & |
---|
1860 | ONLY: coupling_char, message_string |
---|
1861 | |
---|
1862 | USE netcdf_data_input_mod, & |
---|
1863 | ONLY: chem_emis_att, close_input_file, get_attribute, get_dimension_length, get_variable, & |
---|
1864 | open_read_file |
---|
1865 | |
---|
1866 | |
---|
1867 | IMPLICIT NONE |
---|
1868 | |
---|
1869 | |
---|
1870 | INTEGER(iwp) :: att_lod !< lod attribute in chemistry file |
---|
1871 | INTEGER(iwp) :: ncid !< chemistry file netCDF handle |
---|
1872 | |
---|
1873 | IF ( debug_output ) CALL debug_message( 'chem_emissions_header_init_lod2', 'start' ) |
---|
1874 | |
---|
1875 | ! |
---|
1876 | !-- Opens _chemistry input file and obtain header information |
---|
1877 | CALL open_read_file ( TRIM( input_file_chem ) // TRIM( coupling_char ), ncid ) |
---|
1878 | ! |
---|
1879 | !-- Check if LOD in chemistry file matches LOD in namelist |
---|
1880 | CALL get_attribute ( ncid, 'lod', att_lod, .TRUE. ) |
---|
1881 | |
---|
1882 | IF ( att_lod /= emiss_lod ) THEN |
---|
1883 | message_string = '' ! to get around unused variable warning / error |
---|
1884 | WRITE ( message_string, * ) 'LOD mismatch between namelist (emiss_lod) and', & |
---|
1885 | CHAR( 10 ), ' ', 'chemistry input file (global attributes>lod)' |
---|
1886 | CALL message( 'chem_emissions_header_init_lod2', 'CM0468', 1, 2, 0, 6, 0 ) |
---|
1887 | ENDIF |
---|
1888 | ! |
---|
1889 | !-- Obtain unit conversion factor |
---|
1890 | CALL get_attribute ( ncid, 'units', chem_emis_att%units, .FALSE., "emission_values" ) |
---|
1891 | conversion_factor = chem_emissions_convert_base_units ( chem_emis_att%units ) |
---|
1892 | ! |
---|
1893 | !-- Obtain header attributes |
---|
1894 | CALL chem_emissions_init_species ( ncid ) |
---|
1895 | CALL chem_emissions_init_timestamps ( ncid ) |
---|
1896 | ! |
---|
1897 | !-- Done reading file |
---|
1898 | CALL close_input_file (ncid) |
---|
1899 | |
---|
1900 | ! |
---|
1901 | !-- Set previous timestamp index to something different to trigger a read event later on |
---|
1902 | previous_timestamp_index = -1 |
---|
1903 | |
---|
1904 | IF ( debug_output ) CALL debug_message( 'chem_emissions_header_init_lod2', 'end' ) |
---|
1905 | |
---|
1906 | END SUBROUTINE chem_emissions_header_init_lod2 |
---|
1907 | |
---|
1908 | ! |
---|
1909 | !-- 20200203 (ECC) |
---|
1910 | ! |
---|
1911 | !--------------------------------------------------------------------------------------------------! |
---|
1912 | ! Description: |
---|
1913 | ! ------------ |
---|
1914 | !> Reads emission data on demand for LOD2 |
---|
1915 | !--------------------------------------------------------------------------------------------------! |
---|
1916 | |
---|
1917 | SUBROUTINE chem_emissions_update_on_demand_lod2 |
---|
1918 | |
---|
1919 | USE control_parameters, & |
---|
1920 | ONLY: coupling_char, time_since_reference_point |
---|
1921 | |
---|
1922 | USE netcdf_data_input_mod, & |
---|
1923 | ONLY: chem_emis_att, close_input_file, get_variable, open_read_file |
---|
1924 | |
---|
1925 | USE arrays_3d, & |
---|
1926 | ONLY: hyp, pt |
---|
1927 | |
---|
1928 | USE surface_mod, & |
---|
1929 | ONLY: surf_def_h, surf_lsm_h, surf_usm_h |
---|
1930 | |
---|
1931 | |
---|
1932 | IMPLICIT NONE |
---|
1933 | |
---|
1934 | CHARACTER(LEN=80) :: this_timestamp !< writes out timestamp |
---|
1935 | |
---|
1936 | INTEGER(iwp) :: i,j,k,m !< generic counters |
---|
1937 | INTEGER(iwp) :: kmatch !< index of matched species |
---|
1938 | INTEGER(iwp) :: ncid !< netCDF file handle (chemistry file) |
---|
1939 | INTEGER(iwp) :: time_index_location !< location of most recent timestamp |
---|
1940 | |
---|
1941 | REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: cssws_def_h !< dummy default surface array |
---|
1942 | REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: cssws_lsm_h !< dummy LSM surface array |
---|
1943 | REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: cssws_usm_h !< dummy USM surface array |
---|
1944 | REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: mass2mole !< conversion factor mass 2 molar (ppm) flux |
---|
1945 | |
---|
1946 | REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: emis_distrib !< surface emissions |
---|
1947 | |
---|
1948 | REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:,:) :: emissions_raw !< raw emissions data |
---|
1949 | |
---|
1950 | IF ( debug_output ) CALL debug_message ( 'chem_emissions_update_on_demand_lod2', 'start' ) |
---|
1951 | ! |
---|
1952 | !-- Obtain current timestamp and locate index for most recent timestamp element |
---|
1953 | !-- end subroutine (RETURN) if it is still the same index as the existing time index |
---|
1954 | |
---|
1955 | this_timestamp = '' ! string must be initiated before using |
---|
1956 | CALL get_date_time( time_since_reference_point, date_time_str=this_timestamp ) |
---|
1957 | |
---|
1958 | time_index_location = chem_emissions_locate_timestep & |
---|
1959 | ( this_timestamp, timestamps, 1, chem_emis_att%dt_emission ) |
---|
1960 | |
---|
1961 | IF ( time_index_location == previous_timestamp_index ) RETURN |
---|
1962 | |
---|
1963 | ! |
---|
1964 | !-- Begin extract emission data for matched species from netCDF file |
---|
1965 | previous_timestamp_index = time_index_location |
---|
1966 | |
---|
1967 | ALLOCATE ( emis_distrib(n_matched_vars,nys:nyn,nxl:nxr) ) |
---|
1968 | emis_distrib = 0.0_wp |
---|
1969 | |
---|
1970 | ! |
---|
1971 | !-- Open netCDF file and allocate temp memory |
---|
1972 | CALL open_read_file( TRIM( input_file_chem ) // TRIM( coupling_char ), ncid ) |
---|
1973 | ALLOCATE( emissions_raw(1,1,nys:nyn,nxl:nxr,1) ) |
---|
1974 | |
---|
1975 | DO k = 1, n_matched_vars |
---|
1976 | ! |
---|
1977 | !-- Get index for matching species |
---|
1978 | kmatch = chem_emissions_locate_species( spc_names(match_spec_model(k)), & |
---|
1979 | chem_emis_att%species_name ) |
---|
1980 | ! |
---|
1981 | !-- Extract variable as-is |
---|
1982 | !-- Note C index notations for nx and ny due to MPI and reversed index dimension order for netCDF |
---|
1983 | !-- Fortran API) |
---|
1984 | emissions_raw = 0.0_wp |
---|
1985 | |
---|
1986 | CALL get_variable ( ncid, 'emission_values', emissions_raw, & |
---|
1987 | kmatch, nxl+1, nys+1, 1, time_index_location, & |
---|
1988 | 1, nxr-nxl+1, nyn-nys+1, 1, 1, .FALSE. ) |
---|
1989 | ! |
---|
1990 | !-- Transfer emission data |
---|
1991 | DO j = nys,nyn |
---|
1992 | DO i = nxl,nxr |
---|
1993 | emis_distrib(k,j,i) = emissions_raw(1,1,j,i,1) * conversion_factor |
---|
1994 | ENDDO |
---|
1995 | ENDDO |
---|
1996 | |
---|
1997 | ENDDO ! k = n_matched_vars |
---|
1998 | ! |
---|
1999 | !-- netCDF handle and temp memory no longer needed |
---|
2000 | DEALLOCATE( emissions_raw ) |
---|
2001 | CALL close_input_file( ncid ) |
---|
2002 | ! |
---|
2003 | !-- Set emis_dt since chemistry ODEs can be stiff, the option to solve them at every RK substep is |
---|
2004 | !-- present to help improve stability should the need arise |
---|
2005 | dt_emis = dt_3d |
---|
2006 | |
---|
2007 | IF ( call_chem_at_all_substeps ) dt_emis = dt_emis * weight_pres(intermediate_timestep_count) |
---|
2008 | ! |
---|
2009 | !-- Calculate conversion factor from mass flux to molar flux (mixing ratio) |
---|
2010 | ALLOCATE ( mass2mole(nys:nyn,nxl:nxr) ) |
---|
2011 | mass2mole = 0.0_wp |
---|
2012 | |
---|
2013 | DO i = nxl, nxr |
---|
2014 | DO j = nys, nyn |
---|
2015 | mass2mole(j,i) = mass_2_molar_flux ( hyp(nzb), pt(nzb,j,i) ) |
---|
2016 | ENDDO |
---|
2017 | ENDDO |
---|
2018 | |
---|
2019 | ! |
---|
2020 | !-- Calculate surface fluxes |
---|
2021 | !-- NOTE - For some reason I can not pass surf_xxx%cssws as output argument into subroutine |
---|
2022 | !-- assign_surface_flux ( ). The contents got mixed up once the subroutine is finished. I |
---|
2023 | !-- don't know why and I don't have time to investigate. As workaround I declared dummy |
---|
2024 | !-- variables and reassign them one by one (i.e., in a loop) |
---|
2025 | !-- ECC 20200206 |
---|
2026 | |
---|
2027 | ! |
---|
2028 | !-- Allocate and initialize dummy surface fluxes |
---|
2029 | ALLOCATE( cssws_def_h(n_matched_vars,surf_def_h(0)%ns) ) |
---|
2030 | cssws_def_h = 0.0_wp |
---|
2031 | |
---|
2032 | ALLOCATE( cssws_lsm_h(n_matched_vars,surf_lsm_h%ns) ) |
---|
2033 | cssws_lsm_h = 0.0_wp |
---|
2034 | |
---|
2035 | ALLOCATE( cssws_usm_h(n_matched_vars,surf_usm_h%ns) ) |
---|
2036 | cssws_usm_h = 0.0_wp |
---|
2037 | |
---|
2038 | ! |
---|
2039 | !-- Assign and transfer emissions as surface fluxes |
---|
2040 | CALL assign_surface_flux ( cssws_def_h, surf_def_h(0)%ns, & |
---|
2041 | surf_def_h(0)%j, surf_def_h(0)%i, & |
---|
2042 | emis_distrib, mass2mole ) |
---|
2043 | |
---|
2044 | |
---|
2045 | CALL assign_surface_flux ( cssws_lsm_h, surf_lsm_h%ns, & |
---|
2046 | surf_lsm_h%j, surf_lsm_h%i, & |
---|
2047 | emis_distrib, mass2mole ) |
---|
2048 | |
---|
2049 | |
---|
2050 | CALL assign_surface_flux ( cssws_usm_h, surf_usm_h%ns, & |
---|
2051 | surf_usm_h%j, surf_usm_h%i, & |
---|
2052 | emis_distrib, mass2mole ) |
---|
2053 | |
---|
2054 | DO k = 1, n_matched_vars |
---|
2055 | |
---|
2056 | DO m = 1, surf_def_h(0)%ns |
---|
2057 | surf_def_h(0)%cssws(k,m) = cssws_def_h(k,m) |
---|
2058 | ENDDO |
---|
2059 | |
---|
2060 | DO m = 1, surf_lsm_h%ns |
---|
2061 | surf_lsm_h%cssws(k,m) = cssws_lsm_h(k,m) |
---|
2062 | ENDDO |
---|
2063 | |
---|
2064 | DO m = 1, surf_usm_h%ns |
---|
2065 | surf_usm_h%cssws(k,m) = cssws_usm_h(k,m) |
---|
2066 | ENDDO |
---|
2067 | |
---|
2068 | ENDDO |
---|
2069 | |
---|
2070 | ! |
---|
2071 | !-- Cleaning up |
---|
2072 | DEALLOCATE( cssws_def_h ) |
---|
2073 | DEALLOCATE( cssws_lsm_h ) |
---|
2074 | DEALLOCATE( cssws_usm_h ) |
---|
2075 | |
---|
2076 | DEALLOCATE ( emis_distrib ) |
---|
2077 | DEALLOCATE ( mass2mole ) |
---|
2078 | |
---|
2079 | IF ( debug_output ) CALL debug_message ( 'chem_emissions_update_on_demand_lod2', 'end' ) |
---|
2080 | |
---|
2081 | END SUBROUTINE ! chem_emissions_update_on_demand_lod2 |
---|
2082 | |
---|
2083 | |
---|
2084 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2085 | !! |
---|
2086 | !! AUXILIARY SUBROUTINES |
---|
2087 | !! |
---|
2088 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2089 | |
---|
2090 | ! |
---|
2091 | !-- 20200203 (ECC) |
---|
2092 | ! |
---|
2093 | !--------------------------------------------------------------------------------------------------! |
---|
2094 | ! Description: |
---|
2095 | ! ------------ |
---|
2096 | !> Look for matched species between emissions attributes and selected chemical mechanisms and |
---|
2097 | !> determine corresponding molecular weights |
---|
2098 | !--------------------------------------------------------------------------------------------------! |
---|
2099 | |
---|
2100 | SUBROUTINE chem_emissions_init_species ( ncid ) |
---|
2101 | |
---|
2102 | USE netcdf_data_input_mod, & |
---|
2103 | ONLY: chem_emis_att, close_input_file, get_dimension_length, get_variable, open_read_file |
---|
2104 | |
---|
2105 | IMPLICIT NONE |
---|
2106 | |
---|
2107 | INTEGER(iwp) :: ispec !< generic counter 4 species |
---|
2108 | |
---|
2109 | INTEGER(iwp), INTENT(IN) :: ncid !< netcdf file ID |
---|
2110 | |
---|
2111 | IF ( debug_output ) CALL debug_message( 'chem_emissions_init_species', 'start' ) |
---|
2112 | ! |
---|
2113 | !- Assign species |
---|
2114 | CALL get_dimension_length ( ncid, chem_emis_att%n_emiss_species, 'nspecies' ) |
---|
2115 | CALL get_variable ( ncid, 'emission_name', chem_emis_att%species_name, & |
---|
2116 | chem_emis_att%n_emiss_species ) |
---|
2117 | ! |
---|
2118 | !- Backward compatibility for salsa_mod (ECC) |
---|
2119 | chem_emis_att%nspec = chem_emis_att%n_emiss_species |
---|
2120 | ! |
---|
2121 | !-- Get a list of matched species between emission_attributes and selected chemical mechanism |
---|
2122 | emission_output_required = .FALSE. |
---|
2123 | CALL chem_emissions_match( chem_emis_att, n_matched_vars ) |
---|
2124 | |
---|
2125 | ! |
---|
2126 | !-- If matched species found (at least 1), |
---|
2127 | !-- allocate memory for emission attributes, |
---|
2128 | !-- assign molecular masses [kg/mol], |
---|
2129 | !-- see chemistry_model_mod.f90 for reference. |
---|
2130 | IF ( n_matched_vars > 0 ) THEN |
---|
2131 | |
---|
2132 | emission_output_required = .TRUE. |
---|
2133 | |
---|
2134 | ALLOCATE( chem_emis_att%xm(n_matched_vars) ) |
---|
2135 | |
---|
2136 | DO ispec = 1, n_matched_vars |
---|
2137 | chem_emis_att%xm(ispec) = 1.0_wp |
---|
2138 | SELECT CASE ( TRIM( spc_names(match_spec_model(ispec)) ) ) |
---|
2139 | CASE ( 'SO2' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 2 |
---|
2140 | CASE ( 'SO4' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 4 |
---|
2141 | CASE ( 'NO' ); chem_emis_att%xm(ispec) = xm_N + xm_O |
---|
2142 | CASE ( 'NO2' ); chem_emis_att%xm(ispec) = xm_N + xm_O * 2 |
---|
2143 | CASE ( 'NH3' ); chem_emis_att%xm(ispec) = xm_N + xm_H * 3 |
---|
2144 | CASE ( 'CO' ); chem_emis_att%xm(ispec) = xm_C + xm_O |
---|
2145 | CASE ( 'CO2' ); chem_emis_att%xm(ispec) = xm_C + xm_O * 2 |
---|
2146 | CASE ( 'CH4' ); chem_emis_att%xm(ispec) = xm_C + xm_H * 4 |
---|
2147 | CASE ( 'HNO3' ); chem_emis_att%xm(ispec) = xm_H + xm_N + xm_O*3 |
---|
2148 | END SELECT |
---|
2149 | ENDDO |
---|
2150 | |
---|
2151 | ENDIF ! IF ( n_matched_vars > 0 ) |
---|
2152 | |
---|
2153 | IF ( debug_output ) CALL debug_message( 'chem_emissions_init_species', 'end' ) |
---|
2154 | |
---|
2155 | END SUBROUTINE chem_emissions_init_species |
---|
2156 | |
---|
2157 | |
---|
2158 | ! |
---|
2159 | !-- 20200203 (ECC) |
---|
2160 | ! |
---|
2161 | !--------------------------------------------------------------------------------------------------! |
---|
2162 | ! Description: |
---|
2163 | ! ------------ |
---|
2164 | !> Extract timestamps from netCDF input |
---|
2165 | !--------------------------------------------------------------------------------------------------! |
---|
2166 | |
---|
2167 | SUBROUTINE chem_emissions_init_timestamps ( ncid ) |
---|
2168 | |
---|
2169 | USE control_parameters, & |
---|
2170 | ONLY: message_string |
---|
2171 | |
---|
2172 | USE netcdf_data_input_mod, & |
---|
2173 | ONLY: chem_emis_att, close_input_file, get_dimension_length, get_variable, open_read_file |
---|
2174 | |
---|
2175 | IMPLICIT NONE |
---|
2176 | |
---|
2177 | INTEGER(iwp) :: fld_len !< string field length |
---|
2178 | INTEGER(iwp) :: itime !< generic counter (4 species) |
---|
2179 | |
---|
2180 | INTEGER(iwp), INTENT(IN) :: ncid !< netcdf file handle |
---|
2181 | |
---|
2182 | IF ( debug_output ) CALL debug_message( 'chem_emissions_read_timestamps', 'start' ) |
---|
2183 | ! |
---|
2184 | !-- Import timestamps from netCDF input |
---|
2185 | CALL get_dimension_length ( ncid, chem_emis_att%dt_emission, 'time' ) |
---|
2186 | CALL get_dimension_length ( ncid, fld_len, 'field_length' ) |
---|
2187 | CALL get_variable ( ncid, 'timestamp', timestamps, chem_emis_att%dt_emission, fld_len ) |
---|
2188 | ! |
---|
2189 | !-- Throw error at first instance of timestamps not in listed in chronological order. |
---|
2190 | DO itime = 2,chem_emis_att%dt_emission |
---|
2191 | |
---|
2192 | IF ( timestamps(itime-1) > timestamps(itime) ) THEN |
---|
2193 | |
---|
2194 | WRITE( message_string, * ) & |
---|
2195 | 'input timestamps not in chronological order for', & |
---|
2196 | CHAR( 10 ), ' ', & |
---|
2197 | 'index ', (itime-1), ' : ', TRIM( timestamps(itime-1) ), ' and', & |
---|
2198 | CHAR( 10 ), ' ', & |
---|
2199 | 'index ', (itime), ' : ', TRIM( timestamps(itime) ) |
---|
2200 | |
---|
2201 | CALL message( 'chem_emissions_read_timestamps', 'CM0469', 1, 2, 0, 6, 0 ) |
---|
2202 | |
---|
2203 | ENDIF |
---|
2204 | |
---|
2205 | ENDDO |
---|
2206 | |
---|
2207 | IF ( debug_output ) CALL debug_message( 'chem_emissions_read_timestamps', 'end' ) |
---|
2208 | |
---|
2209 | END SUBROUTINE chem_emissions_init_timestamps |
---|
2210 | |
---|
2211 | |
---|
2212 | ! |
---|
2213 | !-- 20200203 (ECC) |
---|
2214 | ! |
---|
2215 | !--------------------------------------------------------------------------------------------------! |
---|
2216 | ! Description: |
---|
2217 | ! ------------ |
---|
2218 | !> Assign emissions as surface fluxes |
---|
2219 | ! |
---|
2220 | !> NOTE: For arguments, I originally wanted to use unspecified dimensions, but I could not get |
---|
2221 | !> this to work properly, hence the dimensioned array arguments. |
---|
2222 | !--------------------------------------------------------------------------------------------------! |
---|
2223 | |
---|
2224 | SUBROUTINE assign_surface_flux ( surf_array, nsurfs, surf_j, surf_i, emis_dist, conv_mole ) |
---|
2225 | |
---|
2226 | USE arrays_3d, & |
---|
2227 | ONLY: rho_air |
---|
2228 | |
---|
2229 | USE netcdf_data_input_mod, & |
---|
2230 | ONLY: chem_emis_att |
---|
2231 | |
---|
2232 | USE surface_mod !< for surf_type |
---|
2233 | |
---|
2234 | IMPLICIT NONE |
---|
2235 | ! |
---|
2236 | !-- Input arguments |
---|
2237 | INTEGER(iwp), INTENT(IN) :: nsurfs !< # surfaces in surf_array |
---|
2238 | |
---|
2239 | INTEGER(iwp), DIMENSION(nsurfs), INTENT(IN) :: surf_i !< i indices 4 surf. elements |
---|
2240 | INTEGER(iwp), DIMENSION(nsurfs), INTENT(IN) :: surf_j !< j indices 4 surf. elements |
---|
2241 | |
---|
2242 | REAL(wp), DIMENSION(nys:nyn,nxl:nxr), INTENT(IN) :: conv_mole !< conv. 2 molar flux |
---|
2243 | REAL(wp), DIMENSION(n_matched_vars,nys:nyn,nxl:nxr), INTENT(IN) :: emis_dist !< surf. emissions |
---|
2244 | |
---|
2245 | REAL(wp), DIMENSION(n_matched_vars,nsurfs), INTENT(INOUT) :: surf_array !< surface listing |
---|
2246 | |
---|
2247 | ! |
---|
2248 | !-- Parameters (magic numbers) |
---|
2249 | CHARACTER(LEN=2), PARAMETER :: sp_PM = 'PM' !< id string for all PMs |
---|
2250 | CHARACTER(LEN=3), PARAMETER :: sp_VOC = 'VOC' !< id string for VOC |
---|
2251 | |
---|
2252 | REAL(wp), PARAMETER :: mol2ppm = 1.0E+06_wp !< conversion from mole 2 ppm |
---|
2253 | ! |
---|
2254 | !-- Local variables |
---|
2255 | CHARACTER(LEN=80) :: this_species_name !< matched species name |
---|
2256 | |
---|
2257 | INTEGER(iwp) :: i,j,k,m !< generic counters |
---|
2258 | |
---|
2259 | REAL(wp) :: flux_conv_factor !< conversion factor |
---|
2260 | |
---|
2261 | IF ( debug_output ) CALL debug_message( 'chem_emissions_header_init_lod2', 'start' ) |
---|
2262 | |
---|
2263 | DO k = 1, n_matched_vars |
---|
2264 | |
---|
2265 | this_species_name = spc_names(k) !< species already matched |
---|
2266 | |
---|
2267 | DO m = 1, nsurfs |
---|
2268 | |
---|
2269 | j = surf_j(m) ! get surface coordinates |
---|
2270 | i = surf_i(m) |
---|
2271 | |
---|
2272 | ! |
---|
2273 | !-- Calculate conversion factor depending on emission species type |
---|
2274 | flux_conv_factor = rho_air(nzb) |
---|
2275 | ! |
---|
2276 | !-- Account for conversion to different types of emisison species |
---|
2277 | IF ( TRIM( this_species_name( 1:LEN( sp_PM ) ) ) == sp_PM ) THEN |
---|
2278 | |
---|
2279 | ! do nothing (use mass flux directly) |
---|
2280 | |
---|
2281 | ELSE IF ( TRIM( this_species_name( 1:LEN( sp_VOC ) ) ) == sp_VOC ) THEN |
---|
2282 | |
---|
2283 | flux_conv_factor = flux_conv_factor * conv_mole(j,i) * mol2ppm |
---|
2284 | |
---|
2285 | ELSE |
---|
2286 | |
---|
2287 | flux_conv_factor = flux_conv_factor * conv_mole(j,i) * mol2ppm / chem_emis_att%xm(k) |
---|
2288 | |
---|
2289 | ENDIF |
---|
2290 | ! |
---|
2291 | !-- Finally assign surface flux |
---|
2292 | surf_array(k,m) = emis_dist(k,j,i) * flux_conv_factor |
---|
2293 | |
---|
2294 | ENDDO ! m = 1, nsurfs |
---|
2295 | |
---|
2296 | ENDDO ! k = 1, n_matched_vars |
---|
2297 | |
---|
2298 | |
---|
2299 | IF ( debug_output ) CALL debug_message( 'chem_emissions_header_init_lod2', 'end' ) |
---|
2300 | |
---|
2301 | END SUBROUTINE assign_surface_flux |
---|
2302 | |
---|
2303 | |
---|
2304 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2305 | !! |
---|
2306 | !! AUXILIARY FUNCTIONS |
---|
2307 | !! |
---|
2308 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2309 | |
---|
2310 | ! |
---|
2311 | !-- 20200203 (ECC) |
---|
2312 | ! |
---|
2313 | !--------------------------------------------------------------------------------------------------! |
---|
2314 | ! Description: |
---|
2315 | ! ------------ |
---|
2316 | !> Given incoming flux units ( mass / area / time ) provide single-valued onversion factor to |
---|
2317 | !> ( kg / m2 / s ) |
---|
2318 | !--------------------------------------------------------------------------------------------------! |
---|
2319 | |
---|
2320 | FUNCTION chem_emissions_convert_base_units ( units_in ) RESULT ( conv_factor ) |
---|
2321 | |
---|
2322 | IMPLICIT NONE |
---|
2323 | ! |
---|
2324 | !-- Function arguments |
---|
2325 | CHARACTER(LEN=*), INTENT(IN) :: units_in !< incoming units (ie emt_att%units) |
---|
2326 | |
---|
2327 | REAL(wp) :: conv_factor !< convertion factor |
---|
2328 | |
---|
2329 | ! |
---|
2330 | !-- Parameters (magic numbers) |
---|
2331 | INTEGER(iwp), PARAMETER :: up2lo = 32 !< convert letter to lower case |
---|
2332 | ! |
---|
2333 | !-- Base unit conversion factors (should be self-explanatory) |
---|
2334 | REAL(wp), PARAMETER :: hour_per_year = 8760.0_wp |
---|
2335 | REAL(wp), PARAMETER :: g_to_kg = 1.0E-03_wp |
---|
2336 | REAL(wp), PARAMETER :: miug_to_kg = 1.0E-09_wp |
---|
2337 | REAL(wp), PARAMETER :: s_per_hour = 3600.0_wp |
---|
2338 | REAL(wp), PARAMETER :: s_per_day = 86400.0_wp |
---|
2339 | REAL(wp), PARAMETER :: tons_to_kg = 100.0_wp |
---|
2340 | |
---|
2341 | REAL(wp), PARAMETER :: day_to_s = 1.0_wp / s_per_day |
---|
2342 | REAL(wp), PARAMETER :: hour_to_s = 1.0_wp / s_per_hour |
---|
2343 | REAL(wp), PARAMETER :: year_to_s = 1.0_wp / s_per_hour / hour_per_year |
---|
2344 | |
---|
2345 | |
---|
2346 | ! |
---|
2347 | !-- Local variables |
---|
2348 | CHARACTER(LEN=LEN(units_in)) :: units_in_lo !< units in lower case |
---|
2349 | |
---|
2350 | INTEGER(iwp) :: j,k !< generic counters |
---|
2351 | INTEGER(iwp) :: str_len !< length of unit string |
---|
2352 | ! |
---|
2353 | !-- Turn units string to lower case |
---|
2354 | units_in_lo = '' |
---|
2355 | str_len = LEN( TRIM( units_in ) ) |
---|
2356 | |
---|
2357 | DO k = 1,str_len |
---|
2358 | j = IACHAR( units_in(k:k) ) |
---|
2359 | units_in_lo(k:k) = ACHAR( j ) |
---|
2360 | IF ( ( j >= IACHAR( "A" ) ) .AND. ( j <= IACHAR( "Z" ) ) ) & |
---|
2361 | units_in_lo(k:k) = ACHAR( j + up2lo ) |
---|
2362 | ENDDO |
---|
2363 | |
---|
2364 | conv_factor = 1.0_wp !< default value (kg/m2/s) |
---|
2365 | |
---|
2366 | SELECT CASE ( TRIM( units_in_lo ) ) |
---|
2367 | CASE ( 'kg/m2/s' ); conv_factor = 1.0_wp |
---|
2368 | CASE ( 'kg/m2/hour' ); conv_factor = hour_to_s |
---|
2369 | CASE ( 'kg/m2/day' ); conv_factor = day_to_s |
---|
2370 | CASE ( 'kg/m2/year' ); conv_factor = year_to_s |
---|
2371 | CASE ( 'ton/m2/s' ); conv_factor = tons_to_kg |
---|
2372 | CASE ( 'ton/m2/hour' ); conv_factor = tons_to_kg * hour_to_s |
---|
2373 | CASE ( 'ton/m2/day' ); conv_factor = tons_to_kg * day_to_s |
---|
2374 | CASE ( 'ton/m2/year' ); conv_factor = tons_to_kg * year_to_s |
---|
2375 | CASE ( 'g/m2/s' ); conv_factor = g_to_kg |
---|
2376 | CASE ( 'g/m2/hour' ); conv_factor = g_to_kg * hour_to_s |
---|
2377 | CASE ( 'g/m2/day' ); conv_factor = g_to_kg * day_to_s |
---|
2378 | CASE ( 'g/m2/year' ); conv_factor = g_to_kg * year_to_s |
---|
2379 | CASE ( 'micrograms/m2/s' ); conv_factor = miug_to_kg |
---|
2380 | CASE ( 'micrograms/m2/hour' ); conv_factor = miug_to_kg * hour_to_s |
---|
2381 | CASE ( 'micrograms/m2/day' ); conv_factor = miug_to_kg * day_to_s |
---|
2382 | CASE ( 'micrograms/m2/year' ); conv_factor = miug_to_kg * year_to_s |
---|
2383 | CASE DEFAULT |
---|
2384 | message_string = '' ! to get around unused variable warning / error |
---|
2385 | WRITE ( message_string, * ) 'Specified emission units (', TRIM( units_in ), & |
---|
2386 | ') not recognized in PALM-4U' |
---|
2387 | CALL message ( 'chem_emission_convert_units', 'CM0446', 2, 2, 0, 6, 0 ) |
---|
2388 | END SELECT |
---|
2389 | |
---|
2390 | END FUNCTION chem_emissions_convert_base_units |
---|
2391 | |
---|
2392 | |
---|
2393 | ! |
---|
2394 | !-- 20200203 (ECC) |
---|
2395 | ! |
---|
2396 | !--------------------------------------------------------------------------------------------------! |
---|
2397 | ! Description: |
---|
2398 | ! ------------ |
---|
2399 | !> Calculates conversion factor from mass flux to ppm (molar flux) |
---|
2400 | !--------------------------------------------------------------------------------------------------! |
---|
2401 | |
---|
2402 | FUNCTION mass_2_molar_flux ( rhogh, theta ) RESULT ( conv_factor ) |
---|
2403 | |
---|
2404 | USE basic_constants_and_equations_mod, & |
---|
2405 | ONLY: p_0, rd_d_cp, rgas_univ |
---|
2406 | |
---|
2407 | IMPLICIT NONE |
---|
2408 | ! |
---|
2409 | !-- Function arguments |
---|
2410 | REAL(wp) :: conv_factor !< conversion factor |
---|
2411 | REAL(wp), INTENT(IN) :: rhogh !< hydrostatic pressure |
---|
2412 | REAL(wp), INTENT(IN) :: theta !< potential temperature |
---|
2413 | |
---|
2414 | conv_factor = ( rgas_univ / rhogh ) * theta * ( ( rhogh / p_0 ) ** rd_d_cp ) |
---|
2415 | |
---|
2416 | END FUNCTION mass_2_molar_flux |
---|
2417 | |
---|
2418 | |
---|
2419 | ! |
---|
2420 | !-- 20200203 (ECC) |
---|
2421 | ! |
---|
2422 | !--------------------------------------------------------------------------------------------------! |
---|
2423 | ! Description: |
---|
2424 | ! ------------ |
---|
2425 | !> Given target sepecies locate index in species array |
---|
2426 | !> returns 0 if none is found |
---|
2427 | !--------------------------------------------------------------------------------------------------! |
---|
2428 | |
---|
2429 | FUNCTION chem_emissions_locate_species ( this_species, species_array ) RESULT ( species_index ) |
---|
2430 | |
---|
2431 | IMPLICIT NONE |
---|
2432 | ! |
---|
2433 | !-- Function arguments |
---|
2434 | INTEGER(iwp) :: species_index !> index matching species |
---|
2435 | |
---|
2436 | CHARACTER(LEN=25), INTENT(IN) :: species_array(:) !> array of species |
---|
2437 | CHARACTER(LEN=*), INTENT(IN) :: this_species !> target species |
---|
2438 | ! |
---|
2439 | !-- Local variables |
---|
2440 | INTEGER(iwp) :: k !> generic counter |
---|
2441 | INTEGER(iwp) :: n_species !> number of species in species_array |
---|
2442 | |
---|
2443 | n_species = SIZE( species_array, 1 ) |
---|
2444 | |
---|
2445 | DO k = 1, n_species |
---|
2446 | IF ( TRIM( species_array(k) ) == TRIM( this_species ) ) EXIT |
---|
2447 | ENDDO |
---|
2448 | |
---|
2449 | species_index = 0 !> assume no matching index is found |
---|
2450 | |
---|
2451 | IF ( TRIM( species_array(k) ) == TRIM( this_species ) ) specieS_index = k |
---|
2452 | |
---|
2453 | END FUNCTION chem_emissions_locate_species |
---|
2454 | |
---|
2455 | |
---|
2456 | ! |
---|
2457 | !-- 20200203 (ECC) |
---|
2458 | ! |
---|
2459 | !--------------------------------------------------------------------------------------------------! |
---|
2460 | ! Description: |
---|
2461 | ! ------------ |
---|
2462 | !> given target timestamp locate most recent timestep in timestamp array |
---|
2463 | !> using bisection search (since array is sorted) |
---|
2464 | !--------------------------------------------------------------------------------------------------! |
---|
2465 | |
---|
2466 | RECURSIVE FUNCTION chem_emissions_locate_timestep & |
---|
2467 | ( this_timestamp, timestamp_array, lower_bound, upper_bound ) & |
---|
2468 | RESULT ( timestamp_index ) |
---|
2469 | |
---|
2470 | ! |
---|
2471 | !-- Function arguments |
---|
2472 | CHARACTER(LEN=*), INTENT(IN) :: this_timestamp !> target timestamp |
---|
2473 | CHARACTER(LEN=512), INTENT(IN) :: timestamp_array(:) !> array of timestamps |
---|
2474 | |
---|
2475 | INTEGER(iwp), INTENT(IN) :: lower_bound, upper_bound !> timestamp_array index bounds |
---|
2476 | |
---|
2477 | INTEGER(iwp) :: timestamp_index !> index for most recent timestamp in timestamp_array |
---|
2478 | |
---|
2479 | ! |
---|
2480 | !-- Local variables |
---|
2481 | INTEGER(iwp) :: k0,km,k1 !> lower, central, and upper index bounds |
---|
2482 | ! |
---|
2483 | !-- Assign bounds |
---|
2484 | k0 = lower_bound |
---|
2485 | k1 = upper_bound |
---|
2486 | ! |
---|
2487 | !-- Make sure k1 is always not smaller than k0 |
---|
2488 | IF ( k0 > k1 ) THEN |
---|
2489 | k0 = upper_bound |
---|
2490 | k1 = lower_bound |
---|
2491 | ENDIF |
---|
2492 | ! |
---|
2493 | !-- Make sure k0 and k1 stay within array bounds by timestamp_array |
---|
2494 | IF ( k0 < 1 ) k0 = 1 |
---|
2495 | IF ( k1 > SIZE( timestamp_array, 1 ) ) k1 = SIZE( timestamp_array, 1 ) |
---|
2496 | ! |
---|
2497 | !-- Terminate if target is contained within 2 consecutive indices otherwise calculate central bound |
---|
2498 | !-- (km) and determine new index bounds for the next iteration |
---|
2499 | |
---|
2500 | IF ( ( k1 - k0 ) > 1 ) THEN |
---|
2501 | km = ( k0 + k1 ) / 2 |
---|
2502 | IF ( TRIM( this_timestamp ) > TRIM( timestamp_array(km) ) ) THEN |
---|
2503 | k0 = km |
---|
2504 | ELSE |
---|
2505 | k1 = km |
---|
2506 | ENDIF |
---|
2507 | timestamp_index = chem_emissions_locate_timestep( this_timestamp, timestamp_array, k0, k1 ) |
---|
2508 | ELSE |
---|
2509 | timestamp_index = k0 |
---|
2510 | ENDIF |
---|
2511 | |
---|
2512 | END FUNCTION chem_emissions_locate_timestep |
---|
2513 | |
---|
2514 | |
---|
2515 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2516 | !! |
---|
2517 | !! END OF MODULE |
---|
2518 | !! |
---|
2519 | !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
---|
2520 | |
---|
2521 | END MODULE chem_emissions_mod |
---|