Go to the documentation of this file.
35 real(kind=
dp),
allocatable :: mix_rat(:)
60 integer,
intent(in) :: n_spec
62 if (
allocated(gas_state%mix_rat))
deallocate(gas_state%mix_rat)
63 allocate(gas_state%mix_rat(n_spec))
77 gas_state%mix_rat = 0d0
90 real(kind=
dp),
intent(in) :: alpha
93 gas_state%mix_rat = gas_state%mix_rat * alpha
110 gas_state%mix_rat = gas_state%mix_rat + gas_state_delta%mix_rat
112 gas_state%mix_rat = gas_state_delta%mix_rat
130 real(kind=
dp),
intent(in) :: alpha
134 gas_state%mix_rat = gas_state%mix_rat &
135 + alpha * gas_state_delta%mix_rat
137 gas_state%mix_rat = alpha * gas_state_delta%mix_rat
155 gas_state%mix_rat = gas_state%mix_rat - gas_state_delta%mix_rat
157 gas_state%mix_rat = - gas_state_delta%mix_rat
172 gas_state%mix_rat = max(gas_state%mix_rat, 0d0)
195 subroutine gas_state_set_camp_conc(gas_state, camp_state, gas_data)
200 type(camp_state_t),
intent(inout) :: camp_state
204 real(kind=
dp),
parameter :: t_steam = 373.15
205 real(kind=
dp) :: a, water_vp
207 camp_state%state_var(1:
size(gas_state%mix_rat)) = gas_state%mix_rat(:) &
214 call assert(590005048,
associated(camp_state%env_states(1)%val))
215 a = 1.0 - t_steam / camp_state%env_states(1)%val%temp
216 a = (((-0.1299 * a - 0.6445) * a - 1.976) * a + 13.3185) * a
217 water_vp = 101325.0 * exp(a)
219 camp_state%state_var(gas_data%i_camp_water) = &
220 camp_state%env_states(1)%val%rel_humid * water_vp * 1.0e6 &
221 / camp_state%env_states(1)%val%pressure
223 end subroutine gas_state_set_camp_conc
228 subroutine gas_state_get_camp_conc(gas_state, camp_state)
233 type(camp_state_t),
intent(in) :: camp_state
235 gas_state%mix_rat(:) = camp_state%state_var(1:
size(gas_state%mix_rat)) &
238 end subroutine gas_state_get_camp_conc
246 rate_list, time, gas_state, rate)
251 real(kind=
dp),
intent(in) :: time_list(
size(gas_state_list))
253 real(kind=
dp),
intent(in) :: rate_list(
size(gas_state_list))
255 real(kind=
dp),
intent(in) :: time
259 real(kind=
dp),
intent(out) :: rate
262 real(kind=
dp) :: y, alpha
264 n =
size(gas_state_list)
265 p =
find_1d(n, time_list, time)
268 gas_state = gas_state_list(1)
272 gas_state = gas_state_list(n)
276 gas_state = gas_state_list(p)
285 subroutine spec_file_read_gas_state(file, gas_data, gas_state)
294 integer :: n_species, species, i
295 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: species_name(:)
296 real(kind=
dp),
allocatable :: species_data(:,:)
325 n_species =
size(species_data, 1)
326 if (.not. ((
size(species_data, 2) == 1) .or. (n_species == 0)))
then
327 call die_msg(686719840,
'each line in ' // trim(file%name) &
328 //
' must contain exactly one data value')
335 if (species == 0)
then
336 call die_msg(129794076,
'unknown species ' // &
337 trim(species_name(i)) //
' in file ' // trim(file%name))
339 gas_state%mix_rat(species) = species_data(i,1)
342 end subroutine spec_file_read_gas_state
348 subroutine spec_file_read_gas_states_times_rates(file, gas_data, &
349 times, rates, gas_states)
356 real(kind=
dp),
allocatable :: times(:)
358 real(kind=
dp),
allocatable :: rates(:)
362 integer :: n_lines, species, i, n_time, i_time
363 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: species_name(:)
364 real(kind=
dp),
allocatable :: species_data(:,:)
435 n_lines =
size(species_data, 1)
436 if (n_lines < 2)
then
437 call die_msg(291542946,
'insufficient data lines in file ' &
440 if (trim(species_name(1)) /=
'time')
then
441 call die_msg(525127793,
'row 1 in file ' &
442 // trim(file%name) //
' must start with: time')
444 if (trim(species_name(2)) /=
'rate')
then
445 call die_msg(506981322,
'row 2 in file ' &
446 // trim(file%name) //
' must start with: rate')
448 n_time =
size(species_data, 2)
450 call die_msg(398532628,
'each line in file ' &
451 // trim(file%name) //
' must contain at least one data value')
455 times = species_data(1,:)
456 rates = species_data(2,:)
457 if (
allocated(gas_states))
deallocate(gas_states)
458 allocate(gas_states(n_time))
464 if (species == 0)
then
465 call die_msg(806500079,
'unknown species ' &
466 // trim(species_name(i)) //
' in file ' &
470 gas_states(i_time)%mix_rat(species) = species_data(i,i_time)
474 end subroutine spec_file_read_gas_states_times_rates
489 val%mix_rat = val_avg%mix_rat
509 val%mix_rat = val_avg%mix_rat
534 character,
intent(inout) :: buffer(:)
536 integer,
intent(inout) :: position
541 integer :: prev_position
543 prev_position = position
557 character,
intent(inout) :: buffer(:)
559 integer,
intent(inout) :: position
564 integer :: prev_position
566 prev_position = position
592 subroutine gas_state_output_netcdf(gas_state, ncid, gas_data)
597 integer,
intent(in) :: ncid
601 integer :: dimid_gas_species
621 "gas_mixing_ratio", (/ dimid_gas_species /), unit=
"ppb", &
622 long_name=
"mixing ratios of gas species")
624 end subroutine gas_state_output_netcdf
634 integer,
intent(in) :: ncid
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.
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.
Wrapper functions for MPI.
subroutine gas_state_mix(val)
Average val over all processes.
subroutine gas_state_ensure_nonnegative(gas_state)
Set any negative values to zero.
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...
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.
subroutine gas_state_mole_dens_to_ppb(gas_state, env_state)
Convert (mol m^{-3}) to (ppb).
subroutine gas_state_set_size(gas_state, n_spec)
Sets the sizes of the gas state.
The gas_data_t structure and associated subroutines.
subroutine pmc_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from a NetCDF file.
subroutine die_msg(code, error_msg)
Error immediately.
integer function pmc_mpi_rank()
Returns the rank of the current process.
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
integer, parameter dp
Kind of a double precision real number.
Reading formatted text input.
subroutine gas_state_add_scaled(gas_state, gas_state_delta, alpha)
Adds the given gas_state_delta scaled by alpha.
Current environment state.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
An input file with extra data for printing messages.
The gas_state_t structure and associated subroutines.
subroutine pmc_mpi_unpack_gas_state(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine gas_state_reduce_avg(val)
Average val over all processes, with the result only on the root process.
real(kind=dp) function env_state_air_molar_den(env_state)
Air molar density (mol m^{-3}).
subroutine gas_state_sub(gas_state, gas_state_delta)
Subtracts the given gas_state_delta.
subroutine gas_state_scale(gas_state, alpha)
Scale a gas state.
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.
The env_state_t structure and associated subroutines.
Current state of the gas mixing ratios in the system.
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...
Common utility subroutines.
elemental integer function gas_data_n_spec(gas_data)
Return the number of gas species.
subroutine pmc_mpi_pack_gas_state(buffer, position, val)
Packs the given value into the buffer, advancing position.
logical function gas_state_is_allocated(gas_state)
Determine whether the gas_state is correctly allocated.
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...
subroutine gas_state_add(gas_state, gas_state_delta)
Adds the given gas_state_delta.
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
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.
integer function pmc_mpi_pack_size_gas_state(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_mpi_unpack_real_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine gas_state_input_netcdf(gas_state, ncid, gas_data)
Read full state.
integer function find_1d(n, x_vals, x)
Find the position of a real number in an arbitrary 1D array.
subroutine gas_state_zero(gas_state)
Zeros the state.