PartMC  2.4.0
gas_data.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2012 Nicole Riemer and Matthew West
2 ! Licensed under the GNU General Public License version 2 or (at your
3 ! option) any later version. See the file COPYING for details.
4 
5 !> \file
6 !> The pmc_gas_data module.
7 
8 !> The gas_data_t structure and associated subroutines.
10 
11  use pmc_spec_file
12  use pmc_mpi
13  use pmc_util
14  use pmc_netcdf
15 #ifdef PMC_USE_MPI
16  use mpi
17 #endif
18 
19  !> Maximum length of the name of a gas.
20  integer, parameter :: gas_name_len = 100
21 
22  !> Constant gas data.
23  !!
24  !! Each gas species is identified by an integer \c i between 1 and
25  !! \c gas_data_n_spec(gas_data). Species \c i has name \c gas_data%%name(i).
26  !! The variable gas data describing the current mixing ratios is stored
27  !! in the gas_state_t structure, so the mixing ratio of species \c i
28  !! is gas_state%%mix_rat(i).
30  !> Species name [length \c gas_data_n_spec(gas_data)].
31  character(len=GAS_NAME_LEN), allocatable :: name(:)
32  !> Index of the corresponding MOSAIC species [length \c
33  !> gas_data_n_spec(gas_data)]. \c to_mosaic(i) is the mosaic index of
34  !> species \c i, or 0 if there is no match.
35  integer, allocatable :: mosaic_index(:)
36  end type gas_data_t
37 
38 contains
39 
40 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
41 
42  !> Return the number of gas species.
43  elemental integer function gas_data_n_spec(gas_data)
44 
45  !> Aero data structure.
46  type(gas_data_t), intent(in) :: gas_data
47 
48  if (allocated(gas_data%name)) then
49  gas_data_n_spec = size(gas_data%name)
50  else
51  gas_data_n_spec = 0
52  end if
53 
54  end function gas_data_n_spec
55 
56 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
57 
58  !> Returns the number of the species in gas with the given name, or
59  !> returns 0 if there is no such species.
60  integer function gas_data_spec_by_name(gas_data, name)
61 
62  !> Gas data.
63  type(gas_data_t), intent(in) :: gas_data
64  !> Name of species to find.
65  character(len=*), intent(in) :: name
66 
67  integer i
68  logical found
69 
70  found = .false.
71  do i = 1,gas_data_n_spec(gas_data)
72  if (name == gas_data%name(i)) then
73  found = .true.
74  exit
75  end if
76  end do
77  if (found) then
79  else
81  end if
82 
83  end function gas_data_spec_by_name
84 
85 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
86 
87  !> Fills in gas_data%%mosaic_index.
88  subroutine gas_data_set_mosaic_map(gas_data)
89 
90  !> Gas data.
91  type(gas_data_t), intent(inout) :: gas_data
92 
93  integer, parameter :: n_mosaic_species = 77
94  character(len=GAS_NAME_LEN), parameter, dimension(n_mosaic_species) &
95  :: mosaic_species = [ &
96  "H2SO4 ", "HNO3 ", "HCl ", &
97  "NH3 ", "NO ", "NO2 ", &
98  "NO3 ", "N2O5 ", "HONO ", &
99  "HNO4 ", "O3 ", "O1D ", &
100  "O3P ", "OH ", "HO2 ", &
101  "H2O2 ", "CO ", "SO2 ", &
102  "CH4 ", "C2H6 ", "CH3O2 ", &
103  "ETHP ", "HCHO ", "CH3OH ", &
104  "ANOL ", "CH3OOH ", "ETHOOH ", &
105  "ALD2 ", "HCOOH ", "RCOOH ", &
106  "C2O3 ", "PAN ", "ARO1 ", &
107  "ARO2 ", "ALK1 ", "OLE1 ", &
108  "API1 ", "API2 ", "LIM1 ", &
109  "LIM2 ", "PAR ", "AONE ", &
110  "MGLY ", "ETH ", "OLET ", &
111  "OLEI ", "TOL ", "XYL ", &
112  "CRES ", "TO2 ", "CRO ", &
113  "OPEN ", "ONIT ", "ROOH ", &
114  "RO2 ", "ANO2 ", "NAP ", &
115  "XO2 ", "XPAR ", "ISOP ", &
116  "ISOPRD ", "ISOPP ", "ISOPN ", &
117  "ISOPO2 ", "API ", "LIM ", &
118  "DMS ", "MSA ", "DMSO ", &
119  "DMSO2 ", "CH3SO2H ", "CH3SCH2OO ", &
120  "CH3SO2 ", "CH3SO3 ", "CH3SO2OO ", &
121  "CH3SO2CH2OO ", "SULFHOX "]
122 
123  integer spec, mosaic_spec, i
124 
125  gas_data%mosaic_index = 0
126  do spec = 1,gas_data_n_spec(gas_data)
127  mosaic_spec = 0
128  do i = 1,n_mosaic_species
129  if (gas_data%name(spec) == mosaic_species(i)) then
130  mosaic_spec = i
131  end if
132  end do
133  if (mosaic_spec > 0) then
134  gas_data%mosaic_index(spec) = mosaic_spec
135  end if
136  end do
137 
138  end subroutine gas_data_set_mosaic_map
139 
140 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
141 
142  !> Read gas data from a .spec file.
143  subroutine spec_file_read_gas_data(file, gas_data)
144 
145  !> Spec file to read data from.
146  type(spec_file_t), intent(inout) :: file
147  !> Gas data.
148  type(gas_data_t), intent(inout) :: gas_data
149 
150  integer :: n_species, species, i
151  character(len=SPEC_LINE_MAX_VAR_LEN), allocatable :: species_name(:)
152  real(kind=dp), allocatable :: species_data(:,:)
153 
154  !> \page input_format_gas_data Input File Format: Gas Material Data
155  !!
156  !! A gas material data file must consist of one line per gas
157  !! species, with each line having the species name. This specifies
158  !! which species are to be recognized as gases. For example, a \c
159  !! gas_data file could contain:
160  !! <pre>
161  !! H2SO4
162  !! HNO3
163  !! HCl
164  !! NH3
165  !! </pre>
166  !!
167  !! See also:
168  !! - \ref spec_file_format --- the input file text format
169  !! - \ref output_format_gas_data --- the corresponding output format
170 
171  ! read the gas data from the specified file
172  call spec_file_read_real_named_array(file, 0, species_name, &
173  species_data)
174 
175  ! check the data size
176  if (size(species_data, 2) /= 0) then
177  call die_msg(614290516, 'each line in ' // trim(file%name) &
178  // ' must only contain the species name')
179  end if
180 
181  ! allocate and copy over the data
182  n_species = size(species_data, 1)
183  call ensure_string_array_size(gas_data%name, n_species)
184  do i = 1,n_species
185  gas_data%name(i) = species_name(i)(1:gas_name_len)
186  end do
187 
188  call ensure_integer_array_size(gas_data%mosaic_index, n_species)
189  call gas_data_set_mosaic_map(gas_data)
190 
191  end subroutine spec_file_read_gas_data
192 
193 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
194 
195  !> Determines the number of bytes required to pack the given value.
196  integer function pmc_mpi_pack_size_gas_data(val)
198  !> Value to pack.
199  type(gas_data_t), intent(in) :: val
200 
203  + pmc_mpi_pack_size_integer_array(val%mosaic_index)
204 
205  end function pmc_mpi_pack_size_gas_data
206 
207 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
208 
209  !> Packs the given value into the buffer, advancing position.
210  subroutine pmc_mpi_pack_gas_data(buffer, position, val)
212  !> Memory buffer.
213  character, intent(inout) :: buffer(:)
214  !> Current buffer position.
215  integer, intent(inout) :: position
216  !> Value to pack.
217  type(gas_data_t), intent(in) :: val
218 
219 #ifdef PMC_USE_MPI
220  integer :: prev_position
221 
222  prev_position = position
223  call pmc_mpi_pack_string_array(buffer, position, val%name)
224  call pmc_mpi_pack_integer_array(buffer, position, val%mosaic_index)
225  call assert(449872094, &
226  position - prev_position <= pmc_mpi_pack_size_gas_data(val))
227 #endif
228 
229  end subroutine pmc_mpi_pack_gas_data
230 
231 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
232 
233  !> Unpacks the given value from the buffer, advancing position.
234  subroutine pmc_mpi_unpack_gas_data(buffer, position, val)
236  !> Memory buffer.
237  character, intent(inout) :: buffer(:)
238  !> Current buffer position.
239  integer, intent(inout) :: position
240  !> Value to pack.
241  type(gas_data_t), intent(inout) :: val
242 
243 #ifdef PMC_USE_MPI
244  integer :: prev_position
245 
246  prev_position = position
247  call pmc_mpi_unpack_string_array(buffer, position, val%name)
248  call pmc_mpi_unpack_integer_array(buffer, position, val%mosaic_index)
249  call assert(137879163, &
250  position - prev_position <= pmc_mpi_pack_size_gas_data(val))
251 #endif
252 
253  end subroutine pmc_mpi_unpack_gas_data
254 
255 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
256 
257  !> Write the gas species dimension to the given NetCDF file if it
258  !> is not already present and in any case return the associated
259  !> dimid.
260  subroutine gas_data_netcdf_dim_gas_species(gas_data, ncid, &
261  dimid_gas_species)
263  !> Gas_data structure.
264  type(gas_data_t), intent(in) :: gas_data
265  !> NetCDF file ID, in data mode.
266  integer, intent(in) :: ncid
267  !> Dimid of the species dimension.
268  integer, intent(out) :: dimid_gas_species
269 
270  integer :: status, i_spec
271  integer :: varid_gas_species
272  integer :: gas_species_centers(gas_data_n_spec(gas_data))
273  character(len=((GAS_NAME_LEN + 2) * gas_data_n_spec(gas_data))) &
274  :: gas_species_names
275 
276  ! try to get the dimension ID
277  status = nf90_inq_dimid(ncid, "gas_species", dimid_gas_species)
278  if (status == nf90_noerr) return
279  if (status /= nf90_ebaddim) call pmc_nc_check(status)
280 
281  ! dimension not defined, so define now define it
282  call pmc_nc_check(nf90_redef(ncid))
283 
284  call pmc_nc_check(nf90_def_dim(ncid, "gas_species", &
285  gas_data_n_spec(gas_data), dimid_gas_species))
286  gas_species_names = ""
287  do i_spec = 1,gas_data_n_spec(gas_data)
288  gas_species_names((len_trim(gas_species_names) + 1):) &
289  = trim(gas_data%name(i_spec))
290  if (i_spec < gas_data_n_spec(gas_data)) then
291  gas_species_names((len_trim(gas_species_names) + 1):) = ","
292  end if
293  end do
294  call pmc_nc_check(nf90_def_var(ncid, "gas_species", nf90_int, &
295  dimid_gas_species, varid_gas_species))
296  call pmc_nc_check(nf90_put_att(ncid, varid_gas_species, "names", &
297  gas_species_names))
298  call pmc_nc_check(nf90_put_att(ncid, varid_gas_species, "description", &
299  "dummy dimension variable (no useful value) - read species names " &
300  // "as comma-separated values from the 'names' attribute"))
301 
302  call pmc_nc_check(nf90_enddef(ncid))
303 
304  do i_spec = 1,gas_data_n_spec(gas_data)
305  gas_species_centers(i_spec) = i_spec
306  end do
307  call pmc_nc_check(nf90_put_var(ncid, varid_gas_species, &
308  gas_species_centers))
309 
310  end subroutine gas_data_netcdf_dim_gas_species
311 
312 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
313 
314  !> Write full state.
315  subroutine gas_data_output_netcdf(gas_data, ncid)
316 
317  !> Gas_data to write.
318  type(gas_data_t), intent(in) :: gas_data
319  !> NetCDF file ID, in data mode.
320  integer, intent(in) :: ncid
321 
322  integer :: dimid_gas_species
323 
324  !> \page output_format_gas_data Output File Format: Gas Material Data
325  !!
326  !! The gas material data NetCDF dimensions are:
327  !! - \b gas_species: number of gas species
328  !!
329  !! The gas material data NetCDF variables are:
330  !! - \b gas_species (dim \c gas_species): dummy dimension variable
331  !! (no useful value) - read species names as comma-separated values
332  !! from the 'names' attribute
333  !! - \b gas_mosaic_index (dim \c gas_species): MOSAIC indices of
334  !! gas species
335  !!
336  !! See also:
337  !! - \ref input_format_gas_data --- the corresponding input format
338 
339  call gas_data_netcdf_dim_gas_species(gas_data, ncid, &
340  dimid_gas_species)
341 
342  call pmc_nc_write_integer_1d(ncid, gas_data%mosaic_index, &
343  "gas_mosaic_index", (/ dimid_gas_species /), &
344  long_name="MOSAIC indices of gas species")
345 
346  end subroutine gas_data_output_netcdf
347 
348 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
349 
350  !> Read full state.
351  subroutine gas_data_input_netcdf(gas_data, ncid)
353  !> Gas_data to read.
354  type(gas_data_t), intent(inout) :: gas_data
355  !> NetCDF file ID, in data mode.
356  integer, intent(in) :: ncid
357 
358  character(len=1000) :: name
359  integer :: dimid_gas_species, n_spec, varid_gas_species, i_spec, i
360  character(len=((GAS_NAME_LEN + 2) * 1000)) :: gas_species_names
361 
362  call pmc_nc_check(nf90_inq_dimid(ncid, "gas_species", dimid_gas_species))
363  call pmc_nc_check(nf90_inquire_dimension(ncid, dimid_gas_species, name, &
364  n_spec))
365  call assert(719237193, n_spec < 1000)
366 
367  if (allocated(gas_data%name)) deallocate(gas_data%name)
368  if (allocated(gas_data%mosaic_index)) deallocate(gas_data%mosaic_index)
369  allocate(gas_data%name(n_spec))
370  allocate(gas_data%mosaic_index(n_spec))
371 
372  call pmc_nc_read_integer_1d(ncid, gas_data%mosaic_index, &
373  "gas_mosaic_index")
374 
375  call pmc_nc_check(nf90_inq_varid(ncid, "gas_species", varid_gas_species))
376  call pmc_nc_check(nf90_get_att(ncid, varid_gas_species, "names", &
377  gas_species_names))
378  ! gas_species_names are comma-separated, so unpack them
379  do i_spec = 1,gas_data_n_spec(gas_data)
380  i = 1
381  do while ((gas_species_names(i:i) /= " ") &
382  .and. (gas_species_names(i:i) /= ","))
383  i = i + 1
384  end do
385  call assert(173021381, i > 1)
386  gas_data%name(i_spec) = gas_species_names(1:(i-1))
387  gas_species_names = gas_species_names((i+1):)
388  end do
389  call assert(469721220, gas_species_names == "")
390 
391  end subroutine gas_data_input_netcdf
392 
393 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
394 
395 end module pmc_gas_data
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:688
An input file with extra data for printing messages.
Definition: spec_file.F90:59
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:447
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
Definition: netcdf.F90:11
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:946
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
Definition: util.F90:1115
integer, parameter gas_name_len
Maximum length of the name of a gas.
Definition: gas_data.F90:20
subroutine pmc_mpi_pack_string_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:752
The gas_data_t structure and associated subroutines.
Definition: gas_data.F90:9
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:103
subroutine pmc_mpi_unpack_gas_data(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: gas_data.F90:235
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
Definition: netcdf.F90:23
subroutine gas_data_input_netcdf(gas_data, ncid)
Read full state.
Definition: gas_data.F90:352
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:134
subroutine gas_data_set_mosaic_map(gas_data)
Fills in gas_data%mosaic_index.
Definition: gas_data.F90:89
subroutine ensure_string_array_size(x, n)
Allocate or reallocate the given array to ensure it is of the given size, without preserving data...
Definition: util.F90:1189
integer function gas_data_spec_by_name(gas_data, name)
Returns the number of the species in gas with the given name, or returns 0 if there is no such specie...
Definition: gas_data.F90:61
elemental integer function gas_data_n_spec(gas_data)
Return the number of gas species.
Definition: gas_data.F90:44
subroutine pmc_mpi_pack_gas_data(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: gas_data.F90:211
Reading formatted text input.
Definition: spec_file.F90:43
subroutine pmc_mpi_unpack_string_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:1010
subroutine spec_file_read_real_named_array(file, max_lines, names, vals)
Read an array of named lines with real data. All lines must have the same number of data elements...
Definition: spec_file.F90:649
Constant gas data.
Definition: gas_data.F90:29
integer function pmc_mpi_pack_size_gas_data(val)
Determines the number of bytes required to pack the given value.
Definition: gas_data.F90:197
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
Definition: netcdf.F90:260
subroutine gas_data_netcdf_dim_gas_species(gas_data, ncid, dimid_gas_species)
Write the gas species dimension to the given NetCDF file if it is not already present and in any case...
Definition: gas_data.F90:262
Common utility subroutines.
Definition: util.F90:9
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine pmc_nc_write_integer_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple integer array to a NetCDF file.
Definition: netcdf.F90:585
integer function pmc_mpi_pack_size_string_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:507