Changeset 4055 for palm/trunk/SOURCE/chem_emissions_mod.f90
- Timestamp:
- Jun 27, 2019 9:47:29 AM (5 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/chem_emissions_mod.f90 (modified) (60 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/chem_emissions_mod.f90
r3968 r4055 27 27 ! ----------------- 28 28 ! $Id$ 29 ! - replaced local thermo. constants w/ module definitions in 30 ! basic_constants_and_equations_mod (rgas_univ, p_0, r_d_cp) 31 ! - initialize array emis_distribution immediately following allocation 32 ! - lots of minor formatting changes based on review sesson in 20190325 33 ! (E.C. Chan) 34 ! 35 ! 3968 2019-05-13 11:04:01Z suehring 29 36 ! - in subroutine chem_emissions_match replace all decision structures relating to 30 37 ! mode_emis to emiss_lod … … 54 61 ! removed unnecessary informative messsages/location 55 62 ! messages and corrected comments on PM units from to kg 56 ! bug fix: spcs_name replaced by nvar in DO loops63 ! bug fix: spcs_name replaced by nvar in DO loops 57 64 ! 58 65 ! 3591 2018-11-30 17:37:32Z suehring … … 118 125 MODULE chem_emissions_mod 119 126 120 USE arrays_3d, 127 USE arrays_3d, & 121 128 ONLY: rho_air 122 129 123 USE control_parameters, & 124 ONLY: debug_output, & 125 end_time, message_string, initializing_actions, & 130 USE basic_constants_and_equations_mod, & 131 ONLY: rgas_univ, p_0, rd_d_cp 132 133 USE control_parameters, & 134 ONLY: debug_output, & 135 end_time, message_string, initializing_actions, & 126 136 intermediate_timestep_count, dt_3d 127 137 … … 134 144 #endif 135 145 136 USE netcdf_data_input_mod, 146 USE netcdf_data_input_mod, & 137 147 ONLY: chem_emis_att_type, chem_emis_val_type 138 148 139 USE date_and_time_mod, 140 ONLY: day_of_month, hour_of_day, 141 index_mm, index_dd, index_hh, 142 month_of_year, hour_of_day, 149 USE date_and_time_mod, & 150 ONLY: day_of_month, hour_of_day, & 151 index_mm, index_dd, index_hh, & 152 month_of_year, hour_of_day, & 143 153 time_default_indices, time_preprocessed_indices 144 154 145 USE chem_gasphase_mod, 155 USE chem_gasphase_mod, & 146 156 ONLY: nvar, spc_names 147 157 148 158 USE chem_modules 149 159 150 USE statistics, 160 USE statistics, & 151 161 ONLY: weight_pres 152 162 … … 159 169 CHARACTER (LEN=80) :: filename_emis !< Variable for the name of the netcdf input file 160 170 171 INTEGER(iwp) :: dt_emis !< Time Step Emissions 161 172 INTEGER(iwp) :: i !< index 1st selected dimension (some dims are not spatial) 162 173 INTEGER(iwp) :: j !< index 2nd selected dimension … … 165 176 INTEGER(iwp) :: j_start !< Index to start read variable from netcdf in additional dims 166 177 INTEGER(iwp) :: j_end !< Index to end read variable from netcdf in additional dims 178 INTEGER(iwp) :: len_index !< length of index (used for several indices) 179 INTEGER(iwp) :: len_index_pm !< length of PMs index 180 INTEGER(iwp) :: len_index_voc !< length of voc index 167 181 INTEGER(iwp) :: z_start !< Index to start read variable from netcdf in additional dims 168 182 INTEGER(iwp) :: z_end !< Index to end read variable from netcdf in additional dims 169 INTEGER(iwp) :: dt_emis !< Time Step Emissions 170 INTEGER(iwp) :: len_index !< length of index (used for several indices) 171 INTEGER(iwp) :: len_index_voc !< length of voc index 172 INTEGER(iwp) :: len_index_pm !< length of PMs index 173 174 REAL(wp) :: con_factor !< Units Conversion Factor 175 176 REAL(wp), PARAMETER :: Rgas = 8.3144 !< gas constant in J/mol/K 177 REAL(wp), PARAMETER :: pref_i = 1.0_wp / 100000.0_wp !< Inverse Reference Pressure (1/Pa) 178 REAL(wp), PARAMETER :: r_cp = 0.286_wp !< R / cp (exponent for potential temperature) 183 184 REAL(wp) :: conversion_factor !< Units Conversion Factor 179 185 180 186 SAVE … … 204 210 ! 205 211 !-- Public Variables 206 PUBLIC con _factor, len_index, len_index_pm, len_index_voc212 PUBLIC conversion_factor, len_index, len_index_pm, len_index_voc 207 213 208 214 CONTAINS … … 249 255 SUBROUTINE chem_emissions_match( emt_att,len_index ) 250 256 251 252 253 257 INTEGER(iwp) :: ind_inp !< Parameters for cycling through chemical input species 254 258 INTEGER(iwp) :: ind_mod !< Parameters for cycling through chemical model species … … 282 286 len_index = 0 283 287 284 IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN 288 ! number of species and number of matched species can be different 289 ! but call is only made if both are greater than zero 290 291 IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN 285 292 286 293 ! 287 294 !-- Cycle over model species 288 295 289 DO ind_mod = 1, nvar296 DO ind_mod = 1, nvar 290 297 ind_inp = 1 291 DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) !< 'novalue' is the default298 DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) !< 'novalue' is the default 292 299 IF ( TRIM( surface_csflux_name(ind_inp) ) == TRIM( spc_names(ind_mod) ) ) THEN 293 300 len_index = len_index + 1 … … 310 317 len_index = 0 311 318 312 DO ind_mod = 1, nvar319 DO ind_mod = 1, nvar 313 320 ind_inp = 1 314 DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' )321 DO WHILE ( TRIM( surface_csflux_name(ind_inp) ) /= 'novalue' ) 315 322 IF ( TRIM(surface_csflux_name(ind_inp)) == & 316 323 TRIM(spc_names(ind_mod)) ) THEN … … 326 333 !-- Check 327 334 328 DO ispec = 1, len_index329 330 IF ( emiss_factor_main(match_spec_input(ispec) ) < 0 .AND.&331 emiss_factor_side(match_spec_input(ispec) ) < 0 ) THEN335 DO ispec = 1, len_index 336 337 IF ( emiss_factor_main(match_spec_input(ispec) ) < 0 .AND. & 338 emiss_factor_side(match_spec_input(ispec) ) < 0 ) THEN 332 339 333 340 message_string = 'PARAMETERIZED emissions mode selected:' // & … … 370 377 CASE (1) 371 378 372 len_index = 0 ! TOTAL number of species (to be accumulated) 373 len_index_voc = 0 ! TOTAL number of VOCs (to be accumulated) 374 len_index_pm = 3 ! TOTAL number of PMs: PM1, PM2.5, PM10. 375 376 IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN 379 len_index = 0 ! total number of species (to be accumulated) 380 len_index_voc = 0 ! total number of VOCs (to be accumulated) 381 len_index_pm = 3 ! total number of PMs: PM1, PM2.5, PM10. 382 383 ! 384 !-- number of model species and input species could be different 385 !-- but process this only when both are non-zero 386 387 IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN 377 388 378 389 ! 379 390 !-- Cycle over model species 380 DO ind_mod = 1, nvar391 DO ind_mod = 1, nvar 381 392 382 393 ! 383 394 !-- Cycle over input species 384 395 385 DO ind_inp = 1, nspec_emis_inp396 DO ind_inp = 1, nspec_emis_inp 386 397 387 398 ! … … 389 400 390 401 IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN 391 DO ind_voc= 1, emt_att%nvoc402 DO ind_voc= 1, emt_att%nvoc 392 403 393 404 IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN … … 483 494 len_index_voc = 0 484 495 485 DO ind_mod = 1, nvar486 DO ind_inp = 1, nspec_emis_inp496 DO ind_mod = 1, nvar 497 DO ind_inp = 1, nspec_emis_inp 487 498 488 499 ! 489 500 !-- VOCs 490 501 491 IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND.&492 ALLOCATED (match_spec_voc_input) ) THEN493 494 DO ind_voc = 1, emt_att%nvoc502 IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" .AND. & 503 ALLOCATED (match_spec_voc_input) ) THEN 504 505 DO ind_voc = 1, emt_att%nvoc 495 506 496 507 IF ( TRIM( emt_att%voc_name(ind_voc) ) == & … … 617 628 len_index_voc = 0 618 629 619 IF ( nvar > 0 .AND. (nspec_emis_inp > 0)) THEN630 IF ( nvar > 0 .AND. nspec_emis_inp > 0 ) THEN 620 631 ! 621 632 !-- Cycle over model species 622 DO ind_mod = 1, nvar633 DO ind_mod = 1, nvar 623 634 624 635 ! 625 636 !-- Cycle over input species 626 DO ind_inp = 1, nspec_emis_inp637 DO ind_inp = 1, nspec_emis_inp 627 638 628 639 ! … … 630 641 631 642 IF ( TRIM( emt_att%species_name(ind_inp) ) == "VOC" ) THEN 632 DO ind_voc = 1, emt_att%nvoc643 DO ind_voc = 1, emt_att%nvoc 633 644 IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN 634 645 len_index = len_index + 1 … … 674 685 !-- Cycle over model species 675 686 676 DO ind_mod = 1, nvar687 DO ind_mod = 1, nvar 677 688 678 689 ! 679 690 !-- Cycle over Input species 680 691 681 DO ind_inp = 1, nspec_emis_inp692 DO ind_inp = 1, nspec_emis_inp 682 693 683 694 ! 684 695 !-- VOCs 685 696 686 IF ( TRIM(emt_att%species_name(ind_inp) ) == "VOC" .AND.&697 IF ( TRIM(emt_att%species_name(ind_inp) ) == "VOC" .AND. & 687 698 ALLOCATED(match_spec_voc_input) ) THEN 688 699 689 DO ind_voc= 1, emt_att%nvoc700 DO ind_voc= 1, emt_att%nvoc 690 701 IF ( TRIM( emt_att%voc_name(ind_voc) ) == TRIM( spc_names(ind_mod) ) ) THEN 691 702 len_index = len_index + 1 … … 716 727 717 728 ! 718 !-- in case there are no species matching 729 !-- in case there are no species matching (just informational message) 719 730 720 731 message_string = 'Non of given emission species' // & … … 767 778 !> fluxes at timestep 0 768 779 !------------------------------------------------------------------------------! 780 769 781 SUBROUTINE chem_emissions_init 770 782 771 USE netcdf_data_input_mod, 783 USE netcdf_data_input_mod, & 772 784 ONLY: chem_emis, chem_emis_att 773 785 774 786 IMPLICIT NONE 775 787 776 INTEGER(iwp) :: ispec!< running index788 INTEGER(iwp) :: ispec !< running index 777 789 778 790 ! … … 810 822 ALLOCATE( chem_emis_att%xm(n_matched_vars) ) 811 823 812 DO ispec = 1, n_matched_vars824 DO ispec = 1, n_matched_vars 813 825 SELECT CASE ( TRIM( spc_names(match_spec_model(ispec)) ) ) 814 826 CASE ( 'SO2' ); chem_emis_att%xm(ispec) = xm_S + xm_O * 2 … … 837 849 !-- is defined, and mode if not) as we can easily take out 838 850 !-- the case structure for mode_emis later on. 839 !-- (ecc 20140424)840 851 841 852 IF ( emiss_lod < 0 ) THEN !-- no LOD defined (not likely) … … 845 856 CASE ( 'PARAMETERIZED' ) ! LOD 0 846 857 847 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN858 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 848 859 ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) ) 849 860 ENDIF … … 853 864 CASE ( 'DEFAULT' ) ! LOD 1 854 865 855 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN866 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 856 867 ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) ) 857 868 ENDIF … … 861 872 CASE ( 'PRE-PROCESSED' ) ! LOD 2 862 873 863 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN874 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 864 875 ALLOCATE( emis_distribution(nzb:nzt+1,0:ny,0:nx,n_matched_vars) ) 865 876 ENDIF … … 875 886 CASE ( 0 ) ! parameterized mode 876 887 877 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN888 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 878 889 ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) ) 879 890 ENDIF … … 883 894 CASE ( 1 ) ! default mode 884 895 885 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN896 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 886 897 ALLOCATE( emis_distribution(1,0:ny,0:nx,n_matched_vars) ) 887 898 ENDIF … … 891 902 CASE ( 2 ) ! pre-processed mode 892 903 893 IF ( .NOT.ALLOCATED( emis_distribution) ) THEN904 IF ( .NOT. ALLOCATED( emis_distribution) ) THEN 894 905 ALLOCATE( emis_distribution(nzb:nzt+1,0:ny,0:nx,n_matched_vars) ) 895 906 ENDIF … … 900 911 901 912 ENDIF 913 914 ! 915 ! -- initialize 916 917 emis_distribution = 0.0_wp 902 918 903 919 ENDIF … … 917 933 SUBROUTINE chem_emissions_setup( emt_att, emt, n_matched_vars ) 918 934 919 USE surface_mod, 935 USE surface_mod, & 920 936 ONLY: surf_def_h, surf_lsm_h, surf_usm_h 921 USE netcdf_data_input_mod, & 937 938 USE netcdf_data_input_mod, & 922 939 ONLY: street_type_f 923 USE arrays_3d, & 940 941 USE arrays_3d, & 924 942 ONLY: hyp, pt 925 943 … … 937 955 938 956 INTEGER(iwp) :: i !< running index for grid in x-direction 957 INTEGER(iwp) :: i_pm_comp !< index for number of PM components 958 INTEGER(iwp) :: icat !< Index for number of categories 959 INTEGER(iwp) :: ispec !< index for number of species 960 INTEGER(iwp) :: ivoc !< Index for number of VOCs 939 961 INTEGER(iwp) :: j !< running index for grid in y-direction 940 962 INTEGER(iwp) :: k !< running index for grid in z-direction 941 963 INTEGER(iwp) :: m !< running index for horizontal surfaces 942 943 INTEGER(iwp) :: icat !< Index for number of categories944 INTEGER(iwp) :: ispec !< index for number of species945 INTEGER(iwp) :: i_pm_comp !< index for number of PM components946 INTEGER(iwp) :: ivoc !< Index for number of VOCs947 948 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: delta_emis949 REAL(wp), ALLOCATABLE, DIMENSION(:) :: time_factor !< factor for time scaling of emissions950 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: emis951 952 REAL(wp), DIMENSION(24) :: par_emis_time_factor !< time factors for the parameterized mode:953 !< fixed houlry profile for example day954 REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: conv_to_ratio !< factor used for converting input955 !< to concentration ratio956 REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: tmp_temp957 964 958 965 ! 959 966 !-- CONVERSION FACTORS: TIME 960 REAL(wp), PARAMETER :: s_per_hour = 3600.0 !< number of sec per hour (s)/(hour)961 REAL(wp), PARAMETER :: s_per_day = 86400.0 !< number of sec per day (s)/(day)962 REAL(wp), PARAMETER :: hour_per_year = 8760.0 !< number of hours in a year of 365 days963 REAL(wp), PARAMETER :: hour_per_day = 24.0 !< number of hours in a day964 965 REAL(wp), PARAMETER :: hour_to_s = 1/s_per_hour !< conversion from hours to seconds (s/hour) ~ 0.2777778966 REAL(wp), PARAMETER :: day_to_s = 1/s_per_day !< conversion from day to seconds (s/day) ~ 1.157407e-05967 REAL(wp), PARAMETER :: year_to_s = 1/(s_per_hour*hour_per_year) !< conversion from year to sec (s/year) ~ 3.170979e-08967 REAL(wp), PARAMETER :: hour_per_year = 8760.0_wp !< number of hours in a year of 365 days 968 REAL(wp), PARAMETER :: hour_per_day = 24.0_wp !< number of hours in a day 969 REAL(wp), PARAMETER :: s_per_hour = 3600.0_wp !< number of sec per hour (s)/(hour) 970 REAL(wp), PARAMETER :: s_per_day = 86400.0_wp !< number of sec per day (s)/(day) 971 972 REAL(wp), PARAMETER :: day_to_s = 1.0_wp/s_per_day !< conversion day -> sec 973 REAL(wp), PARAMETER :: hour_to_s = 1.0_wp/s_per_hour !< conversion hours -> sec 974 REAL(wp), PARAMETER :: year_to_s = 1.0_wp/(s_per_hour*hour_per_year) !< conversion year -> sec 968 975 ! 969 !-- CONVERSION FACTORS: WEIGHT970 REAL(wp), PARAMETER :: tons_to_kg = 100 !< Conversion from tons to kg (kg/tons)971 REAL(wp), PARAMETER :: g_to_kg = 0.001 !< Conversion from g to kg (kg/g)972 REAL(wp), PARAMETER :: miug_to_kg = 0.000000001 !< Conversion from g to kg (kg/g)976 !-- CONVERSION FACTORS: MASS 977 REAL(wp), PARAMETER :: g_to_kg = 1.0E-03_wp !< Conversion from g to kg (kg/g) 978 REAL(wp), PARAMETER :: miug_to_kg = 1.0E-09_wp !< Conversion from g to kg (kg/g) 979 REAL(wp), PARAMETER :: tons_to_kg = 100.0_wp !< Conversion from tons to kg (kg/tons) 973 980 ! 974 !-- CONVERSION FACTORS: fraction to ppm 975 REAL(wp), PARAMETER :: ratio2ppm = 1.0e06 976 !------------------------------------------------------ 981 !-- CONVERSION FACTORS: PPM 982 REAL(wp), PARAMETER :: ratio2ppm = 1.0E+06_wp 983 984 REAL(wp), DIMENSION(24) :: par_emis_time_factor !< time factors for the parameterized mode: 985 !< fixed houlry profile for example day 986 REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: conv_to_ratio !< factor used for converting input 987 !< to concentration ratio 988 REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: tmp_temp !< temporary variable for abs. temperature 989 990 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: delta_emis !< incremental emission factor 991 REAL(wp), DIMENSION(:), ALLOCATABLE :: time_factor !< factor for time scaling of emissions 992 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: emis !< emission factor 977 993 978 994 IF ( emission_output_required ) THEN 979 995 980 ! 981 !-- Set emis_dt 996 ! 997 !-- Set emis_dt to be used - since chemistry ODEs can be stiff, the option 998 !-- to solve them at every RK substep is present to help improve stability 999 !-- should the need arises 1000 982 1001 IF ( call_chem_at_all_substeps ) THEN 983 1002 … … 994 1013 !-- In PARAMETERIZED mode no conversion is performed: in this case input units are fixed 995 1014 996 IF ( TRIM( mode_emis ) == "DEFAULT" .OR.TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN1015 IF ( TRIM( mode_emis ) == "DEFAULT" .OR. TRIM( mode_emis ) == "PRE-PROCESSED" ) THEN 997 1016 998 1017 SELECT CASE ( TRIM( emt_att%units ) ) 999 1018 1000 CASE ( 'kg/m2/s', 'KG/M2/S' ); con_factor = 1.0_wp ! kg1001 CASE ( 'kg/m2/hour', 'KG/M2/HOUR' ); con _factor = hour_to_s1002 CASE ( 'kg/m2/day', 'KG/M2/DAY' ); con_factor = day_to_s1003 CASE ( 'kg/m2/year', 'KG/M2/YEAR' ); con _factor = year_to_s1004 1005 CASE ( 'ton/m2/s', 'TON/M2/S' ); con_factor = tons_to_kg ! tonnes1006 CASE ( 'ton/m2/hour', 'TON/M2/HOUR' ); con _factor = tons_to_kg*hour_to_s1007 CASE ( 'ton/m2/year', 'TON/M2/YEAR' ); con _factor = tons_to_kg*year_to_s1008 1009 CASE ( 'g/m2/s', 'G/M2/S' ); con_factor = g_to_kg ! grams1010 CASE ( 'g/m2/hour', 'G/M2/HOUR' ); con _factor = g_to_kg*hour_to_s1011 CASE ( 'g/m2/year', 'G/M2/YEAR' ); con _factor = g_to_kg*year_to_s1012 1013 CASE ( 'micrograms/m2/s', 'MICROGRAMS/M2/S' ); con_factor = miug_to_kg ! ug1014 CASE ( 'micrograms/m2/hour', 'MICROGRAMS/M2/HOUR' ); con _factor = miug_to_kg*hour_to_s1015 CASE ( 'micrograms/m2/year', 'MICROGRAMS/M2/YEAR' ); con _factor = miug_to_kg*year_to_s1019 CASE ( 'kg/m2/s', 'KG/M2/S' ); conversion_factor = 1.0_wp ! kg 1020 CASE ( 'kg/m2/hour', 'KG/M2/HOUR' ); conversion_factor = hour_to_s 1021 CASE ( 'kg/m2/day', 'KG/M2/DAY' ); conversion_factor = day_to_s 1022 CASE ( 'kg/m2/year', 'KG/M2/YEAR' ); conversion_factor = year_to_s 1023 1024 CASE ( 'ton/m2/s', 'TON/M2/S' ); conversion_factor = tons_to_kg ! tonnes 1025 CASE ( 'ton/m2/hour', 'TON/M2/HOUR' ); conversion_factor = tons_to_kg*hour_to_s 1026 CASE ( 'ton/m2/year', 'TON/M2/YEAR' ); conversion_factor = tons_to_kg*year_to_s 1027 1028 CASE ( 'g/m2/s', 'G/M2/S' ); conversion_factor = g_to_kg ! grams 1029 CASE ( 'g/m2/hour', 'G/M2/HOUR' ); conversion_factor = g_to_kg*hour_to_s 1030 CASE ( 'g/m2/year', 'G/M2/YEAR' ); conversion_factor = g_to_kg*year_to_s 1031 1032 CASE ( 'micrograms/m2/s', 'MICROGRAMS/M2/S' ); conversion_factor = miug_to_kg ! ug 1033 CASE ( 'micrograms/m2/hour', 'MICROGRAMS/M2/HOUR' ); conversion_factor = miug_to_kg*hour_to_s 1034 CASE ( 'micrograms/m2/year', 'MICROGRAMS/M2/YEAR' ); conversion_factor = miug_to_kg*year_to_s 1016 1035 1017 1036 ! … … 1031 1050 !-- Conversion factor to convert kg/m**2/s to ppm/s 1032 1051 1033 DO i = nxl, nxr1034 DO j = nys, nyn1052 DO i = nxl, nxr 1053 DO j = nys, nyn 1035 1054 1036 1055 ! 1037 1056 !-- Derive Temperature from Potential Temperature 1038 1057 1039 tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * ( hyp(nzb:nzt+1) * pref_i )**r_cp 1058 tmp_temp(nzb:nzt+1,j,i) = pt(nzb:nzt+1,j,i) * & 1059 ( hyp(nzb:nzt+1) / p_0 )**rd_d_cp 1040 1060 1041 1061 ! … … 1046 1066 !-- V/N = RT/P 1047 1067 1068 conv_to_ratio(nzb:nzt+1,j,i) = rgas_univ * & ! J K-1 mol-1 1069 tmp_temp(nzb:nzt+1,j,i) / & ! K 1070 hyp(nzb:nzt+1) ! Pa 1071 1072 ! (ecc) for reference 1048 1073 ! m**3/Nmole (J/mol)*K^-1 K Pa 1049 conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) )1074 ! conv_to_ratio(nzb:nzt+1,j,i) = ( (Rgas * tmp_temp(nzb:nzt+1,j,i)) / ((hyp(nzb:nzt+1))) ) 1050 1075 1051 1076 ENDDO … … 1053 1078 1054 1079 1080 ! (ecc) moved initialization immediately after allocation 1055 1081 ! 1056 1082 !-- Initialize 1057 1083 1058 emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp1084 ! emis_distribution(:,nys:nyn,nxl:nxr,:) = 0.0_wp 1059 1085 1060 1086 … … 1084 1110 !-- Allocate array where to store temporary emission values 1085 1111 1086 IF ( .NOT.ALLOCATED(emis) ) ALLOCATE( emis(nys:nyn,nxl:nxr) )1112 IF ( .NOT. ALLOCATED(emis) ) ALLOCATE( emis(nys:nyn,nxl:nxr) ) 1087 1113 1088 1114 ! … … 1129 1155 ! 1130 1156 !-- Update time indices 1131 CALL time_default_indices( daytype_mdh, month_of_year, day_of_month, 1157 CALL time_default_indices( daytype_mdh, month_of_year, day_of_month, & 1132 1158 hour_of_day, index_mm, index_dd,index_hh ) 1133 1159 … … 1178 1204 !-- assign constant values of time factors, diurnal profile for traffic sector 1179 1205 1180 par_emis_time_factor( : ) = (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, 1181 0.056, 0.053, 0.051, 0.051, 0.052, 0.055, 1182 0.059, 0.061, 0.064, 0.067, 0.069, 0.069, 1206 par_emis_time_factor( : ) = (/ 0.009, 0.004, 0.004, 0.009, 0.029, 0.039, & 1207 0.056, 0.053, 0.051, 0.051, 0.052, 0.055, & 1208 0.059, 0.061, 0.064, 0.067, 0.069, 0.069, & 1183 1209 0.049, 0.039, 0.039, 0.029, 0.024, 0.019 /) 1184 1210 1185 IF ( .NOT.ALLOCATED (time_factor) ) ALLOCATE (time_factor(1))1211 IF ( .NOT. ALLOCATED (time_factor) ) ALLOCATE (time_factor(1)) 1186 1212 1187 1213 ! … … 1205 1231 ! IF ( TRIM( mode_emis ) == "PARAMETERIZED" ) THEN 1206 1232 1207 DO ispec = 1, n_matched_vars1233 DO ispec = 1, n_matched_vars 1208 1234 1209 1235 ! … … 1214 1240 time_factor(1) * hour_to_s 1215 1241 1216 ENDDO 1242 ENDDO 1217 1243 1218 1244 … … 1234 1260 !-- Cycle over categories 1235 1261 1236 DO icat = 1, emt_att%ncat1262 DO icat = 1, emt_att%ncat 1237 1263 1238 1264 ! … … 1240 1266 !-- in common between the emission input data and the chemistry mechanism used 1241 1267 1242 DO ispec = 1, n_matched_vars1268 DO ispec = 1, n_matched_vars 1243 1269 1244 1270 emis(nys:nyn,nxl:nxr) = & … … 1254 1280 time_factor(icat) * & 1255 1281 emt_att%nox_comp(icat,1) * & 1256 con _factor * hour_per_day1282 conversion_factor * hour_per_day 1257 1283 1258 1284 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1267 1293 time_factor(icat) * & 1268 1294 emt_att%nox_comp(icat,2) * & 1269 con _factor * hour_per_day1295 conversion_factor * hour_per_day 1270 1296 1271 1297 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1280 1306 time_factor(icat) * & 1281 1307 emt_att%sox_comp(icat,1) * & 1282 con _factor * hour_per_day1308 conversion_factor * hour_per_day 1283 1309 1284 1310 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1295 1321 time_factor(icat) * & 1296 1322 emt_att%sox_comp(icat,2) * & 1297 con _factor * hour_per_day1323 conversion_factor * hour_per_day 1298 1324 1299 1325 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1312 1338 time_factor(icat) * & 1313 1339 emt_att%pm_comp(icat,i_pm_comp,1) * & 1314 con _factor * hour_per_day1340 conversion_factor * hour_per_day 1315 1341 1316 1342 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1330 1356 time_factor(icat) * & 1331 1357 emt_att%pm_comp(icat,i_pm_comp,2) * & 1332 con _factor * hour_per_day1358 conversion_factor * hour_per_day 1333 1359 1334 1360 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1348 1374 time_factor(icat) * & 1349 1375 emt_att%pm_comp(icat,i_pm_comp,3) * & 1350 con _factor * hour_per_day1376 conversion_factor * hour_per_day 1351 1377 1352 1378 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1369 1395 time_factor(icat) * & 1370 1396 emt_att%voc_comp(icat,match_spec_voc_input(ivoc)) * & 1371 con _factor * hour_per_day1397 conversion_factor * hour_per_day 1372 1398 1373 1399 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1386 1412 delta_emis(nys:nyn,nxl:nxr) = emis(nys:nyn,nxl:nxr) * & 1387 1413 time_factor(icat) * & 1388 con _factor * hour_per_day1414 conversion_factor * hour_per_day 1389 1415 1390 1416 emis_distribution(1,nys:nyn,nxl:nxr,ispec) = & … … 1394 1420 ENDIF ! TRIM spc_names 1395 1421 1396 emis (:,:)= 01422 emis = 0 1397 1423 1398 1424 ENDDO 1399 1425 1400 delta_emis (:,:)=01426 delta_emis = 0 1401 1427 1402 1428 ENDDO … … 1414 1440 !-- in common between the emission input data and the chemistry mechanism used 1415 1441 1416 DO ispec = 1, n_matched_vars1442 DO ispec = 1, n_matched_vars 1417 1443 1418 1444 ! (ecc) … … 1420 1446 emt(match_spec_input(ispec))% & 1421 1447 preproc_emission_data(index_hh,1,nys+1:nyn+1,nxl+1:nxr+1) * & 1422 con _factor1448 conversion_factor 1423 1449 1424 1450 … … 1426 1452 ! emt(match_spec_input(ispec))% & 1427 1453 ! preproc_emission_data(index_hh,1,:,:) * & 1428 ! con _factor1454 ! conversion_factor 1429 1455 ENDDO 1430 1456 … … 1464 1490 k = surf_lsm_h%k(m) 1465 1491 1492 ! 1493 !-- set everything to zero then reassign according to street type 1494 1495 surf_lsm_h%cssws(:,m) = 0.0_wp 1496 1466 1497 IF ( street_type_f%var(j,i) >= main_street_id .AND. & 1467 1498 street_type_f%var(j,i) < max_street_id ) THEN … … 1475 1506 !-- PMs are already in kilograms 1476 1507 1477 IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. &1478 TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. &1508 IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. & 1509 TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. & 1479 1510 TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN 1480 1511 … … 1516 1547 !-- PMs are already in kilograms 1517 1548 1518 IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR.&1519 TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR.&1549 IF ( TRIM(spc_names(match_spec_model(ispec))) == "PM1" .OR. & 1550 TRIM(spc_names(match_spec_model(ispec))) == "PM25" .OR. & 1520 1551 TRIM(spc_names(match_spec_model(ispec))) == "PM10" ) THEN 1521 1552 … … 1548 1579 !-- If no street type is defined, then assign zero emission to all the species 1549 1580 1550 ELSE 1551 1552 surf_lsm_h%cssws(:,m) = 0.0_wp 1581 ! (ecc) moved to front (for reference) 1582 ! ELSE 1583 ! 1584 ! surf_lsm_h%cssws(:,m) = 0.0_wp 1553 1585 1554 1586 ENDIF ! street type … … 1566 1598 1567 1599 1568 DO ispec = 1, n_matched_vars1600 DO ispec = 1, n_matched_vars 1569 1601 1570 1602 ! … … 1702 1734 ! 1703 1735 !-- VOCs 1704 IF ( len_index_voc > 0 1736 IF ( len_index_voc > 0 .AND. & 1705 1737 emt_att%species_name(match_spec_input(ispec)) == "VOC" ) THEN 1706 1738
Note: See TracChangeset
for help on using the changeset viewer.
