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