PartMC  2.2.0
gas_state.F90
Go to the documentation of this file.
00001 ! Copyright (C) 2007-2011 Matthew West
00002 ! Licensed under the GNU General Public License version 2 or (at your
00003 ! option) any later version. See the file COPYING for details.
00004 
00005 !> \file
00006 !> The pmc_gas_state module.
00007 
00008 !> The gas_state_t structure and associated subroutines.
00009 module pmc_gas_state
00010 
00011   use pmc_util
00012   use pmc_spec_file
00013   use pmc_gas_data
00014   use pmc_mpi
00015   use pmc_netcdf
00016 #ifdef PMC_USE_MPI
00017   use mpi
00018 #endif
00019 
00020   !> Current state of the gas mixing ratios in the system.
00021   !!
00022   !! The gas species are defined by the gas_data_t structure, so that
00023   !! \c gas_state%%mix_rat(i) is the current mixing ratio of the gas
00024   !! with name \c gas_data%%name(i), etc.
00025   type gas_state_t
00026      !> Length n_spec, mixing ratio (ppb).
00027      real(kind=dp), pointer :: mix_rat(:)
00028   end type gas_state_t
00029 
00030 contains
00031 
00032 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00033 
00034   !> Allocate storage for gas species.
00035   subroutine gas_state_allocate(gas_state)
00036 
00037     !> Gas state to be allocated.
00038     type(gas_state_t), intent(out) :: gas_state
00039 
00040     allocate(gas_state%mix_rat(0))
00041 
00042   end subroutine gas_state_allocate
00043 
00044 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00045 
00046   !> Allocate storage for gas species of the given size.
00047   subroutine gas_state_allocate_size(gas_state, n_spec)
00048 
00049     !> Gas state to be allocated.
00050     type(gas_state_t), intent(out) :: gas_state
00051     !> Number of species.
00052     integer, intent(in) :: n_spec
00053 
00054     allocate(gas_state%mix_rat(n_spec))
00055     call gas_state_zero(gas_state)
00056 
00057   end subroutine gas_state_allocate_size
00058 
00059 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00060 
00061   !> Free all storage.
00062   subroutine gas_state_deallocate(gas_state)
00063 
00064     !> Gas state to be freed.
00065     type(gas_state_t), intent(inout) :: gas_state
00066 
00067     deallocate(gas_state%mix_rat)
00068 
00069   end subroutine gas_state_deallocate
00070 
00071 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00072 
00073   !> Zeros the state.
00074   subroutine gas_state_zero(gas_state)
00075 
00076     !> Gas state.
00077     type(gas_state_t), intent(inout) :: gas_state
00078 
00079     gas_state%mix_rat = 0d0
00080 
00081   end subroutine gas_state_zero
00082 
00083 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00084 
00085   !> Copy to an already allocated to_state.
00086   subroutine gas_state_copy(from_state, to_state)
00087 
00088     !> Existing gas state.
00089     type(gas_state_t), intent(in) :: from_state
00090     !> Must be allocated already.
00091     type(gas_state_t), intent(inout) :: to_state
00092 
00093     integer :: n_spec
00094 
00095     n_spec = size(from_state%mix_rat)
00096     deallocate(to_state%mix_rat)
00097     allocate(to_state%mix_rat(n_spec))
00098     to_state%mix_rat = from_state%mix_rat
00099 
00100   end subroutine gas_state_copy
00101 
00102 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00103 
00104   !> Scale a gas state.
00105   subroutine gas_state_scale(gas_state, alpha)
00106 
00107     !> Existing gas state.
00108     type(gas_state_t), intent(inout) :: gas_state
00109     !> Scale factor.
00110     real(kind=dp), intent(in) :: alpha
00111 
00112     gas_state%mix_rat = gas_state%mix_rat * alpha
00113 
00114   end subroutine gas_state_scale
00115 
00116 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00117 
00118   !> Adds the given gas_state_delta.
00119   subroutine gas_state_add(gas_state, gas_state_delta)
00120 
00121     !> Existing gas state.
00122     type(gas_state_t), intent(inout) :: gas_state
00123     !> Incremental state.
00124     type(gas_state_t), intent(in) :: gas_state_delta
00125 
00126     gas_state%mix_rat = gas_state%mix_rat + gas_state_delta%mix_rat
00127 
00128   end subroutine gas_state_add
00129 
00130 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00131 
00132   !> Subtracts the given gas_state_delta.
00133   subroutine gas_state_sub(gas_state, gas_state_delta)
00134 
00135     !> Existing gas state.
00136     type(gas_state_t), intent(inout) :: gas_state
00137     !> Incremental state.
00138     type(gas_state_t), intent(in) :: gas_state_delta
00139 
00140     gas_state%mix_rat = gas_state%mix_rat - gas_state_delta%mix_rat
00141 
00142   end subroutine gas_state_sub
00143 
00144 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00145 
00146   !> Sets y = a * x + y.
00147   subroutine gas_state_axpy(alpha, gas_state_x, gas_state_y)
00148 
00149     !> Coefficient.
00150     real(kind=dp), intent(in) :: alpha
00151     !> Incremental state.
00152     type(gas_state_t), intent(in) :: gas_state_x
00153     !> Existing gas state.
00154     type(gas_state_t), intent(inout) :: gas_state_y
00155 
00156     gas_state_y%mix_rat = alpha * gas_state_x%mix_rat + gas_state_y%mix_rat
00157 
00158   end subroutine gas_state_axpy
00159 
00160 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00161 
00162   !> Set any negative values to zero.
00163   subroutine gas_state_ensure_nonnegative(gas_state)
00164 
00165     !> Gas state.
00166     type(gas_state_t) :: gas_state
00167 
00168     gas_state%mix_rat = max(gas_state%mix_rat, 0d0)
00169 
00170   end subroutine gas_state_ensure_nonnegative
00171 
00172 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00173 
00174   !> Determine the current gas_state and rate by interpolating at the
00175   !> current time with the lists of gas_states and rates.
00176   subroutine gas_state_interp_1d(gas_state_list, time_list, &
00177          rate_list, time, gas_state, rate)
00178 
00179     !> Gas states.
00180     type(gas_state_t), intent(in) :: gas_state_list(:)
00181     !> Times (s).
00182     real(kind=dp), intent(in) :: time_list(size(gas_state_list))
00183     !> Rates (s^{-1}).
00184     real(kind=dp), intent(in) :: rate_list(size(gas_state_list))
00185     !> Current time (s).
00186     real(kind=dp), intent(in) :: time
00187     !> Current gas state.
00188     type(gas_state_t), intent(inout) :: gas_state
00189     !> Current rate (s^{-1}).
00190     real(kind=dp), intent(out) :: rate
00191 
00192     integer :: n, p
00193     real(kind=dp) :: y, alpha
00194 
00195     n = size(gas_state_list)
00196     p = find_1d(n, time_list, time)
00197     if (p == 0) then
00198        ! before the start, just use the first state and rate
00199        call gas_state_copy(gas_state_list(1), gas_state)
00200        rate = rate_list(1)
00201     elseif (p == n) then
00202        ! after the end, just use the last state and rate
00203        call gas_state_copy(gas_state_list(n), gas_state)
00204        rate = rate_list(n)
00205     else
00206        ! in the middle, use the previous state
00207        call gas_state_copy(gas_state_list(p), gas_state)
00208        rate = rate_list(p)
00209     end if
00210 
00211   end subroutine gas_state_interp_1d
00212 
00213 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00214 
00215   !> Read gas state from the file named on the line read from file.
00216   subroutine spec_file_read_gas_state(file, gas_data, gas_state)
00217 
00218     !> File to read gas state from.
00219     type(spec_file_t), intent(inout) :: file
00220     !> Gas data.
00221     type(gas_data_t), intent(in) :: gas_data
00222     !> Gas data to read.
00223     type(gas_state_t), intent(inout) :: gas_state
00224 
00225     integer :: n_species, species, i
00226     character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: species_name(:)
00227     real(kind=dp), pointer :: species_data(:,:)
00228 
00229     !> \page input_format_gas_state Input File Format: Gas State
00230     !!
00231     !! A gas state input file must consist of one line per gas
00232     !! species, with each line having the species name followed by the
00233     !! species mixing ratio in ppb (parts per billion). The valid
00234     !! species names are those specfied by the \ref
00235     !! input_format_gas_data file, but not all species have to be
00236     !! listed. Any missing species will have mixing ratios of
00237     !! zero. For example, a gas state file could contain:
00238     !! <pre>
00239     !! # gas  mixing ratio (ppb)
00240     !! H2SO4  0
00241     !! HNO3   1
00242     !! HCl    0.7
00243     !! NH3    0.5
00244     !! </pre>
00245     !!
00246     !! See also:
00247     !!   - \ref spec_file_format --- the input file text format
00248     !!   - \ref output_format_gas_state --- the corresponding output format
00249     !!   - \ref input_format_gas_data --- the gas species list and
00250     !!     material data
00251 
00252     allocate(species_name(0))
00253     allocate(species_data(0,0))
00254     call spec_file_read_real_named_array(file, 0, species_name, &
00255          species_data)
00256 
00257     ! check the data size
00258     n_species = size(species_data, 1)
00259     if (.not. ((size(species_data, 2) == 1) .or. (n_species == 0))) then
00260        call die_msg(686719840, 'each line in ' // trim(file%name) &
00261             // ' must contain exactly one data value')
00262     end if
00263 
00264     ! copy over the data
00265     call gas_state_deallocate(gas_state)
00266     call gas_state_allocate_size(gas_state, gas_data%n_spec)
00267     gas_state%mix_rat = 0d0
00268     do i = 1,n_species
00269        species = gas_data_spec_by_name(gas_data, species_name(i))
00270        if (species == 0) then
00271           call die_msg(129794076, 'unknown species ' // &
00272                trim(species_name(i)) // ' in file ' // trim(file%name))
00273        end if
00274        gas_state%mix_rat(species) = species_data(i,1)
00275     end do
00276     deallocate(species_name)
00277     deallocate(species_data)
00278 
00279   end subroutine spec_file_read_gas_state
00280 
00281 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00282 
00283   !> Read an array of gas states with associated times and rates from
00284   !> the file named on the line read from the given file.
00285   subroutine spec_file_read_gas_states_times_rates(file, gas_data, &
00286        times, rates, gas_states)
00287 
00288     !> Spec file.
00289     type(spec_file_t), intent(inout) :: file
00290     !> Gas data.
00291     type(gas_data_t), intent(in) :: gas_data
00292     !> Times (s).
00293     real(kind=dp), pointer :: times(:)
00294     !> Rates (s^{-1}).
00295     real(kind=dp), pointer :: rates(:)
00296     !> Gas states.
00297     type(gas_state_t), pointer :: gas_states(:)
00298 
00299     integer :: n_lines, species, i, n_time, i_time
00300     character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: species_name(:)
00301     real(kind=dp), pointer :: species_data(:,:)
00302 
00303     !> \page input_format_gas_profile Input File Format: Gas Profile
00304     !!
00305     !! A gas profile input file must consist of three or more
00306     !! lines, consisting of:
00307     !! - the first line must begin with \c time and should be followed
00308     !!   by \f$N\f$ space-separated real scalars, giving the times (in
00309     !!   s after the start of the simulation) of the gas set
00310     !!   points --- the times must be in increasing order
00311     !! - the second line must begin with \c rate and should be
00312     !!   followed by \f$N\f$ space-separated real scalars, giving the
00313     !!   scalings (dimensionless) at the corresponding times
00314     !! - the third and subsequent lines specify gas species, one
00315     !!   species per line, with each line beginning with the species
00316     !!   name and followed by \f$N\f$ space-separated scalars giving
00317     !!   the gas state of that species at the corresponding times
00318     !!
00319     !! The units of the species lines depends on the type of gas
00320     !! profile:
00321     !! - emissions gas profiles have units of mol/(m^2 s) --- the
00322     !!   emission rate is divided by the current mixing layer
00323     !!   height to give a per-volume emission rate
00324     !! - background gas profiles have units of ppb (parts per billion)
00325     !!
00326     !! The species names must be those specified by the \ref
00327     !! input_format_gas_data. Any species not listed are taken to be
00328     !! zero.
00329     !!
00330     !! The scaling parameter is used to multiply all the gas state
00331     !! values at the corresponding time, giving a simple way of
00332     !! scaling the overall gas state.
00333     !!
00334     !! Between the specified times the gas profile is interpolated
00335     !! step-wise and kept constant at its last value. That is, if the
00336     !! times are \f$t_i\f$, the rates are \f$r_i\f$, and the gas
00337     !! states are \f$g_i\f$ (all with \f$i = 1,\ldots,n\f$), then
00338     !! between times \f$t_i\f$ and \f$t_{i+1}\f$ the gas state is
00339     !! constant at \f$r_i g_i\f$. Before time \f$t_1\f$ the gas state
00340     !! is \f$r_1 g_1\f$, while after time \f$t_n\f$ it is \f$r_n
00341     !! g_n\f$.
00342     !!
00343     !! Example: an emissions gas profile could be:
00344     !! <pre>
00345     !! time   0       600     1800    # time (in s) after simulation start
00346     !! rate   1       0.5     1       # scaling factor
00347     !! H2SO4  0       0       0       # emission rate in mol/(m^2 s)
00348     !! SO2    4e-9    5.6e-9  5e-9    # emission rate in mol/(m^2 s)
00349     !! </pre>
00350     !! Here there are no emissions of \f$\rm H_2SO_4\f$, while \f$\rm
00351     !! SO_2\f$ starts out being emitted at \f$4\times
00352     !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ at the start of the simulation,
00353     !! before falling to a rate of \f$2.8\times 
00354     !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ at 10&nbsp;min (note the scaling
00355     !! of 0.5), and then rising again to \f$5\times
00356     !! 10^{-9}\rm\ mol\ m^{-2}\ s^{-1}\f$ after 30&nbsp;min. Between
00357     !! 0&nbsp;min and 10&nbsp;min the emissions are the same as at
00358     !! 0&nbsp;min, while between 10&nbsp;min and 30&nbsp;min they are
00359     !! the same as at 10&nbsp;min. After 30&nbsp;min they are held
00360     !! constant at their final value.
00361     !!
00362     !! See also:
00363     !!   - \ref spec_file_format --- the input file text format
00364     !!   - \ref input_format_gas_data --- the gas species list and
00365     !!     material data
00366 
00367     ! read the data from the file
00368     allocate(species_name(0))
00369     allocate(species_data(0,0))
00370     call spec_file_read_real_named_array(file, 0, species_name, &
00371          species_data)
00372 
00373     ! check the data size
00374     n_lines = size(species_data, 1)
00375     if (n_lines < 2) then
00376        call die_msg(291542946, 'insufficient data lines in file ' &
00377             // trim(file%name))
00378     end if
00379     if (trim(species_name(1)) /= 'time') then
00380        call die_msg(525127793, 'row 1 in file ' &
00381             // trim(file%name) // ' must start with: time')
00382     end if
00383     if (trim(species_name(2)) /= 'rate') then
00384        call die_msg(506981322, 'row 2 in file ' &
00385             // trim(file%name) // ' must start with: rate')
00386     end if
00387     n_time = size(species_data, 2)
00388     if (n_time < 1) then
00389        call die_msg(398532628, 'each line in file ' &
00390             // trim(file%name) // ' must contain at least one data value')
00391     end if
00392 
00393     ! copy over the data
00394     do i_time = 1,size(gas_states)
00395        call gas_state_deallocate(gas_states(i_time))
00396     end do
00397     deallocate(gas_states)
00398     deallocate(times)
00399     deallocate(rates)
00400     allocate(gas_states(n_time))
00401     allocate(times(n_time))
00402     allocate(rates(n_time))
00403     do i_time = 1,n_time
00404        call gas_state_allocate_size(gas_states(i_time), gas_data%n_spec)
00405        times(i_time) = species_data(1,i_time)
00406        rates(i_time) = species_data(2,i_time)
00407     end do
00408     do i = 3,n_lines
00409        species = gas_data_spec_by_name(gas_data, species_name(i))
00410        if (species == 0) then
00411           call die_msg(806500079, 'unknown species ' &
00412                // trim(species_name(i)) // ' in file ' &
00413                // trim(file%name))
00414        end if
00415        do i_time = 1,n_time
00416           gas_states(i_time)%mix_rat(species) = species_data(i,i_time)
00417        end do
00418     end do
00419     deallocate(species_name)
00420     deallocate(species_data)
00421 
00422   end subroutine spec_file_read_gas_states_times_rates
00423 
00424 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00425 
00426   !> Average val over all processes.
00427   subroutine gas_state_mix(val)
00428 
00429     !> Value to average.
00430     type(gas_state_t), intent(inout) :: val
00431 
00432 #ifdef PMC_USE_MPI
00433     type(gas_state_t) :: val_avg
00434 
00435     call gas_state_allocate_size(val_avg, size(val%mix_rat))
00436     call pmc_mpi_allreduce_average_real_array(val%mix_rat, val_avg%mix_rat)
00437     val%mix_rat = val_avg%mix_rat
00438     call gas_state_deallocate(val_avg)
00439 #endif
00440 
00441   end subroutine gas_state_mix
00442 
00443 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00444 
00445   !> Average val over all processes, with the result only on the root
00446   !> process.
00447   subroutine gas_state_reduce_avg(val)
00448 
00449     !> Value to average.
00450     type(gas_state_t), intent(inout) :: val
00451 
00452 #ifdef PMC_USE_MPI
00453     type(gas_state_t) :: val_avg
00454 
00455     call gas_state_allocate_size(val_avg, size(val%mix_rat))
00456     call pmc_mpi_reduce_avg_real_array(val%mix_rat, val_avg%mix_rat)
00457     if (pmc_mpi_rank() == 0) then
00458        val%mix_rat = val_avg%mix_rat
00459     end if
00460     call gas_state_deallocate(val_avg)
00461 #endif
00462 
00463   end subroutine gas_state_reduce_avg
00464 
00465 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00466 
00467   !> Determines the number of bytes required to pack the given value.
00468   integer function pmc_mpi_pack_size_gas_state(val)
00469 
00470     !> Value to pack.
00471     type(gas_state_t), intent(in) :: val
00472 
00473     pmc_mpi_pack_size_gas_state = &
00474          + pmc_mpi_pack_size_real_array(val%mix_rat)
00475 
00476   end function pmc_mpi_pack_size_gas_state
00477 
00478 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00479 
00480   !> Packs the given value into the buffer, advancing position.
00481   subroutine pmc_mpi_pack_gas_state(buffer, position, val)
00482 
00483     !> Memory buffer.
00484     character, intent(inout) :: buffer(:)
00485     !> Current buffer position.
00486     integer, intent(inout) :: position
00487     !> Value to pack.
00488     type(gas_state_t), intent(in) :: val
00489 
00490 #ifdef PMC_USE_MPI
00491     integer :: prev_position
00492 
00493     prev_position = position
00494     call pmc_mpi_pack_real_array(buffer, position, val%mix_rat)
00495     call assert(655827004, &
00496          position - prev_position <= pmc_mpi_pack_size_gas_state(val))
00497 #endif
00498 
00499   end subroutine pmc_mpi_pack_gas_state
00500 
00501 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00502 
00503   !> Unpacks the given value from the buffer, advancing position.
00504   subroutine pmc_mpi_unpack_gas_state(buffer, position, val)
00505 
00506     !> Memory buffer.
00507     character, intent(inout) :: buffer(:)
00508     !> Current buffer position.
00509     integer, intent(inout) :: position
00510     !> Value to pack.
00511     type(gas_state_t), intent(inout) :: val
00512 
00513 #ifdef PMC_USE_MPI
00514     integer :: prev_position
00515 
00516     prev_position = position
00517     call pmc_mpi_unpack_real_array(buffer, position, val%mix_rat)
00518     call assert(520815247, &
00519          position - prev_position <= pmc_mpi_pack_size_gas_state(val))
00520 #endif
00521 
00522   end subroutine pmc_mpi_unpack_gas_state
00523 
00524 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00525 
00526   !> Computes the average of val across all processes, storing the
00527   !> result in val_avg on the root process.
00528   subroutine pmc_mpi_reduce_avg_gas_state(val, val_avg)
00529 
00530     !> Value to average.
00531     type(gas_state_t), intent(in) :: val
00532     !> Result.
00533     type(gas_state_t), intent(inout) :: val_avg
00534 
00535     call pmc_mpi_reduce_avg_real_array(val%mix_rat, val_avg%mix_rat)
00536 
00537   end subroutine pmc_mpi_reduce_avg_gas_state
00538 
00539 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00540 
00541   !> Write full state.
00542   subroutine gas_state_output_netcdf(gas_state, ncid, gas_data)
00543     
00544     !> Gas state to write.
00545     type(gas_state_t), intent(in) :: gas_state
00546     !> NetCDF file ID, in data mode.
00547     integer, intent(in) :: ncid
00548     !> Gas data.
00549     type(gas_data_t), intent(in) :: gas_data
00550 
00551     integer :: dimid_gas_species
00552 
00553     !> \page output_format_gas_state Output File Format: Gas State
00554     !!
00555     !! The gas state uses the \c gas_species NetCDF dimension as specified
00556     !! in the \ref output_format_gas_data section.
00557     !!
00558     !! The gas state NetCDF variables are:
00559     !!   - \b gas_mixing_ratio (unit ppb, dim \c gas_species): current mixing
00560     !!     ratios of each gas species
00561     !!
00562     !! See also:
00563     !!   - \ref output_format_gas_data --- the gas species list and material
00564     !!     data output format
00565     !!   - \ref input_format_gas_state --- the corresponding input format
00566 
00567     call gas_data_netcdf_dim_gas_species(gas_data, ncid, &
00568          dimid_gas_species)
00569 
00570     call pmc_nc_write_real_1d(ncid, gas_state%mix_rat, &
00571          "gas_mixing_ratio", (/ dimid_gas_species /), unit="ppb", &
00572          long_name="mixing ratios of gas species")
00573 
00574   end subroutine gas_state_output_netcdf
00575 
00576 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00577 
00578   !> Read full state.
00579   subroutine gas_state_input_netcdf(gas_state, ncid, gas_data)
00580     
00581     !> Gas state to read.
00582     type(gas_state_t), intent(inout) :: gas_state
00583     !> NetCDF file ID, in data mode.
00584     integer, intent(in) :: ncid
00585     !> Gas data.
00586     type(gas_data_t), intent(in) :: gas_data
00587 
00588     call gas_state_deallocate(gas_state)
00589     call gas_state_allocate_size(gas_state, gas_data%n_spec)
00590     call pmc_nc_read_real_1d(ncid, gas_state%mix_rat, &
00591          "gas_mixing_ratio")
00592 
00593   end subroutine gas_state_input_netcdf
00594 
00595 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00596 
00597 end module pmc_gas_state