Home

Context Navigation

source: palm/trunk/SOURCE/chem_emissions_mod.f90 @ 4144

Last change on this file since 4144 was 4144, checked in by raasch, 5 years ago
relational operators .EQ., .NE., etc. replaced by ==, /=, etc.
Property svn:keywords set to `Id`
File size: 64.5 KB

Line
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
6	! terms of the GNU General Public License as published by the Free Software
7	! Foundation, either version 3 of the License, or (at your option) any later
8	! version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 2018-2019 Leibniz Universitaet Hannover
18	! Copyright 2018-2019 Freie Universitaet Berlin
19	! Copyright 2018-2019 Karlsruhe Institute of Technology
20	!--------------------------------------------------------------------------------!
21	!
22	! Current revisions:
23	! ------------------
24	!
25	!
26	! Former revisions:
27	! -----------------
28	! $Id: chem_emissions_mod.f90 4144 2019-08-06 09:11:47Z raasch $
29	! relational operators .EQ., .NE., etc. replaced by ==, /=, etc.
30	!
31	! 4055 2019-06-27 09:47:29Z suehring
32	! - replaced local thermo. constants w/ module definitions in
33	! basic_constants_and_equations_mod (rgas_univ, p_0, r_d_cp)
34	! - initialize array emis_distribution immediately following allocation
35	! - lots of minor formatting changes based on review sesson in 20190325
36	! (E.C. Chan)
37	!
38	! 3968 2019-05-13 11:04:01Z suehring
39	! - in subroutine chem_emissions_match replace all decision structures relating to
40	! mode_emis to emiss_lod
41	! - in subroutine chem_check_parameters replace emt%nspec with emt%n_emiss_species
42	! - spring cleaning (E.C. Chan)
43	!
44	! 3885 2019-04-11 11:29:34Z kanani
45	! Changes related to global restructuring of location messages and introduction
46	! of additional debug messages
47	!
48	! 3831 2019-03-28 09:11:22Z forkel
49	! added nvar to USE chem_gasphase_mod (chem_modules will not include nvar anymore)
50	!
51	! 3788 2019-03-07 11:40:09Z banzhafs
52	! Removed unused variables from chem_emissions_mod
53	!
54	! 3772 2019-02-28 15:51:57Z suehring
55	! - In case of parametrized emissions, assure that emissions are only on natural
56	! surfaces (i.e. streets) and not on urban surfaces.
57	! - some unnecessary if clauses removed
58	!
59	! 3685 2019 -01-21 01:02:11Z knoop
60	! Some interface calls moved to module_interface + cleanup
61	!
62	! 3611 2018-12-07 14:14:11Z banzhafs
63	! Code update to comply PALM coding rules
64	! removed unnecessary informative messsages/location
65	! messages and corrected comments on PM units from to kg
66	! bug fix: spcs_name replaced by nvar in DO loops
67	!
68	! 3591 2018-11-30 17:37:32Z suehring
69	! - Removed salsa dependency.
70	! - Enabled PARAMETRIZED mode for default surfaces when LSM is not applied but
71	! salsa is (M. Kurppa)
72	!
73	! 3582 2018-11-29 19:16:36Z suehring
74	! resler:
75	! Break lines at 132 characters
76	!
77	! 3483 2018-11-02 14:19:26Z raasch
78	! bugfix: wrong locations of netCDF directives fixed
79	!
80	! 3467 2018-10-30 19:05:21Z suehring
81	! Enabled PARAMETRIZED mode for default surfaces when LSM is not applied but
82	! salsa is used
83	!
84	! 3458 2018-10-30 14:51:23Z kanani
85	! from chemistry branch r3443, banzhafs, Russo:
86	! Additional correction for index of input file of pre-processed mode
87	! Removed atomic and molecular weights as now available in chem_modules and
88	! added par_emis_time_factor (formerly in netcdf_data_input_mod)
89	! Added correction for index of input file of pre-processed mode
90	! Added a condition for default mode necessary for information read from ncdf file
91	! in pre-processed and default mode
92	! Correction of multiplying factor necessary for scaling emission values in time
93	! Introduced correction for scaling units in the case of DEFAULT emission mode
94	!
95	! 3373 2018-10-18 15:25:56Z kanani
96	! Fix wrong location of __netcdf directive
97	!
98	! 3337 2018-10-12 15:17:09Z kanani
99	! (from branch resler)
100	! Formatting
101	!
102	! 3312 2018-10-06 14:15:46Z knoop
103	! Initial revision
104	!
105	! 3286 2018-09-28 07:41:39Z forkel
106	!
107	! Authors:
108	! --------
109	! @author Emmanuele Russo (FU-Berlin)
110	! @author Sabine Banzhaf (FU-Berlin)
111	! @author Martijn Schaap (FU-Berlin, TNO Utrecht)
112	!
113	! Description:
114	! ------------
115	!> MODULE for reading-in Chemistry Emissions data
116	!>
117	!> @todo Rename nspec to n_emis to avoid inteferece with nspec from chem_gasphase_mod
118	!> @todo Check_parameters routine should be developed: for now it includes just one condition
119	!> @todo Use of Restart files not contempled at the moment
120	!> @todo revise indices of files read from the netcdf: preproc_emission_data and expert_emission_data
121	!> @todo for now emission data may be passed on a singular vertical level: need to be more flexible
122	!> @todo fill/activate restart option in chem_emissions_init
123	!> @todo discuss dt_emis
124	!> @note <Enter notes on the module>
125	!> @bug <Enter known bugs here>
126	!------------------------------------------------------------------------------!
127
128	MODULE chem_emissions_mod
129
130	USE arrays_3d, &
131	ONLY: rho_air
132
133	USE basic_constants_and_equations_mod, &
134	ONLY: rgas_univ, p_0, rd_d_cp
135
136	USE control_parameters, &
137	ONLY: debug_output, &
138	end_time, message_string, initializing_actions, &
139	intermediate_timestep_count, dt_3d
140
141	USE indices
142
143	USE kinds
144
145	#if defined( __netcdf )
146	USE netcdf
147	#endif
148
149	USE netcdf_data_input_mod, &
150	ONLY: chem_emis_att_type, chem_emis_val_type
151
152	USE date_and_time_mod, &
153	ONLY: day_of_month, hour_of_day, &
154	index_mm, index_dd, index_hh, &
155	month_of_year, hour_of_day, &
156	time_default_indices, time_preprocessed_indices
157
158	USE chem_gasphase_mod, &
159	ONLY: nvar, spc_names
160
161	USE chem_modules
162
163	USE statistics, &
164	ONLY: weight_pres
165
166
167	IMPLICIT NONE
168
169	!
170	!-- Declare all global variables within the module
171
172	CHARACTER (LEN=80) :: filename_emis !< Variable for the name of the netcdf input file
173
174	INTEGER(iwp) :: dt_emis !< Time Step Emissions
175	INTEGER(iwp) :: i !< index 1st selected dimension (some dims are not spatial)
176	INTEGER(iwp) :: j !< index 2nd selected dimension
177	INTEGER(iwp) :: i_start !< Index to start read variable from netcdf along one dims
178	INTEGER(iwp) :: i_end !< Index to end read variable from netcdf in one dims
179	INTEGER(iwp) :: j_start !< Index to start read variable from netcdf in additional dims
180	INTEGER(iwp) :: j_end !< Index to end read variable from netcdf in additional dims
181	INTEGER(iwp) :: len_index !< length of index (used for several indices)
182	INTEGER(iwp) :: len_index_pm !< length of PMs index
183	INTEGER(iwp) :: len_index_voc !< length of voc index
184	INTEGER(iwp) :: z_start !< Index to start read variable from netcdf in additional dims
185	INTEGER(iwp) :: z_end !< Index to end read variable from netcdf in additional dims
186
187	REAL(wp) :: conversion_factor !< Units Conversion Factor
188
189	SAVE
190
191	!
192	!-- Checks Input parameters
193	INTERFACE chem_emissions_check_parameters
194	MODULE PROCEDURE chem_emissions_check_parameters
195	END INTERFACE chem_emissions_check_parameters
196	!
197	!-- Matching Emissions actions
198	INTERFACE chem_emissions_match
199	MODULE PROCEDURE chem_emissions_match
200	END INTERFACE chem_emissions_match
201	!
202	!-- Initialization actions
203	INTERFACE chem_emissions_init
204	MODULE PROCEDURE chem_emissions_init
205	END INTERFACE chem_emissions_init
206	!
207	!-- Setup of Emissions
208	INTERFACE chem_emissions_setup
209	MODULE PROCEDURE chem_emissions_setup
210	END INTERFACE chem_emissions_setup
211
212	PUBLIC chem_emissions_init, chem_emissions_match, chem_emissions_setup
213	!
214	!-- Public Variables
215	PUBLIC conversion_factor, len_index, len_index_pm, len_index_voc
216
217	CONTAINS
218
219	!------------------------------------------------------------------------------!
220	! Description:
221	! ------------
222	!> Routine for checking input parameters
223	!------------------------------------------------------------------------------!
224	SUBROUTINE chem_emissions_check_parameters
225
226	IMPLICIT NONE
227
228	TYPE(chem_emis_att_type) :: emt
229
230	!
231	!-- Check if species count matches the number of names
232	!-- passed for the chemiscal species
233
234	IF ( SIZE(emt%species_name) /= emt%n_emiss_species ) THEN
235	! IF ( SIZE(emt%species_name) /= emt%nspec ) THEN
236
237	message_string = 'Numbers of input emission species names and number of species' // &
238	'for which emission values are given do not match'
239	CALL message( 'chem_emissions_check_parameters', 'CM0437', 2, 2, 0, 6, 0 )
240
241	ENDIF
242
243	END SUBROUTINE chem_emissions_check_parameters
244
245
246	!------------------------------------------------------------------------------!
247	! Description:
248	! ------------
249	!> Matching the chemical species indices. The routine checks what are the
250	!> indices of the emission input species and the corresponding ones of the
251	!> model species. The routine gives as output a vector containing the number
252	!> of common species: it is important to note that while the model species
253	!> are distinct, their values could be given to a single species in input.
254	!> For example, in the case of NO2 and NO, values may be passed in input as
255	!> NOX values.
256	!------------------------------------------------------------------------------!
257
258	SUBROUTINE chem_emissions_match( emt_att,len_index )
259
260	INTEGER(iwp) :: ind_inp !< Parameters for cycling through chemical input species
261	INTEGER(iwp) :: ind_mod !< Parameters for cycling through chemical model species
262	INTEGER(iwp) :: ind_voc !< Indices to check whether a split for voc should be done
263	INTEGER(iwp) :: ispec !< index for cycle over effective number of emission species
264	INTEGER(iwp) :: nspec_emis_inp !< Variable where to store # of emission species in input
265
266	INTEGER(iwp), INTENT(INOUT) :: len_index !< number of common species between input dataset & model species
267
268	TYPE(chem_emis_att_type), INTENT(INOUT) :: emt_att !< Chemistry Emission Array (decl. netcdf_data_input.f90)
269
270
271	IF ( debug_output ) CALL debug_message( 'chem_emissions_match', 'start' )
272
273	!
274	!-- Number of input emission species
275
276	nspec_emis_inp = emt_att%n_emiss_species
277	! nspec_emis_inp=emt_att%nspec
278
279	!
280	!-- Check the emission LOD: 0 (PARAMETERIZED), 1 (DEFAULT), 2 (PRE-PROCESSED)
281	!
282	SELECT CASE (emiss_lod)
283
284	!
285	!-- LOD 0 (PARAMETERIZED mode)
286
287	CASE (0)
288
289	len_index = 0
290
291	! number of species and number of matched species can be different
292	! but call is only made if both are greater than zero
293
294	IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN
295
296	!
297	!-- Cycle over model species
298
299	DO ind_mod = 1, nvar
300	ind_inp = 1
301	DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) !< 'novalue' is the default
302	IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN
303	len_index = len_index + 1
304	ENDIF
305	ind_inp = ind_inp + 1
306	ENDDO
307	ENDDO
308
309	IF ( len_index > 0 ) THEN
310
311	!
312	!-- Allocation of Arrays of the matched species
313
314	ALLOCATE ( match_spec_input(len_index) )
315	ALLOCATE ( match_spec_model(len_index) )
316
317	!
318	!-- Pass species indices to declared arrays
319
320	len_index = 0
321
322	DO ind_mod = 1, nvar
323	ind_inp = 1
324	DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' )
325	IF ( TRIM(surface_csflux_name(ind_inp)) == &
326	TRIM(spc_names(ind_mod)) ) THEN
327	len_index = len_index + 1
328	match_spec_input(len_index) = ind_inp
329	match_spec_model(len_index) = ind_mod
330	ENDIF
331	ind_inp = ind_inp + 1
332	END DO
333	END DO
334
335	!
336	!-- Check
337
338	DO ispec = 1, len_index
339
340	IF ( emiss_factor_main(match_spec_input(ispec) ) < 0 .AND. &
341	emiss_factor_side(match_spec_input(ispec) ) < 0 ) THEN
342
343	message_string = 'PARAMETERIZED emissions mode selected:' // &
344	' EMISSIONS POSSIBLE ONLY ON STREET SURFACES' // &
345	' but values of scaling factors for street types' // &
346	' emiss_factor_main AND emiss_factor_side' // &
347	' not provided for each of the emissions species' // &
348	' or not provided at all: PLEASE set a finite value' // &
349	' for these parameters in the chemistry namelist'
350	CALL message( 'chem_emissions_matching', 'CM0442', 2, 2, 0, 6, 0 )
351
352	ENDIF
353
354	END DO
355
356
357	ELSE
358
359	message_string = 'Non of given Emission Species' // &
360	' matches' // &
361	' model chemical species' // &
362	' Emission routine is not called'
363	CALL message( 'chem_emissions_matching', 'CM0443', 0, 0, 0, 6, 0 )
364
365	ENDIF
366
367	ELSE
368
369	message_string = 'Array of Emission species not allocated: ' // &
370	' Either no emission species are provided as input or' // &
371	' no chemical species are used by PALM.' // &
372	' Emission routine is not called'
373	CALL message( 'chem_emissions_matching', 'CM0444', 0, 2, 0, 6, 0 )
374
375	ENDIF
376
377	!
378	!-- LOD 1 (DEFAULT mode)
379
380	CASE (1)
381
382	len_index = 0 ! total number of species (to be accumulated)
383	len_index_voc = 0 ! total number of VOCs (to be accumulated)
384	len_index_pm = 3 ! total number of PMs: PM1, PM2.5, PM10.
385
386	!
387	!-- number of model species and input species could be different
388	!-- but process this only when both are non-zero
389
390	IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN
391
392	!
393	!-- Cycle over model species
394	DO ind_mod = 1, nvar
395
396	!
397	!-- Cycle over input species
398
399	DO ind_inp = 1, nspec_emis_inp
400
401	!
402	!-- Check for VOC Species
403
404	IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN
405	DO ind_voc= 1, emt_att%nvoc
406
407	IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN
408	len_index = len_index + 1
409	len_index_voc = len_index_voc + 1
410	ENDIF
411
412	END DO
413	ENDIF
414
415	!
416	!-- PMs: There is one input species name for all PM
417	!-- This variable has 3 dimensions, one for PM1, PM2.5 and PM10
418
419	IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" ) THEN
420
421	IF ( TRIM( spc_names(ind_mod) ) == "PM1" ) THEN
422	len_index = len_index + 1
423	ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" ) THEN
424	len_index = len_index + 1
425	ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" ) THEN
426	len_index = len_index + 1
427	ENDIF
428
429	ENDIF
430
431	!
432	!-- NOX: NO2 and NO
433
434	IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOX" ) THEN
435
436	IF ( TRIM( spc_names(ind_mod) ) == "NO" ) THEN
437	len_index = len_index + 1
438	ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" ) THEN
439	len_index = len_index + 1
440	ENDIF
441
442	ENDIF
443
444	!
445	!-- SOX: SO2 and SO4
446
447	IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOX" ) THEN
448
449	IF ( TRIM( spc_names(ind_mod) ) == "SO2" ) THEN
450	len_index = len_index + 1
451	ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" ) THEN
452	len_index = len_index + 1
453	ENDIF
454
455	ENDIF
456
457	!
458	!-- Other Species
459
460	IF ( TRIM( emt_att%species_name(ind_inp) ) == &
461	TRIM( spc_names(ind_mod) ) ) THEN
462	len_index = len_index + 1
463	ENDIF
464
465	END DO ! ind_inp ...
466
467	END DO ! ind_mod ...
468
469
470	!
471	!-- Allocate arrays
472
473	IF ( len_index > 0 ) THEN
474
475	ALLOCATE ( match_spec_input(len_index) )
476	ALLOCATE ( match_spec_model(len_index) )
477
478	IF ( len_index_voc > 0 ) THEN
479
480	!
481	!-- Contains indices of the VOC model species
482
483	ALLOCATE( match_spec_voc_model(len_index_voc) )
484
485	!
486	!-- Contains the indices of different values of VOC composition
487	!-- of input variable VOC_composition
488
489	ALLOCATE( match_spec_voc_input(len_index_voc) )
490
491	ENDIF
492
493	!
494	!-- Pass the species indices to declared arrays
495
496	len_index = 0
497	len_index_voc = 0
498
499	DO ind_mod = 1, nvar
500	DO ind_inp = 1, nspec_emis_inp
501
502	!
503	!-- VOCs
504
505	IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND. &
506	ALLOCATED (match_spec_voc_input) ) THEN
507
508	DO ind_voc = 1, emt_att%nvoc
509
510	IF ( TRIM( emt_att%voc_name(ind_voc) ) == &
511	TRIM( spc_names(ind_mod) ) ) THEN
512
513	len_index = len_index + 1
514	len_index_voc = len_index_voc + 1
515
516	match_spec_input(len_index) = ind_inp
517	match_spec_model(len_index) = ind_mod
518
519	match_spec_voc_input(len_index_voc) = ind_voc
520	match_spec_voc_model(len_index_voc) = ind_mod
521
522	ENDIF
523
524	END DO
525
526	ENDIF
527
528	!
529	!-- PMs
530
531	IF ( TRIM( emt_att%species_name(ind_inp) ) == "PM" ) THEN
532
533	IF ( TRIM( spc_names(ind_mod) ) == "PM1" ) THEN
534	len_index = len_index + 1
535	match_spec_input(len_index) = ind_inp
536	match_spec_model(len_index) = ind_mod
537	ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM25" ) THEN
538	len_index = len_index + 1
539	match_spec_input(len_index) = ind_inp
540	match_spec_model(len_index) = ind_mod
541	ELSEIF ( TRIM( spc_names(ind_mod) ) == "PM10" ) THEN
542	len_index = len_index + 1
543	match_spec_input(len_index) = ind_inp
544	match_spec_model(len_index) = ind_mod
545	ENDIF
546
547	ENDIF
548
549	!
550	!-- NOX
551	IF ( TRIM( emt_att%species_name(ind_inp) ) == "NOX" ) THEN
552
553	IF ( TRIM( spc_names(ind_mod) ) == "NO" ) THEN
554	len_index = len_index + 1
555
556	match_spec_input(len_index) = ind_inp
557	match_spec_model(len_index) = ind_mod
558
559	ELSEIF ( TRIM( spc_names(ind_mod) ) == "NO2" ) THEN
560	len_index = len_index + 1
561
562	match_spec_input(len_index) = ind_inp
563	match_spec_model(len_index) = ind_mod
564
565	ENDIF
566
567	ENDIF
568
569
570	!
571	!-- SOX
572
573	IF ( TRIM( emt_att%species_name(ind_inp) ) == "SOX" ) THEN
574
575	IF ( TRIM( spc_names(ind_mod) ) == "SO2" ) THEN
576	len_index = len_index + 1
577	match_spec_input(len_index) = ind_inp
578	match_spec_model(len_index) = ind_mod
579	ELSEIF ( TRIM( spc_names(ind_mod) ) == "SO4" ) THEN
580	len_index = len_index + 1
581	match_spec_input(len_index) = ind_inp
582	match_spec_model(len_index) = ind_mod
583	ENDIF
584
585	ENDIF
586
587	!
588	!-- Other Species
589
590	IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN
591	len_index = len_index + 1
592	match_spec_input(len_index) = ind_inp
593	match_spec_model(len_index) = ind_mod
594	ENDIF
595
596	END DO ! inp_ind
597
598	END DO ! inp_mod
599
600	!
601	!-- Error reporting (no matching)
602
603	ELSE
604
605	message_string = 'None of given Emission Species matches' // &
606	' model chemical species' // &
607	' Emission routine is not called'
608	CALL message( 'chem_emissions_matching', 'CM0440', 0, 0, 0, 6, 0 )
609
610	ENDIF
611
612	!
613	!-- Error reporting (no species)
614
615	ELSE
616
617	message_string = 'Array of Emission species not allocated: ' // &
618	' Either no emission species are provided as input or' // &
619	' no chemical species are used by PALM:' // &
620	' Emission routine is not called'
621	CALL message( 'chem_emissions_matching', 'CM0441', 0, 2, 0, 6, 0 )
622
623	ENDIF
624
625	!
626	!-- LOD 2 (PRE-PROCESSED mode)
627
628	CASE (2)
629
630	len_index = 0
631	len_index_voc = 0
632
633	IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN
634	!
635	!-- Cycle over model species
636	DO ind_mod = 1, nvar
637
638	!
639	!-- Cycle over input species
640	DO ind_inp = 1, nspec_emis_inp
641
642	!
643	!-- Check for VOC Species
644
645	IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN
646	DO ind_voc = 1, emt_att%nvoc
647	IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN
648	len_index = len_index + 1
649	len_index_voc = len_index_voc + 1
650	ENDIF
651	END DO
652	ENDIF
653
654	!
655	!-- Other Species
656
657	IF ( TRIM(emt_att%species_name(ind_inp)) == TRIM(spc_names(ind_mod)) ) THEN
658	len_index = len_index + 1
659	ENDIF
660	ENDDO
661	ENDDO
662
663	!
664	!-- Allocate array for storing the indices of the matched species
665
666	IF ( len_index > 0 ) THEN
667
668	ALLOCATE ( match_spec_input(len_index) )
669
670	ALLOCATE ( match_spec_model(len_index) )
671
672	IF ( len_index_voc > 0 ) THEN
673	!
674	!-- contains indices of the VOC model species
675	ALLOCATE( match_spec_voc_model(len_index_voc) )
676	!
677	!-- contains the indices of different values of VOC composition of input variable VOC_composition
678	ALLOCATE( match_spec_voc_input(len_index_voc) )
679
680	ENDIF
681
682	!
683	!-- pass the species indices to declared arrays
684
685	len_index = 0
686
687	!
688	!-- Cycle over model species
689
690	DO ind_mod = 1, nvar
691
692	!
693	!-- Cycle over Input species
694
695	DO ind_inp = 1, nspec_emis_inp
696
697	!
698	!-- VOCs
699
700	IF ( TRIM(emt_att%species_name(ind_inp) ) == "VOC" .AND. &
701	ALLOCATED(match_spec_voc_input) ) THEN
702
703	DO ind_voc= 1, emt_att%nvoc
704	IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN
705	len_index = len_index + 1
706	len_index_voc = len_index_voc + 1
707
708	match_spec_input(len_index) = ind_inp
709	match_spec_model(len_index) = ind_mod
710
711	match_spec_voc_input(len_index_voc) = ind_voc
712	match_spec_voc_model(len_index_voc) = ind_mod
713	ENDIF
714	END DO
715	ENDIF
716
717	!
718	!-- Other Species
719
720	IF ( TRIM( emt_att%species_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN
721	len_index = len_index + 1
722	match_spec_input(len_index) = ind_inp
723	match_spec_model(len_index) = ind_mod
724	ENDIF
725
726	END DO ! ind_inp
727	END DO ! ind_mod
728
729	ELSE ! if len_index_voc <= 0
730
731	!
732	!-- in case there are no species matching (just informational message)
733
734	message_string = 'Non of given emission species' // &
735	' matches' // &
736	' model chemical species:' // &
737	' Emission routine is not called'
738	CALL message( 'chem_emissions_matching', 'CM0438', 0, 0, 0, 6, 0 )
739	ENDIF
740
741	!
742	!-- Error check (no matching)
743
744	ELSE
745
746	!
747	!-- either spc_names is zero or nspec_emis_inp is not allocated
748	message_string = 'Array of Emission species not allocated:' // &
749	' Either no emission species are provided as input or' // &
750	' no chemical species are used by PALM:' // &
751	' Emission routine is not called'
752	CALL message( 'chem_emissions_matching', 'CM0439', 0, 2, 0, 6, 0 )
753
754	ENDIF
755
756	!
757	!-- If emission module is switched on but mode_emis is not specified or it is given the wrong name
758
759	!
760	!-- Error check (no species)
761
762	CASE DEFAULT
763
764	message_string = 'Emission Module switched ON, but' // &
765	' either no emission mode specified or incorrectly given :' // &
766	' please, pass the correct value to the namelist parameter "mode_emis"'
767	CALL message( 'chem_emissions_matching', 'CM0445', 2, 2, 0, 6, 0 )
768
769	END SELECT
770
771	IF ( debug_output ) CALL debug_message( 'chem_emissions_match', 'end' )
772
773	END SUBROUTINE chem_emissions_match
774
775
776	!------------------------------------------------------------------------------!
777	! Description:
778	! ------------
779	!> Initialization:
780	!> Netcdf reading, arrays allocation and first calculation of cssws
781	!> fluxes at timestep 0
782	!------------------------------------------------------------------------------!
783
784	SUBROUTINE chem_emissions_init
785
786	USE netcdf_data_input_mod, &
787	ONLY: chem_emis, chem_emis_att
788
789	IMPLICIT NONE
790
791	INTEGER(iwp) :: ispec !< running index
792
793	!
794	!-- Actions for initial runs
795	! IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN
796	!-- ...
797	!
798	!
799	!-- Actions for restart runs
800	! ELSE
801	!-- ...
802	!
803	! ENDIF
804
805
806	IF ( debug_output ) CALL debug_message( 'chem_emissions_init', 'start' )
807
808	!
809	!-- Matching
810
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
825	ALLOCATE( chem_emis_att%xm(n_matched_vars) )
826
827	DO ispec = 1, n_matched_vars
828	SELECT CASE ( TRIM( spc_names(match_spec_model(ispec)) ) )
829	CASE ( 'SO2' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 2
830	CASE ( 'SO4' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 4
831	CASE ( 'NO' ); chem_emis_att%xm(ispec) = xm_N + xm_O
832	CASE ( 'NO2' ); chem_emis_att%xm(ispec) = xm_N + xm_O * 2
833	CASE ( 'NH3' ); chem_emis_att%xm(ispec) = xm_N + xm_H * 3
834	CASE ( 'CO' ); chem_emis_att%xm(ispec) = xm_C + xm_O
835	CASE ( 'CO2' ); chem_emis_att%xm(ispec) = xm_C + xm_O * 2
836	CASE ( 'CH4' ); chem_emis_att%xm(ispec) = xm_C + xm_H * 4
837	CASE ( 'HNO3' ); chem_emis_att%xm(ispec) = xm_H + xm_N + xm_O*3
838	CASE DEFAULT
839	chem_emis_att%xm(ispec) = 1.0_wp
840	END SELECT
841	ENDDO
842
843
844	!
845	!-- Get emissions for the first time step base on LOD (if defined)
846	!-- or emission mode (if no LOD defined)
847
848	!
849	!-- NOTE - I could use a combined if ( lod = xxx .or. mode = 'XXX' )
850	!-- type of decision structure but I think it is much better
851	!-- to implement it this way (i.e., conditional on lod if it
852	!-- is defined, and mode if not) as we can easily take out
853	!-- the case structure for mode_emis later on.
854
855	IF ( emiss_lod < 0 ) THEN !-- no LOD defined (not likely)
856
857	SELECT CASE ( TRIM( mode_emis ) )
858
859	CASE ( 'PARAMETERIZED' ) ! LOD 0
860
861	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
862	ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) )
863	ENDIF
864
865	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars)
866
867	CASE ( 'DEFAULT' ) ! LOD 1
868
869	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
870	ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) )
871	ENDIF
872
873	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars )
874
875	CASE ( 'PRE-PROCESSED' ) ! LOD 2
876
877	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
878	ALLOCATE( emis_distribution(nzb:nzt+1,0:ny,0:nx,n_matched_vars) )
879	ENDIF
880
881	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars )
882
883	END SELECT
884
885	ELSE ! if LOD is defined
886
887	SELECT CASE ( emiss_lod )
888
889	CASE ( 0 ) ! parameterized mode
890
891	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
892	ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) )
893	ENDIF
894
895	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars)
896
897	CASE ( 1 ) ! default mode
898
899	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
900	ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) )
901	ENDIF
902
903	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars )
904
905	CASE ( 2 ) ! pre-processed mode
906
907	IF ( .NOT. ALLOCATED( emis_distribution) ) THEN
908	ALLOCATE( emis_distribution(nzb:nzt+1,0:ny,0:nx,n_matched_vars) )
909	ENDIF
910
911	CALL chem_emissions_setup( chem_emis_att, chem_emis, n_matched_vars )
912
913	END SELECT
914
915	ENDIF
916
917	!
918	! -- initialize
919
920	emis_distribution = 0.0_wp
921
922	ENDIF
923
924	IF ( debug_output ) CALL debug_message( 'chem_emissions_init', 'end' )
925
926	END SUBROUTINE chem_emissions_init
927
928
929
930	!------------------------------------------------------------------------------!
931	! Description:
932	! ------------
933	!> Routine for Update of Emission values at each timestep
934	!-------------------------------------------------------------------------------!
935
936	SUBROUTINE chem_emissions_setup( emt_att, emt, n_matched_vars )
937
938	USE surface_mod, &
939	ONLY: surf_def_h, surf_lsm_h, surf_usm_h
940
941	USE netcdf_data_input_mod, &
942	ONLY: street_type_f
943
944	USE arrays_3d, &
945	ONLY: hyp, pt
946
947
948	IMPLICIT NONE
949
950
951	TYPE(chem_emis_att_type), INTENT(INOUT) :: emt_att !< variable to store emission information
952
953	TYPE(chem_emis_val_type), INTENT(INOUT), ALLOCATABLE, DIMENSION(:) :: emt !< variable to store emission input values,
954	!< depending on the considered species
955
956	INTEGER,INTENT(IN) :: n_matched_vars !< Output of matching routine with number
957	!< of matched species
958
959	INTEGER(iwp) :: i !< running index for grid in x-direction
960	INTEGER(iwp) :: i_pm_comp !< index for number of PM components
961	INTEGER(iwp) :: icat !< Index for number of categories
962	INTEGER(iwp) :: ispec !< index for number of species
963	INTEGER(iwp) :: ivoc !< Index for number of VOCs
964	INTEGER(iwp) :: j !< running index for grid in y-direction
965	INTEGER(iwp) :: k !< running index for grid in z-direction
966	INTEGER(iwp) :: m !< running index for horizontal surfaces
967
968	!
969	!-- CONVERSION FACTORS: TIME
970	REAL(wp), PARAMETER :: hour_per_year = 8760.0_wp !< number of hours in a year of 365 days
971	REAL(wp), PARAMETER :: hour_per_day = 24.0_wp !< number of hours in a day
972	REAL(wp), PARAMETER :: s_per_hour = 3600.0_wp !< number of sec per hour (s)/(hour)
973	REAL(wp), PARAMETER :: s_per_day = 86400.0_wp !< number of sec per day (s)/(day)
974
975	REAL(wp), PARAMETER :: day_to_s = 1.0_wp/s_per_day !< conversion day -> sec
976	REAL(wp), PARAMETER :: hour_to_s = 1.0_wp/s_per_hour !< conversion hours -> sec
977	REAL(wp), PARAMETER :: year_to_s = 1.0_wp/(s_per_hour*hour_per_year) !< conversion year -> sec
978	!
979	!-- CONVERSION FACTORS: MASS
980	REAL(wp), PARAMETER :: g_to_kg = 1.0E-03_wp !< Conversion from g to kg (kg/g)
981	REAL(wp), PARAMETER :: miug_to_kg = 1.0E-09_wp !< Conversion from g to kg (kg/g)
982	REAL(wp), PARAMETER :: tons_to_kg = 100.0_wp !< Conversion from tons to kg (kg/tons)
983	!
984	!-- CONVERSION FACTORS: PPM
985	REAL(wp), PARAMETER :: ratio2ppm = 1.0E+06_wp
986
987	REAL(wp), DIMENSION(24) :: par_emis_time_factor !< time factors for the parameterized mode:
988	!< fixed houlry profile for example day
989	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: conv_to_ratio !< factor used for converting input
990	!< to concentration ratio
991	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: tmp_temp !< temporary variable for abs. temperature
992
993	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: delta_emis !< incremental emission factor
994	REAL(wp), DIMENSION(:), ALLOCATABLE :: time_factor !< factor for time scaling of emissions
995	REAL(wp), DIMENSION(:,:), ALLOCATABLE :: emis !< emission factor
996
997	IF ( emission_output_required ) THEN
998
999	!
1000	!-- Set emis_dt to be used - since chemistry ODEs can be stiff, the option
1001	!-- to solve them at every RK substep is present to help improve stability
1002	!-- should the need arises
1003
1004	IF ( call_chem_at_all_substeps ) THEN
1005
1006	dt_emis = dt_3d * weight_pres(intermediate_timestep_count)
1007
1008	ELSE
1009
1010	dt_emis = dt_3d
1011
1012	ENDIF
1013
1014	!
1015	!-- Conversion of units to the ones employed in PALM
1016	!-- In PARAMETERIZED mode no conversion is performed: in this case input units are fixed
1017
1018	IF ( TRIM( mode_emis ) == "DEFAULT" .OR. TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN
1019
1020	SELECT CASE ( TRIM( emt_att%units ) )
1021
1022	CASE ( 'kg/m2/s', 'KG/M2/S' ); conversion_factor = 1.0_wp ! kg
1023	CASE ( 'kg/m2/hour', 'KG/M2/HOUR' ); conversion_factor = hour_to_s
1024	CASE ( 'kg/m2/day', 'KG/M2/DAY' ); conversion_factor = day_to_s
1025	CASE ( 'kg/m2/year', 'KG/M2/YEAR' ); conversion_factor = year_to_s
1026
1027	CASE ( 'ton/m2/s', 'TON/M2/S' ); conversion_factor = tons_to_kg ! tonnes
1028	CASE ( 'ton/m2/hour', 'TON/M2/HOUR' ); conversion_factor = tons_to_kg*hour_to_s
1029	CASE ( 'ton/m2/year', 'TON/M2/YEAR' ); conversion_factor = tons_to_kg*year_to_s
1030
1031	CASE ( 'g/m2/s', 'G/M2/S' ); conversion_factor = g_to_kg ! grams
1032	CASE ( 'g/m2/hour', 'G/M2/HOUR' ); conversion_factor = g_to_kg*hour_to_s
1033	CASE ( 'g/m2/year', 'G/M2/YEAR' ); conversion_factor = g_to_kg*year_to_s
1034
1035	CASE ( 'micrograms/m2/s', 'MICROGRAMS/M2/S' ); conversion_factor = miug_to_kg ! ug
1036	CASE ( 'micrograms/m2/hour', 'MICROGRAMS/M2/HOUR' ); conversion_factor = miug_to_kg*hour_to_s
1037	CASE ( 'micrograms/m2/year', 'MICROGRAMS/M2/YEAR' ); conversion_factor = miug_to_kg*year_to_s
1038
1039	!
1040	!-- Error check (need units)
1041
1042	CASE DEFAULT
1043	message_string = 'The Units of the provided emission input' // &
1044	' are not the ones required by PALM-4U: please check ' // &
1045	' emission module documentation.'
1046	CALL message( 'chem_emissions_setup', 'CM0446', 2, 2, 0, 6, 0 )
1047
1048	END SELECT
1049
1050	ENDIF
1051
1052	!
1053	!-- Conversion factor to convert kg/m**2/s to ppm/s
1054
1055	DO i = nxl, nxr
1056	DO j = nys, nyn
1057
1058	!
1059	!-- Derive Temperature from Potential Temperature
1060
1061	tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * &
1062	( hyp(nzb:nzt+1) / p_0 )**rd_d_cp
1063
1064	!
1065	!-- We need to pass to cssws <- (ppm/s) * dz
1066	!-- Input is Nmole/(m^2*s)
1067	!-- To go to ppmdz multiply the input by (m2/N)dz
1068	!-- (m*2/N)dz == V/N
1069	!-- V/N = RT/P
1070
1071	conv_to_ratio(nzb:nzt+1,j,i) = rgas_univ * & ! J K-1 mol-1
1072	tmp_temp(nzb:nzt+1,j,i) / & ! K
1073	hyp(nzb:nzt+1) ! Pa
1074
1075	! (ecc) for reference
1076	! m*3/Nmole (J/mol)K^-1 K Pa
1077	! conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) )
1078
1079	ENDDO
1080	ENDDO
1081
1082
1083	! (ecc) moved initialization immediately after allocation
1084	!
1085	!-- Initialize
1086
1087	! emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp
1088
1089
1090	!
1091	!-- LOD 2 (PRE-PROCESSED MODE)
1092
1093	IF ( emiss_lod == 2 ) THEN
1094
1095	! for reference (ecc)
1096	! IF ( TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN
1097
1098	!
1099	!-- Update time indices
1100
1101	CALL time_preprocessed_indices( index_hh )
1102
1103
1104	!
1105	!-- LOD 1 (DEFAULT MODE)
1106
1107	ELSEIF ( emiss_lod == 1 ) THEN
1108
1109	! for reference (ecc)
1110	! ELSEIF ( TRIM( mode_emis ) == "DEFAULT" ) THEN
1111
1112	!
1113	!-- Allocate array where to store temporary emission values
1114
1115	IF ( .NOT. ALLOCATED(emis) ) ALLOCATE( emis(nys:nyn,nxl:nxr) )
1116
1117	!
1118	!-- Allocate time factor per category
1119
1120	ALLOCATE( time_factor(emt_att%ncat) )
1121
1122	!
1123	!-- Read-in hourly emission time factor
1124
1125	IF ( TRIM(time_fac_type) == "HOUR" ) THEN
1126
1127	!
1128	!-- Update time indices
1129
1130	CALL time_default_indices( month_of_year, day_of_month, hour_of_day, index_hh )
1131
1132	!
1133	!-- Check if the index is less or equal to the temporal dimension of HOURLY emission files
1134
1135	IF ( index_hh <= SIZE( emt_att%hourly_emis_time_factor(1,:) ) ) THEN
1136
1137	!
1138	!-- Read-in the correspondant time factor
1139
1140	time_factor(:) = emt_att%hourly_emis_time_factor(:,index_hh)
1141
1142	!
1143	!-- Error check (time out of range)
1144
1145	ELSE
1146
1147	message_string = 'The "HOUR" time-factors in the DEFAULT mode ' // &
1148	' are not provided for each hour of the total simulation time'
1149	CALL message( 'chem_emissions_setup', 'CM0448', 2, 2, 0, 6, 0 )
1150
1151	ENDIF
1152
1153	!
1154	!-- Read-in MDH emissions time factors
1155
1156	ELSEIF ( TRIM( time_fac_type ) == "MDH" ) THEN
1157
1158	!
1159	!-- Update time indices
1160	CALL time_default_indices( daytype_mdh, month_of_year, day_of_month, &
1161	hour_of_day, index_mm, index_dd,index_hh )
1162
1163	!
1164	!-- Check if the index is less or equal to the temporal dimension of MDH emission files
1165
1166	IF ( ( index_hh + index_dd + index_mm) <= SIZE( emt_att%mdh_emis_time_factor(1,:) ) ) THEN
1167
1168	!
1169	!-- Read in corresponding time factor
1170
1171	time_factor(:) = emt_att%mdh_emis_time_factor(:,index_mm) * &
1172	emt_att%mdh_emis_time_factor(:,index_dd) * &
1173	emt_att%mdh_emis_time_factor(:,index_hh)
1174
1175	!
1176	!-- Error check (MDH time factor not provided)
1177
1178	ELSE
1179
1180	message_string = 'The "MDH" time-factors in the DEFAULT mode ' // &
1181	' are not provided for each hour/day/month of the total simulation time'
1182	CALL message( 'chem_emissions_setup', 'CM0449', 2, 2, 0, 6, 0 )
1183
1184	ENDIF
1185
1186	!
1187	!-- Error check (no time factor defined)
1188
1189	ELSE
1190
1191	message_string = 'In the DEFAULT mode the time factor' // &
1192	' has to be defined in the NAMELIST'
1193	CALL message( 'chem_emissions_setup', 'CM0450', 2, 2, 0, 6, 0 )
1194
1195	ENDIF
1196
1197	!
1198	!-- PARAMETERIZED MODE
1199
1200	ELSEIF ( emiss_lod == 0 ) THEN
1201
1202
1203	! for reference (ecc)
1204	! ELSEIF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN
1205
1206	!
1207	!-- assign constant values of time factors, diurnal profile for traffic sector
1208
1209	par_emis_time_factor( : ) = (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, &
1210	0.056, 0.053, 0.051, 0.051, 0.052, 0.055, &
1211	0.059, 0.061, 0.064, 0.067, 0.069, 0.069, &
1212	0.049, 0.039, 0.039, 0.029, 0.024, 0.019 /)
1213
1214	IF ( .NOT. ALLOCATED (time_factor) ) ALLOCATE (time_factor(1))
1215
1216	!
1217	!-- Get time-factor for specific hour
1218
1219	index_hh = hour_of_day
1220	time_factor(1) = par_emis_time_factor(index_hh)
1221
1222	ENDIF ! emiss_lod
1223
1224
1225	!
1226	!-- Emission distribution calculation
1227
1228	!
1229	!-- LOD 0 (PARAMETERIZED mode)
1230
1231	IF ( emiss_lod == 0 ) THEN
1232
1233	! for reference (ecc)
1234	! IF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN
1235
1236	DO ispec = 1, n_matched_vars
1237
1238	!
1239	!-- Units are micromoles/m*2day (or kilograms/m*2day for PMs)
1240
1241	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1242	surface_csflux(match_spec_input(ispec)) * &
1243	time_factor(1) * hour_to_s
1244
1245	ENDDO
1246
1247
1248	!
1249	!-- LOD 1 (DEFAULT mode)
1250
1251
1252	ELSEIF ( emiss_lod == 1 ) THEN
1253
1254	! for referene (ecc)
1255	! ELSEIF ( TRIM( mode_emis ) == "DEFAULT" ) THEN
1256
1257	!
1258	!-- Allocate array for the emission value corresponding to a specific category and time factor
1259
1260	ALLOCATE (delta_emis(nys:nyn,nxl:nxr))
1261
1262	!
1263	!-- Cycle over categories
1264
1265	DO icat = 1, emt_att%ncat
1266
1267	!
1268	!-- Cycle over Species: n_matched_vars represents the number of species
1269	!-- in common between the emission input data and the chemistry mechanism used
1270
1271	DO ispec = 1, n_matched_vars
1272
1273	emis(nys:nyn,nxl:nxr) = &
1274	emt(match_spec_input(ispec))% &
1275	default_emission_data(icat,nys+1:nyn+1,nxl+1:nxr+1)
1276
1277	!
1278	!-- NO
1279
1280	IF ( TRIM(spc_names(match_spec_model(ispec))) == "NO" ) THEN
1281
1282	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1283	time_factor(icat) * &
1284	emt_att%nox_comp(icat,1) * &
1285	conversion_factor * hour_per_day
1286
1287	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1288	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1289	delta_emis(nys:nyn,nxl:nxr)
1290	!
1291	!-- NO2
1292
1293	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "NO2" ) THEN
1294
1295	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1296	time_factor(icat) * &
1297	emt_att%nox_comp(icat,2) * &
1298	conversion_factor * hour_per_day
1299
1300	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1301	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1302	delta_emis(nys:nyn,nxl:nxr)
1303
1304	!
1305	!-- SO2
1306	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "SO2" ) THEN
1307
1308	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1309	time_factor(icat) * &
1310	emt_att%sox_comp(icat,1) * &
1311	conversion_factor * hour_per_day
1312
1313	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1314	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1315	delta_emis(nys:nyn,nxl:nxr)
1316
1317	!
1318	!-- SO4
1319
1320	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "SO4" ) THEN
1321
1322
1323	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1324	time_factor(icat) * &
1325	emt_att%sox_comp(icat,2) * &
1326	conversion_factor * hour_per_day
1327
1328	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1329	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1330	delta_emis(nys:nyn,nxl:nxr)
1331
1332
1333	!
1334	!-- PM1
1335
1336	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" ) THEN
1337
1338	DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,1) ) ! cycle through components
1339
1340	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1341	time_factor(icat) * &
1342	emt_att%pm_comp(icat,i_pm_comp,1) * &
1343	conversion_factor * hour_per_day
1344
1345	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1346	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1347	delta_emis(nys:nyn,nxl:nxr)
1348
1349	ENDDO
1350
1351	!
1352	!-- PM2.5
1353
1354	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "PM25" ) THEN
1355
1356	DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,2) ) ! cycle through components
1357
1358	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1359	time_factor(icat) * &
1360	emt_att%pm_comp(icat,i_pm_comp,2) * &
1361	conversion_factor * hour_per_day
1362
1363	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1364	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1365	delta_emis(nys:nyn,nxl:nxr)
1366
1367	ENDDO
1368
1369	!
1370	!-- PM10
1371
1372	ELSEIF ( TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN
1373
1374	DO i_pm_comp = 1, SIZE( emt_att%pm_comp(1,:,3) ) ! cycle through components
1375
1376	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & ! kg/m2/s
1377	time_factor(icat) * &
1378	emt_att%pm_comp(icat,i_pm_comp,3) * &
1379	conversion_factor * hour_per_day
1380
1381	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1382	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1383	delta_emis(nys:nyn,nxl:nxr)
1384
1385	ENDDO
1386
1387	!
1388	!-- VOCs
1389
1390	ELSEIF ( SIZE( match_spec_voc_input ) > 0 ) THEN
1391
1392	DO ivoc = 1, SIZE( match_spec_voc_input ) ! cycle through components
1393
1394	IF ( TRIM(spc_names(match_spec_model(ispec))) == &
1395	TRIM(emt_att%voc_name(ivoc)) ) THEN
1396
1397	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * &
1398	time_factor(icat) * &
1399	emt_att%voc_comp(icat,match_spec_voc_input(ivoc)) * &
1400	conversion_factor * hour_per_day
1401
1402	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1403	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1404	delta_emis(nys:nyn,nxl:nxr)
1405
1406	ENDIF
1407
1408	ENDDO
1409
1410	!
1411	!-- any other species
1412
1413	ELSE
1414
1415	delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * &
1416	time_factor(icat) * &
1417	conversion_factor * hour_per_day
1418
1419	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1420	emis_distribution(1,nys:nyn,nxl:nxr,ispec) + &
1421	delta_emis(nys:nyn,nxl:nxr)
1422
1423	ENDIF ! TRIM spc_names
1424
1425	emis = 0
1426
1427	ENDDO
1428
1429	delta_emis = 0
1430
1431	ENDDO
1432
1433	!
1434	!-- LOD 2 (PRE-PROCESSED mode)
1435
1436	ELSEIF ( emiss_lod == 2 ) THEN
1437
1438	! for reference (ecc)
1439	! ELSEIF ( TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN
1440
1441	!
1442	!-- Cycle over species: n_matched_vars represents the number of species
1443	!-- in common between the emission input data and the chemistry mechanism used
1444
1445	DO ispec = 1, n_matched_vars
1446
1447	! (ecc)
1448	emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1449	emt(match_spec_input(ispec))% &
1450	preproc_emission_data(index_hh,1,nys+1:nyn+1,nxl+1:nxr+1) * &
1451	conversion_factor
1452
1453
1454	! emis_distribution(1,nys:nyn,nxl:nxr,ispec) = &
1455	! emt(match_spec_input(ispec))% &
1456	! preproc_emission_data(index_hh,1,:,:) * &
1457	! conversion_factor
1458	ENDDO
1459
1460	ENDIF ! emiss_lod
1461
1462
1463	!
1464	!-- Cycle to transform x,y coordinates to the one of surface_mod and to assign emission values to cssws
1465
1466	!
1467	!-- LOD 0 (PARAMETERIZED mode)
1468	!-- Units of inputs are micromoles/m2/s
1469
1470	IF ( emiss_lod == 0 ) THEN
1471	! for reference (ecc)
1472	! IF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN
1473
1474	IF (street_type_f%from_file) THEN
1475
1476	!
1477	!-- Streets are lsm surfaces, hence, no usm surface treatment required.
1478	!-- However, urban surface may be initialized via default initialization
1479	!-- in surface_mod, e.g. at horizontal urban walls that are at k == 0
1480	!-- (building is lower than the first grid point). Hence, in order to
1481	!-- have only emissions at streets, set the surfaces emissions to zero
1482	!-- at urban walls.
1483
1484	IF ( surf_usm_h%ns >=1 ) surf_usm_h%cssws = 0.0_wp
1485
1486	!
1487	!-- Treat land-surfaces.
1488
1489	DO m = 1, surf_lsm_h%ns
1490
1491	i = surf_lsm_h%i(m)
1492	j = surf_lsm_h%j(m)
1493	k = surf_lsm_h%k(m)
1494
1495	!
1496	!-- set everything to zero then reassign according to street type
1497
1498	surf_lsm_h%cssws(:,m) = 0.0_wp
1499
1500	IF ( street_type_f%var(j,i) >= main_street_id .AND. &
1501	street_type_f%var(j,i) < max_street_id ) THEN
1502
1503	!
1504	!-- Cycle over matched species
1505
1506	DO ispec = 1, n_matched_vars
1507
1508	!
1509	!-- PMs are already in kilograms
1510
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^2s) 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
1526	ELSE
1527
1528	!
1529	!-- ppm/s m kg/m^3
1530	surf_lsm_h%cssws(match_spec_model(ispec),m) = &
1531	emiss_factor_main(match_spec_input(ispec)) * &
1532	emis_distribution(1,j,i,ispec) * & ! micromoles/(m^2*s)
1533	conv_to_ratio(k,j,i) * & ! m^3/Nmole
1534	rho_air(k) ! kg/m^3
1535
1536	ENDIF
1537
1538	ENDDO ! ispec
1539
1540
1541	ELSEIF ( street_type_f%var(j,i) >= side_street_id .AND. &
1542	street_type_f%var(j,i) < main_street_id ) THEN
1543
1544	!
1545	!-- Cycle over matched species
1546
1547	DO ispec = 1, n_matched_vars
1548
1549	!
1550	!-- PMs are already in kilograms
1551
1552	IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. &
1553	TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. &
1554	TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN
1555
1556	!
1557	!-- kg/(m^2s) kg/m^3
1558	surf_lsm_h%cssws(match_spec_model(ispec),m) = &
1559	emiss_factor_side(match_spec_input(ispec)) * &
1560	emis_distribution(1,j,i,ispec) * & ! kg/(m^2*s)
1561	rho_air(k) ! kg/m^3
1562	!
1563	!-- Other species
1564	!-- Inputs are micromoles
1565
1566	ELSE
1567
1568	!
1569	!-- ppm/s m kg/m^3
1570
1571	surf_lsm_h%cssws(match_spec_model(ispec),m) = &
1572	emiss_factor_side(match_spec_input(ispec)) * &
1573	emis_distribution(1,j,i,ispec) * & ! micromoles/(m^2*s)
1574	conv_to_ratio(k,j,i) * & ! m^3/Nmole
1575	rho_air(k) ! kg/m^3
1576
1577	ENDIF
1578
1579	ENDDO ! ispec
1580
1581	!
1582	!-- If no street type is defined, then assign zero emission to all the species
1583
1584	! (ecc) moved to front (for reference)
1585	! ELSE
1586	!
1587	! surf_lsm_h%cssws(:,m) = 0.0_wp
1588
1589	ENDIF ! street type
1590
1591	ENDDO ! m
1592
1593	ENDIF ! street_type_f%from_file
1594
1595
1596	!
1597	!-- LOD 1 (DEFAULT) and LOD 2 (PRE-PROCESSED)
1598
1599
1600	ELSE
1601
1602
1603	DO ispec = 1, n_matched_vars
1604
1605	!
1606	!-- Default surfaces
1607
1608	DO m = 1, surf_def_h(0)%ns
1609
1610	i = surf_def_h(0)%i(m)
1611	j = surf_def_h(0)%j(m)
1612
1613	IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN
1614
1615	!
1616	!-- PMs
1617	IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. &
1618	TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. &
1619	TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN
1620
1621	surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3
1622	emis_distribution(1,j,i,ispec)* & ! kg/m2/s
1623	rho_air(nzb) ! kg/m^3
1624
1625	ELSE
1626
1627	!
1628	!-- VOCs
1629	IF ( len_index_voc > 0 .AND. &
1630	emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN
1631
1632	surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1633	emis_distribution(1,j,i,ispec) * & ! mole/m2/s
1634	conv_to_ratio(nzb,j,i) * & ! m^3/mole
1635	ratio2ppm * & ! ppm
1636	rho_air(nzb) ! kg/m^3
1637
1638
1639	!
1640	!-- Other species
1641
1642	ELSE
1643
1644	surf_def_h(0)%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1645	emis_distribution(1,j,i,ispec) * & ! kg/m2/s
1646	( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg
1647	conv_to_ratio(nzb,j,i) * & ! m^3/mole
1648	ratio2ppm * & ! ppm
1649	rho_air(nzb) ! kg/m^3
1650
1651	ENDIF ! VOC
1652
1653	ENDIF ! PM
1654
1655	ENDIF ! emis_distribution > 0
1656
1657	ENDDO ! m
1658
1659	!
1660	!-- LSM surfaces
1661
1662	DO m = 1, surf_lsm_h%ns
1663
1664	i = surf_lsm_h%i(m)
1665	j = surf_lsm_h%j(m)
1666	k = surf_lsm_h%k(m)
1667
1668	IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN
1669
1670	!
1671	!-- PMs
1672	IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. &
1673	TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. &
1674	TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN
1675
1676	surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3
1677	emis_distribution(1,j,i,ispec) * & ! kg/m2/s
1678	rho_air(k) ! kg/m^3
1679
1680	ELSE
1681
1682	!
1683	!-- VOCs
1684
1685	IF ( len_index_voc > 0 .AND. &
1686	emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN
1687
1688	surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1689	emis_distribution(1,j,i,ispec) * & ! mole/m2/s
1690	conv_to_ratio(k,j,i) * & ! m^3/mole
1691	ratio2ppm * & ! ppm
1692	rho_air(k) ! kg/m^3
1693
1694	!
1695	!-- Other species
1696
1697	ELSE
1698
1699	surf_lsm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1700	emis_distribution(1,j,i,ispec) * & ! kg/m2/s
1701	( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg
1702	conv_to_ratio(k,j,i) * & ! m^3/mole
1703	ratio2ppm * & ! ppm
1704	rho_air(k) ! kg/m^3
1705
1706	ENDIF ! VOC
1707
1708	ENDIF ! PM
1709
1710	ENDIF ! emis_distribution
1711
1712	ENDDO ! m
1713
1714	!
1715	!-- USM surfaces
1716
1717	DO m = 1, surf_usm_h%ns
1718
1719	i = surf_usm_h%i(m)
1720	j = surf_usm_h%j(m)
1721	k = surf_usm_h%k(m)
1722
1723	IF ( emis_distribution(1,j,i,ispec) > 0.0_wp ) THEN
1724
1725	!
1726	!-- PMs
1727	IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. &
1728	TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. &
1729	TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN
1730
1731	surf_usm_h%cssws(match_spec_model(ispec),m) = & ! kg/m2/s * kg/m3
1732	emis_distribution(1,j,i,ispec)* & ! kg/m2/s
1733	rho_air(k) ! kg/m^3
1734
1735	ELSE
1736
1737	!
1738	!-- VOCs
1739	IF ( len_index_voc > 0 .AND. &
1740	emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN
1741
1742	surf_usm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1743	emis_distribution(1,j,i,ispec) * & ! m2/s
1744	conv_to_ratio(k,j,i) * & ! m^3/mole
1745	ratio2ppm * & ! ppm
1746	rho_air(k) ! kg/m^3
1747
1748	!
1749	!-- Other species
1750	ELSE
1751
1752	surf_usm_h%cssws(match_spec_model(ispec),m) = & ! ppm/s * m * kg/m3
1753	emis_distribution(1,j,i,ispec) * & ! kg/m2/s
1754	( 1.0_wp / emt_att%xm(ispec) ) * & ! mole/kg
1755	conv_to_ratio(k,j,i) * & ! m^3/mole
1756	ratio2ppm* & ! ppm
1757	rho_air(k) ! kg/m^3
1758
1759
1760	ENDIF ! VOC
1761
1762	ENDIF ! PM
1763
1764	ENDIF ! emis_distribution
1765
1766	ENDDO ! m
1767
1768	ENDDO
1769
1770	ENDIF
1771
1772	!
1773	!-- Deallocate time_factor in case of DEFAULT mode)
1774
1775	IF ( ALLOCATED (time_factor) ) DEALLOCATE (time_factor)
1776
1777	ENDIF
1778
1779	END SUBROUTINE chem_emissions_setup
1780
1781	END MODULE chem_emissions_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |