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