Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |