PartMC  2.4.0
gas_state.F90
Go to the documentation of this file.
1 ! Copyright (C) 2007-2012 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_state module.
7 
8 !> The gas_state_t structure and associated subroutines.
10 
11  use pmc_util
12  use pmc_spec_file
13  use pmc_gas_data
14  use pmc_env_state
15  use pmc_mpi
16  use pmc_netcdf
17 #ifdef PMC_USE_MPI
18  use mpi
19 #endif
20 
21  !> Current state of the gas mixing ratios in the system.
22  !!
23  !! The gas species are defined by the gas_data_t structure, so that
24  !! \c gas_state%%mix_rat(i) is the current mixing ratio of the gas
25  !! with name \c gas_data%%name(i), etc.
26  !!
27  !! By convention, if gas_state_is_allocated() return \c .false.,
28  !! then the gas_state is treated as zero for all operations on
29  !! it. This will be the case for new \c gas_state_t structures.
31  !> Length n_spec, mixing ratio (ppb).
32  real(kind=dp), allocatable :: mix_rat(:)
33  end type gas_state_t
34 
35 contains
36 
37 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
38 
39  !> Determine whether the \c gas_state is correctly allocated.
40  logical function gas_state_is_allocated(gas_state)
41 
42  !> Gas state to check.
43  type(gas_state_t), intent(in) :: gas_state
44 
45  gas_state_is_allocated = allocated(gas_state%mix_rat)
46 
47  end function gas_state_is_allocated
48 
49 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
50 
51  !> Sets the sizes of the gas state.
52  subroutine gas_state_set_size(gas_state, n_spec)
53 
54  !> Gas state to be allocated.
55  type(gas_state_t), intent(inout) :: gas_state
56  !> Number of species.
57  integer, intent(in) :: n_spec
58 
59  if (allocated(gas_state%mix_rat)) deallocate(gas_state%mix_rat)
60  allocate(gas_state%mix_rat(n_spec))
61  call gas_state_zero(gas_state)
62 
63  end subroutine gas_state_set_size
64 
65 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
66 
67  !> Zeros the state.
68  subroutine gas_state_zero(gas_state)
69 
70  !> Gas state.
71  type(gas_state_t), intent(inout) :: gas_state
72 
73  if (gas_state_is_allocated(gas_state)) then
74  gas_state%mix_rat = 0d0
75  end if
76 
77  end subroutine gas_state_zero
78 
79 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
80 
81  !> Scale a gas state.
82  subroutine gas_state_scale(gas_state, alpha)
83 
84  !> Existing gas state.
85  type(gas_state_t), intent(inout) :: gas_state
86  !> Scale factor.
87  real(kind=dp), intent(in) :: alpha
88 
89  if (gas_state_is_allocated(gas_state)) then
90  gas_state%mix_rat = gas_state%mix_rat * alpha
91  end if
92 
93  end subroutine gas_state_scale
94 
95 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
96 
97  !> Adds the given gas_state_delta.
98  subroutine gas_state_add(gas_state, gas_state_delta)
99 
100  !> Existing gas state.
101  type(gas_state_t), intent(inout) :: gas_state
102  !> Incremental state.
103  type(gas_state_t), intent(in) :: gas_state_delta
104 
105  if (gas_state_is_allocated(gas_state_delta)) then
106  if (gas_state_is_allocated(gas_state)) then
107  gas_state%mix_rat = gas_state%mix_rat + gas_state_delta%mix_rat
108  else
109  gas_state%mix_rat = gas_state_delta%mix_rat
110  end if
111  end if
112 
113  end subroutine gas_state_add
114 
115 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
116 
117  !> Adds the given \c gas_state_delta scaled by \c alpha.
118  !!
119  !! Does gas_state += alpha * gas_state_delta.
120  subroutine gas_state_add_scaled(gas_state, gas_state_delta, alpha)
122  !> Existing gas state.
123  type(gas_state_t), intent(inout) :: gas_state
124  !> Incremental state.
125  type(gas_state_t), intent(in) :: gas_state_delta
126  !> Scale factor.
127  real(kind=dp), intent(in) :: alpha
128 
129  if (gas_state_is_allocated(gas_state_delta)) then
130  if (gas_state_is_allocated(gas_state)) then
131  gas_state%mix_rat = gas_state%mix_rat &
132  + alpha * gas_state_delta%mix_rat
133  else
134  gas_state%mix_rat = alpha * gas_state_delta%mix_rat
135  end if
136  end if
137 
138  end subroutine gas_state_add_scaled
139 
140 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
141 
142  !> Subtracts the given gas_state_delta.
143  subroutine gas_state_sub(gas_state, gas_state_delta)
145  !> Existing gas state.
146  type(gas_state_t), intent(inout) :: gas_state
147  !> Incremental state.
148  type(gas_state_t), intent(in) :: gas_state_delta
149 
150  if (gas_state_is_allocated(gas_state_delta)) then
151  if (gas_state_is_allocated(gas_state)) then
152  gas_state%mix_rat = gas_state%mix_rat - gas_state_delta%mix_rat
153  else
154  gas_state%mix_rat = - gas_state_delta%mix_rat
155  end if
156  end if
157 
158  end subroutine gas_state_sub
159 
160 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
161 
162  !> Set any negative values to zero.
163  subroutine gas_state_ensure_nonnegative(gas_state)
165  !> Gas state.
166  type(gas_state_t) :: gas_state
167 
168  if (gas_state_is_allocated(gas_state)) then
169  gas_state%mix_rat = max(gas_state%mix_rat, 0d0)
170  end if
171 
172  end subroutine gas_state_ensure_nonnegative
173 
174 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
175 
176  !> Convert (mol m^{-3}) to (ppb).
177  subroutine gas_state_mole_dens_to_ppb(gas_state, env_state)
179  !> Gas state.
180  type(gas_state_t), intent(inout) :: gas_state
181  !> Environment state.
182  type(env_state_t), intent(in) :: env_state
183 
184  call gas_state_scale(gas_state, 1d9 / env_state_air_molar_den(env_state))
185 
186  end subroutine gas_state_mole_dens_to_ppb
187 
188 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
189 
190  !> Determine the current gas_state and rate by interpolating at the
191  !> current time with the lists of gas_states and rates.
192  subroutine gas_state_interp_1d(gas_state_list, time_list, &
193  rate_list, time, gas_state, rate)
195  !> Gas states.
196  type(gas_state_t), intent(in) :: gas_state_list(:)
197  !> Times (s).
198  real(kind=dp), intent(in) :: time_list(size(gas_state_list))
199  !> Rates (s^{-1}).
200  real(kind=dp), intent(in) :: rate_list(size(gas_state_list))
201  !> Current time (s).
202  real(kind=dp), intent(in) :: time
203  !> Current gas state.
204  type(gas_state_t), intent(inout) :: gas_state
205  !> Current rate (s^{-1}).
206  real(kind=dp), intent(out) :: rate
207 
208  integer :: n, p
209  real(kind=dp) :: y, alpha
210 
211  n = size(gas_state_list)
212  p = find_1d(n, time_list, time)
213  if (p == 0) then
214  ! before the start, just use the first state and rate
215  gas_state = gas_state_list(1)
216  rate = rate_list(1)
217  elseif (p == n) then
218  ! after the end, just use the last state and rate
219  gas_state = gas_state_list(n)
220  rate = rate_list(n)
221  else
222  ! in the middle, use the previous state
223  gas_state = gas_state_list(p)
224  rate = rate_list(p)
225  end if
226 
227  end subroutine gas_state_interp_1d
228 
229 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
230 
231  !> Read gas state from the file named on the line read from file.
232  subroutine spec_file_read_gas_state(file, gas_data, gas_state)
233 
234  !> File to read gas state from.
235  type(spec_file_t), intent(inout) :: file
236  !> Gas data.
237  type(gas_data_t), intent(in) :: gas_data
238  !> Gas data to read.
239  type(gas_state_t), intent(inout) :: gas_state
240 
241  integer :: n_species, species, i
242  character(len=SPEC_LINE_MAX_VAR_LEN), allocatable :: species_name(:)
243  real(kind=dp), allocatable :: species_data(:,:)
244 
245  !> \page input_format_gas_state Input File Format: Gas State
246  !!
247  !! A gas state input file must consist of one line per gas
248  !! species, with each line having the species name followed by the
249  !! species mixing ratio in ppb (parts per billion). The valid
250  !! species names are those specfied by the \ref
251  !! input_format_gas_data file, but not all species have to be
252  !! listed. Any missing species will have mixing ratios of
253  !! zero. For example, a gas state file could contain:
254  !! <pre>
255  !! # gas mixing ratio (ppb)
256  !! H2SO4 0
257  !! HNO3 1
258  !! HCl 0.7
259  !! NH3 0.5
260  !! </pre>
261  !!
262  !! See also:
263  !! - \ref spec_file_format --- the input file text format
264  !! - \ref output_format_gas_state --- the corresponding output format
265  !! - \ref input_format_gas_data --- the gas species list and
266  !! material data
267 
268  call spec_file_read_real_named_array(file, 0, species_name, &
269  species_data)
270 
271  ! check the data size
272  n_species = size(species_data, 1)
273  if (.not. ((size(species_data, 2) == 1) .or. (n_species == 0))) then
274  call die_msg(686719840, 'each line in ' // trim(file%name) &
275  // ' must contain exactly one data value')
276  end if
277 
278  ! copy over the data
279  call gas_state_set_size(gas_state, gas_data_n_spec(gas_data))
280  do i = 1,n_species
281  species = gas_data_spec_by_name(gas_data, species_name(i))
282  if (species == 0) then
283  call die_msg(129794076, 'unknown species ' // &
284  trim(species_name(i)) // ' in file ' // trim(file%name))
285  end if
286  gas_state%mix_rat(species) = species_data(i,1)
287  end do
288 
289  end subroutine spec_file_read_gas_state
290 
291 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
292 
293  !> Read an array of gas states with associated times and rates from
294  !> the file named on the line read from the given file.
295  subroutine spec_file_read_gas_states_times_rates(file, gas_data, &
296  times, rates, gas_states)
297 
298  !> Spec file.
299  type(spec_file_t), intent(inout) :: file
300  !> Gas data.
301  type(gas_data_t), intent(in) :: gas_data
302  !> Times (s).
303  real(kind=dp), allocatable :: times(:)
304  !> Rates (s^{-1}).
305  real(kind=dp), allocatable :: rates(:)
306  !> Gas states.
307  type(gas_state_t), allocatable :: gas_states(:)
308 
309  integer :: n_lines, species, i, n_time, i_time
310  character(len=SPEC_LINE_MAX_VAR_LEN), allocatable :: species_name(:)
311  real(kind=dp), allocatable :: species_data(:,:)
312 
313  !> \page input_format_gas_profile Input File Format: Gas Profile
314  !!
315  !! A gas profile input file must consist of three or more
316  !! lines, consisting of:
317  !! - the first line must begin with \c time and should be followed
318  !! by \f$N\f$ space-separated real scalars, giving the times (in
319  !! s after the start of the simulation) of the gas set
320  !! points --- the times must be in increasing order
321  !! - the second line must begin with \c rate and should be
322  !! followed by \f$N\f$ space-separated real scalars, giving the
323  !! values at the corresponding times
324  !! - the third and subsequent lines specify gas species, one
325  !! species per line, with each line beginning with the species
326  !! name and followed by \f$N\f$ space-separated scalars giving
327  !! the gas state of that species at the corresponding times
328  !!
329  !! The units and meanings of the rate and species lines depends on
330  !! the type of gas profile:
331  !! - emissions gas profiles have dimensionless rates that are used
332  !! to scale the species rates and species giving emission rates
333  !! with units of mol/(m^2 s) --- the emission rate is divided by
334  !! the current mixing layer height to give a per-volume emission
335  !! rate
336  !! - background gas profiles have rates with units s^{-1} that are
337  !! dilution rates and species with units of ppb (parts per
338  !! billion) that are the background mixing ratios
339  !!
340  !! The species names must be those specified by the \ref
341  !! input_format_gas_data. Any species not listed are taken to be
342  !! zero.
343  !!
344  !! Between the specified times the gas profile is interpolated
345  !! step-wise and kept constant at its last value. That is, if the
346  !! times are \f$t_i\f$, the rates are \f$r_i\f$, and the gas
347  !! states are \f$g_i\f$ (all with \f$i = 1,\ldots,n\f$), then
348  !! between times \f$t_i\f$ and \f$t_{i+1}\f$ the gas state is
349  !! constant at \f$r_i g_i\f$. Before time \f$t_1\f$ the gas state
350  !! is \f$r_1 g_1\f$, while after time \f$t_n\f$ it is \f$r_n
351  !! g_n\f$.
352  !!
353  !! Example: an emissions gas profile could be:
354  !! <pre>
355  !! time 0 600 1800 # time (in s) after simulation start
356  !! rate 1 0.5 1 # scaling factor
357  !! H2SO4 0 0 0 # emission rate in mol/(m^2 s)
358  !! SO2 4e-9 5.6e-9 5e-9 # emission rate in mol/(m^2 s)
359  !! </pre>
360  !! Here there are no emissions of \f$\rm H_2SO_4\f$, while \f$\rm
361  !! SO_2\f$ starts out being emitted at \f$4\times
362  !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ at the start of the simulation,
363  !! before falling to a rate of \f$2.8\times
364  !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ at 10&nbsp;min (note the scaling
365  !! of 0.5), and then rising again to \f$5\times
366  !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ after 30&nbsp;min. Between
367  !! 0&nbsp;min and 10&nbsp;min the emissions are the same as at
368  !! 0&nbsp;min, while between 10&nbsp;min and 30&nbsp;min they are
369  !! the same as at 10&nbsp;min. After 30&nbsp;min they are held
370  !! constant at their final value.
371  !!
372  !! See also:
373  !! - \ref spec_file_format --- the input file text format
374  !! - \ref input_format_gas_data --- the gas species list and
375  !! material data
376 
377  ! read the data from the file
378  call spec_file_read_real_named_array(file, 0, species_name, &
379  species_data)
380 
381  ! check the data size
382  n_lines = size(species_data, 1)
383  if (n_lines < 2) then
384  call die_msg(291542946, 'insufficient data lines in file ' &
385  // trim(file%name))
386  end if
387  if (trim(species_name(1)) /= 'time') then
388  call die_msg(525127793, 'row 1 in file ' &
389  // trim(file%name) // ' must start with: time')
390  end if
391  if (trim(species_name(2)) /= 'rate') then
392  call die_msg(506981322, 'row 2 in file ' &
393  // trim(file%name) // ' must start with: rate')
394  end if
395  n_time = size(species_data, 2)
396  if (n_time < 1) then
397  call die_msg(398532628, 'each line in file ' &
398  // trim(file%name) // ' must contain at least one data value')
399  end if
400 
401  ! copy over the data
402  times = species_data(1,:)
403  rates = species_data(2,:)
404  if (allocated(gas_states)) deallocate(gas_states)
405  allocate(gas_states(n_time))
406  do i_time = 1,n_time
407  call gas_state_set_size(gas_states(i_time), gas_data_n_spec(gas_data))
408  end do
409  do i = 3,n_lines
410  species = gas_data_spec_by_name(gas_data, species_name(i))
411  if (species == 0) then
412  call die_msg(806500079, 'unknown species ' &
413  // trim(species_name(i)) // ' in file ' &
414  // trim(file%name))
415  end if
416  do i_time = 1,n_time
417  gas_states(i_time)%mix_rat(species) = species_data(i,i_time)
418  end do
419  end do
420 
421  end subroutine spec_file_read_gas_states_times_rates
422 
423 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
424 
425  !> Average val over all processes.
426  subroutine gas_state_mix(val)
428  !> Value to average.
429  type(gas_state_t), intent(inout) :: val
430 
431 #ifdef PMC_USE_MPI
432  type(gas_state_t) :: val_avg
433 
434  call gas_state_set_size(val_avg, size(val%mix_rat))
435  call pmc_mpi_allreduce_average_real_array(val%mix_rat, val_avg%mix_rat)
436  val%mix_rat = val_avg%mix_rat
437 #endif
438 
439  end subroutine gas_state_mix
440 
441 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
442 
443  !> Average val over all processes, with the result only on the root
444  !> process.
445  subroutine gas_state_reduce_avg(val)
447  !> Value to average.
448  type(gas_state_t), intent(inout) :: val
449 
450 #ifdef PMC_USE_MPI
451  type(gas_state_t) :: val_avg
452 
453  call gas_state_set_size(val_avg, size(val%mix_rat))
454  call pmc_mpi_reduce_avg_real_array(val%mix_rat, val_avg%mix_rat)
455  if (pmc_mpi_rank() == 0) then
456  val%mix_rat = val_avg%mix_rat
457  end if
458 #endif
459 
460  end subroutine gas_state_reduce_avg
461 
462 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
463 
464  !> Determines the number of bytes required to pack the given value.
465  integer function pmc_mpi_pack_size_gas_state(val)
467  !> Value to pack.
468  type(gas_state_t), intent(in) :: val
469 
471  + pmc_mpi_pack_size_real_array(val%mix_rat)
472 
473  end function pmc_mpi_pack_size_gas_state
474 
475 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
476 
477  !> Packs the given value into the buffer, advancing position.
478  subroutine pmc_mpi_pack_gas_state(buffer, position, val)
480  !> Memory buffer.
481  character, intent(inout) :: buffer(:)
482  !> Current buffer position.
483  integer, intent(inout) :: position
484  !> Value to pack.
485  type(gas_state_t), intent(in) :: val
486 
487 #ifdef PMC_USE_MPI
488  integer :: prev_position
489 
490  prev_position = position
491  call pmc_mpi_pack_real_array(buffer, position, val%mix_rat)
492  call assert(655827004, &
493  position - prev_position <= pmc_mpi_pack_size_gas_state(val))
494 #endif
495 
496  end subroutine pmc_mpi_pack_gas_state
497 
498 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
499 
500  !> Unpacks the given value from the buffer, advancing position.
501  subroutine pmc_mpi_unpack_gas_state(buffer, position, val)
503  !> Memory buffer.
504  character, intent(inout) :: buffer(:)
505  !> Current buffer position.
506  integer, intent(inout) :: position
507  !> Value to pack.
508  type(gas_state_t), intent(inout) :: val
509 
510 #ifdef PMC_USE_MPI
511  integer :: prev_position
512 
513  prev_position = position
514  call pmc_mpi_unpack_real_array(buffer, position, val%mix_rat)
515  call assert(520815247, &
516  position - prev_position <= pmc_mpi_pack_size_gas_state(val))
517 #endif
518 
519  end subroutine pmc_mpi_unpack_gas_state
520 
521 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
522 
523  !> Computes the average of val across all processes, storing the
524  !> result in val_avg on the root process.
525  subroutine pmc_mpi_reduce_avg_gas_state(val, val_avg)
527  !> Value to average.
528  type(gas_state_t), intent(in) :: val
529  !> Result.
530  type(gas_state_t), intent(inout) :: val_avg
531 
532  call pmc_mpi_reduce_avg_real_array(val%mix_rat, val_avg%mix_rat)
533 
534  end subroutine pmc_mpi_reduce_avg_gas_state
535 
536 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
537 
538  !> Write full state.
539  subroutine gas_state_output_netcdf(gas_state, ncid, gas_data)
540 
541  !> Gas state to write.
542  type(gas_state_t), intent(in) :: gas_state
543  !> NetCDF file ID, in data mode.
544  integer, intent(in) :: ncid
545  !> Gas data.
546  type(gas_data_t), intent(in) :: gas_data
547 
548  integer :: dimid_gas_species
549 
550  !> \page output_format_gas_state Output File Format: Gas State
551  !!
552  !! The gas state uses the \c gas_species NetCDF dimension as specified
553  !! in the \ref output_format_gas_data section.
554  !!
555  !! The gas state NetCDF variables are:
556  !! - \b gas_mixing_ratio (unit ppb, dim \c gas_species): current mixing
557  !! ratios of each gas species
558  !!
559  !! See also:
560  !! - \ref output_format_gas_data --- the gas species list and material
561  !! data output format
562  !! - \ref input_format_gas_state --- the corresponding input format
563 
564  call gas_data_netcdf_dim_gas_species(gas_data, ncid, &
565  dimid_gas_species)
566 
567  call pmc_nc_write_real_1d(ncid, gas_state%mix_rat, &
568  "gas_mixing_ratio", (/ dimid_gas_species /), unit="ppb", &
569  long_name="mixing ratios of gas species")
570 
571  end subroutine gas_state_output_netcdf
572 
573 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
574 
575  !> Read full state.
576  subroutine gas_state_input_netcdf(gas_state, ncid, gas_data)
578  !> Gas state to read.
579  type(gas_state_t), intent(inout) :: gas_state
580  !> NetCDF file ID, in data mode.
581  integer, intent(in) :: ncid
582  !> Gas data.
583  type(gas_data_t), intent(in) :: gas_data
584 
585  call pmc_nc_read_real_1d(ncid, gas_state%mix_rat, "gas_mixing_ratio")
586 
587  end subroutine gas_state_input_netcdf
588 
589 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
590 
591 end module pmc_gas_state
real(kind=dp) function env_state_air_molar_den(env_state)
Air molar density (mol m^{-3}).
Definition: env_state.F90:164
An input file with extra data for printing messages.
Definition: spec_file.F90:59
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_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from a NetCDF file.
Definition: netcdf.F90:219
subroutine gas_state_add_scaled(gas_state, gas_state_delta, alpha)
Adds the given gas_state_delta scaled by alpha.
Definition: gas_state.F90:121
integer function pmc_mpi_rank()
Returns the rank of the current process.
Definition: mpi.F90:117
subroutine gas_state_input_netcdf(gas_state, ncid, gas_data)
Read full state.
Definition: gas_state.F90:577
subroutine gas_state_mole_dens_to_ppb(gas_state, env_state)
Convert (mol m^{-3}) to (ppb).
Definition: gas_state.F90:178
The env_state_t structure and associated subroutines.
Definition: env_state.F90:9
subroutine gas_state_interp_1d(gas_state_list, time_list, rate_list, time, gas_state, rate)
Determine the current gas_state and rate by interpolating at the current time with the lists of gas_s...
Definition: gas_state.F90:194
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:477
The gas_data_t structure and associated subroutines.
Definition: gas_data.F90:9
subroutine gas_state_ensure_nonnegative(gas_state)
Set any negative values to zero.
Definition: gas_state.F90:164
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:103
subroutine pmc_mpi_allreduce_average_real_array(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on all processes...
Definition: mpi.F90:1305
Current environment state.
Definition: env_state.F90:26
subroutine pmc_mpi_reduce_avg_real_array(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: mpi.F90:1226
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:134
subroutine pmc_mpi_reduce_avg_gas_state(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: gas_state.F90:526
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 gas_state_set_size(gas_state, n_spec)
Sets the sizes of the gas state.
Definition: gas_state.F90:53
The gas_state_t structure and associated subroutines.
Definition: gas_state.F90:9
Reading formatted text input.
Definition: spec_file.F90:43
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:720
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:978
integer function find_1d(n, x_vals, x)
Find the position of a real number in an arbitrary 1D array.
Definition: util.F90:592
subroutine gas_state_mix(val)
Average val over all processes.
Definition: gas_state.F90:427
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
subroutine gas_state_zero(gas_state)
Zeros the state.
Definition: gas_state.F90:69
Current state of the gas mixing ratios in the system.
Definition: gas_state.F90:30
subroutine pmc_nc_write_real_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple real array to a NetCDF file.
Definition: netcdf.F90:536
logical function gas_state_is_allocated(gas_state)
Determine whether the gas_state is correctly allocated.
Definition: gas_state.F90:41
subroutine gas_state_reduce_avg(val)
Average val over all processes, with the result only on the root process.
Definition: gas_state.F90:446
Constant gas data.
Definition: gas_data.F90:29
subroutine pmc_mpi_unpack_gas_state(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: gas_state.F90:502
subroutine gas_state_add(gas_state, gas_state_delta)
Adds the given gas_state_delta.
Definition: gas_state.F90:99
subroutine gas_state_scale(gas_state, alpha)
Scale a gas state.
Definition: gas_state.F90:83
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
subroutine pmc_mpi_pack_gas_state(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: gas_state.F90:479
subroutine gas_state_sub(gas_state, gas_state_delta)
Subtracts the given gas_state_delta.
Definition: gas_state.F90:144
Common utility subroutines.
Definition: util.F90:9
Wrapper functions for MPI.
Definition: mpi.F90:13
integer function pmc_mpi_pack_size_gas_state(val)
Determines the number of bytes required to pack the given value.
Definition: gas_state.F90:466