PartMC  2.5.0
env_state.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_env_state module.
7 
8 !> The env_state_t structure and associated subroutines.
10 
11  use pmc_constants
12  use pmc_util
13  use pmc_spec_file
14  use pmc_mpi
15  use pmc_netcdf
16 #ifdef PMC_USE_MPI
17  use mpi
18 #endif
19 
20  !> Current environment state.
21  !!
22  !! All quantities are instantaneous, describing the state at a
23  !! particular instant of time. Constant data and other data not
24  !! associated with the current environment state is stored in
25  !! scenario_t.
27  !> Temperature (K).
28  real(kind=dp) :: temp
29  !> Relative humidity (1).
30  real(kind=dp) :: rel_humid
31  !> Ambient pressure (Pa).
32  real(kind=dp) :: pressure
33  !> Longitude (degrees).
34  real(kind=dp) :: longitude
35  !> Latitude (degrees).
36  real(kind=dp) :: latitude
37  !> Altitude (m).
38  real(kind=dp) :: altitude
39  !> Start time (s since 00:00 UTC on \c start_day).
40  real(kind=dp) :: start_time
41  !> Start day of year (UTC).
42  integer :: start_day
43  !> Time since \c start_time (s).
44  real(kind=dp) :: elapsed_time
45  !> Solar zenith angle (radians from zenith).
46  real(kind=dp) :: solar_zenith_angle
47  !> Box height (m).
48  real(kind=dp) :: height
49  end type env_state_t
50 
51 contains
52 
53 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
54 
55  !> env_state += env_state_delta
56  subroutine env_state_add(env_state, env_state_delta)
57 
58  !> Environment.
59  type(env_state_t), intent(inout) :: env_state
60  !> Increment.
61  type(env_state_t), intent(in) :: env_state_delta
62 
63  env_state%temp = env_state%temp + env_state_delta%temp
64  env_state%rel_humid = env_state%rel_humid + env_state_delta%rel_humid
65  env_state%pressure = env_state%pressure + env_state_delta%pressure
66  env_state%longitude = env_state%longitude + env_state_delta%longitude
67  env_state%latitude = env_state%latitude + env_state_delta%latitude
68  env_state%altitude = env_state%altitude + env_state_delta%altitude
69  env_state%start_time = env_state%start_time + env_state_delta%start_time
70  env_state%start_day = env_state%start_day + env_state_delta%start_day
71  env_state%elapsed_time = env_state%elapsed_time &
72  + env_state_delta%elapsed_time
73  env_state%solar_zenith_angle = env_state%solar_zenith_angle &
74  + env_state_delta%solar_zenith_angle
75  env_state%height = env_state%height + env_state_delta%height
76 
77  end subroutine env_state_add
78 
79 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
80 
81  !> env_state *= alpha
82  subroutine env_state_scale(env_state, alpha)
83 
84  !> Environment.
85  type(env_state_t), intent(inout) :: env_state
86  !> Scale factor.
87  real(kind=dp), intent(in) :: alpha
88 
89  env_state%temp = env_state%temp * alpha
90  env_state%rel_humid = env_state%rel_humid * alpha
91  env_state%pressure = env_state%pressure * alpha
92  env_state%longitude = env_state%longitude * alpha
93  env_state%latitude = env_state%latitude * alpha
94  env_state%altitude = env_state%altitude * alpha
95  env_state%start_time = env_state%start_time * alpha
96  env_state%start_day = nint(real(env_state%start_day, kind=dp) * alpha)
97  env_state%elapsed_time = env_state%elapsed_time * alpha
98  env_state%solar_zenith_angle = env_state%solar_zenith_angle * alpha
99  env_state%height = env_state%height * alpha
100 
101  end subroutine env_state_scale
102 
103 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
104 
105  !> Adds the given water volume to the water vapor and updates all
106  !> environment quantities.
107  subroutine env_state_change_water_volume(env_state, dv)
109  !> Environment state to update.
110  type(env_state_t), intent(inout) :: env_state
111  !> Volume concentration of water added (m^3/m^3).
112  real(kind=dp), intent(in) :: dv
113 
114  real(kind=dp) pmv ! ambient water vapor pressure (Pa)
115  real(kind=dp) mv ! ambient water vapor density (kg m^{-3})
116  ! pmv and mv are related by the factor molec_weight/(R*T)
117  real(kind=dp) dmv ! change of water density (kg m^{-3})
118 
119  dmv = dv * const%water_density
120  pmv = env_state_sat_vapor_pressure(env_state) * env_state%rel_humid
121  mv = const%water_molec_weight / (const%univ_gas_const*env_state%temp) * pmv
122  mv = mv - dmv
123  if (mv < 0d0) then
124  call warn_msg(980320483, "relative humidity tried to go negative")
125  mv = 0d0
126  end if
127  env_state%rel_humid = const%univ_gas_const * env_state%temp &
128  / const%water_molec_weight * mv &
129  / env_state_sat_vapor_pressure(env_state)
130 
131  end subroutine env_state_change_water_volume
132 
133 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
134 
135  !> Computes the current saturation vapor pressure (Pa).
136  real(kind=dp) function env_state_sat_vapor_pressure(env_state)
138  !> Environment state.
139  type(env_state_t), intent(in) :: env_state
140 
141  env_state_sat_vapor_pressure = const%water_eq_vap_press &
142  * 10d0**(7.45d0 * (env_state%temp - const%water_freeze_temp) &
143  / (env_state%temp - 38d0))
144 
145  end function env_state_sat_vapor_pressure
146 
147 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
148 
149  !> Air density (kg m^{-3}).
150  real(kind=dp) function env_state_air_den(env_state)
152  !> Environment state.
153  type(env_state_t), intent(in) :: env_state
154 
155  env_state_air_den = const%air_molec_weight &
156  * env_state_air_molar_den(env_state)
157 
158  end function env_state_air_den
159 
160 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
161 
162  !> Air molar density (mol m^{-3}).
163  real(kind=dp) function env_state_air_molar_den(env_state)
165  !> Environment state.
166  type(env_state_t), intent(in) :: env_state
167 
168  env_state_air_molar_den = env_state%pressure &
169  / (const%univ_gas_const * env_state%temp)
170 
171  end function env_state_air_molar_den
172 
173 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
174 
175  !> Condensation \f$A\f$ parameter.
176  real(kind=dp) function env_state_a(env_state)
178  !> Environment state.
179  type(env_state_t), intent(in) :: env_state
180 
181  env_state_a = 4d0 * const%water_surf_eng * const%water_molec_weight &
182  / (const%univ_gas_const * env_state%temp * const%water_density)
183 
184  end function env_state_a
185 
186 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
187 
188  !> Convert (ppb) to (molecules m^{-3}).
189  real(kind=dp) function env_state_ppb_to_conc(env_state, ppb)
191  !> Environment state.
192  type(env_state_t), intent(in) :: env_state
193  !> Mixing ratio (ppb).
194  real(kind=dp), intent(in) :: ppb
195 
196  env_state_ppb_to_conc = ppb / 1d9 * env_state_air_molar_den(env_state) &
197  * const%avagadro
198 
199  end function env_state_ppb_to_conc
200 
201 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
202 
203  !> Convert (molecules m^{-3}) to (ppb).
204  real(kind=dp) function env_state_conc_to_ppb(env_state, conc)
206  !> Environment state.
207  type(env_state_t), intent(in) :: env_state
208  !> Concentration (molecules m^{-3}).
209  real(kind=dp), intent(in) :: conc
210 
211  env_state_conc_to_ppb = conc * 1d9 / env_state_air_molar_den(env_state) &
212  / const%avagadro
213 
214  end function env_state_conc_to_ppb
215 
216 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
217 
218  !> Read environment specification from a spec file.
219  subroutine spec_file_read_env_state(file, env_state)
220 
221  !> Spec file.
222  type(spec_file_t), intent(inout) :: file
223  !> Environment data.
224  type(env_state_t), intent(inout) :: env_state
225 
226  !> \page input_format_env_state Input File Format: Environment State
227  !!
228  !! The environment parameters are divided into those specified at
229  !! the start of the simulation and then either held constant or
230  !! computed for the rest of the simulation, and those parameters
231  !! given as prescribed profiles for the entire simulation
232  !! duration. The variables below are for the first type --- for
233  !! the prescribed profiles see \ref input_format_scenario.
234  !!
235  !! The environment state is specified by the parameters:
236  !! - \b rel_humidity (real, dimensionless): the relative humidity
237  !! (0 is completely unsaturated and 1 is fully saturated)
238  !! - \b latitude (real, unit degrees_north): the latitude of the
239  !! simulation location
240  !! - \b longitude (real, unit degrees_east): the longitude of the
241  !! simulation location
242  !! - \b altitude (real, unit m): the altitude of the simulation
243  !! location
244  !! - \b start_time (real, unit s): the time-of-day of the start of
245  !! the simulation (in seconds past midnight)
246  !! - \b start_day (integer): the day-of-year of the start of the
247  !! simulation (starting from 1 on the first day of the year)
248  !!
249  !! See also:
250  !! - \ref spec_file_format --- the input file text format
251  !! - \ref output_format_env_state --- the corresponding output
252  !! format
253  !! - \ref input_format_scenario --- the prescribed profiles of
254  !! other environment data
255 
256  call spec_file_read_real(file, 'rel_humidity', env_state%rel_humid)
257  call spec_file_read_real(file, 'latitude', env_state%latitude)
258  call spec_file_read_real(file, 'longitude', env_state%longitude)
259  call spec_file_read_real(file, 'altitude', env_state%altitude)
260  call spec_file_read_real(file, 'start_time', env_state%start_time)
261  call spec_file_read_integer(file, 'start_day', env_state%start_day)
262 
263  end subroutine spec_file_read_env_state
264 
265 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
266 
267  !> Average val over all processes.
268  subroutine env_state_mix(val)
270  !> Value to average.
271  type(env_state_t), intent(inout) :: val
272 
273 #ifdef PMC_USE_MPI
274  type(env_state_t) :: val_avg
275 
276  call pmc_mpi_allreduce_average_real(val%temp, val_avg%temp)
277  call pmc_mpi_allreduce_average_real(val%rel_humid, val_avg%rel_humid)
278  call pmc_mpi_allreduce_average_real(val%pressure, val_avg%pressure)
279  val%temp = val_avg%temp
280  val%rel_humid = val_avg%rel_humid
281  val%pressure = val_avg%pressure
282 #endif
283 
284  end subroutine env_state_mix
285 
286 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
287 
288  !> Average val over all processes, with the result only on the root
289  !> process.
290  subroutine env_state_reduce_avg(val)
292  !> Value to average.
293  type(env_state_t), intent(inout) :: val
294 
295 #ifdef PMC_USE_MPI
296  type(env_state_t) :: val_avg
297 
298  call pmc_mpi_reduce_avg_real(val%temp, val_avg%temp)
299  call pmc_mpi_reduce_avg_real(val%rel_humid, val_avg%rel_humid)
300  call pmc_mpi_reduce_avg_real(val%pressure, val_avg%pressure)
301  if (pmc_mpi_rank() == 0) then
302  val%temp = val_avg%temp
303  val%rel_humid = val_avg%rel_humid
304  val%pressure = val_avg%pressure
305  end if
306 #endif
307 
308  end subroutine env_state_reduce_avg
309 
310 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
311 
312  !> Determines the number of bytes required to pack the given value.
313  integer function pmc_mpi_pack_size_env_state(val)
315  !> Value to pack.
316  type(env_state_t), intent(in) :: val
317 
319  pmc_mpi_pack_size_real(val%temp) &
320  + pmc_mpi_pack_size_real(val%rel_humid) &
321  + pmc_mpi_pack_size_real(val%pressure) &
322  + pmc_mpi_pack_size_real(val%longitude) &
323  + pmc_mpi_pack_size_real(val%latitude) &
324  + pmc_mpi_pack_size_real(val%altitude) &
325  + pmc_mpi_pack_size_real(val%start_time) &
326  + pmc_mpi_pack_size_integer(val%start_day) &
327  + pmc_mpi_pack_size_real(val%elapsed_time) &
328  + pmc_mpi_pack_size_real(val%solar_zenith_angle) &
329  + pmc_mpi_pack_size_real(val%height)
330 
331  end function pmc_mpi_pack_size_env_state
332 
333 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
334 
335  !> Packs the given value into the buffer, advancing position.
336  subroutine pmc_mpi_pack_env_state(buffer, position, val)
338  !> Memory buffer.
339  character, intent(inout) :: buffer(:)
340  !> Current buffer position.
341  integer, intent(inout) :: position
342  !> Value to pack.
343  type(env_state_t), intent(in) :: val
344 
345 #ifdef PMC_USE_MPI
346  integer :: prev_position
347 
348  prev_position = position
349  call pmc_mpi_pack_real(buffer, position, val%temp)
350  call pmc_mpi_pack_real(buffer, position, val%rel_humid)
351  call pmc_mpi_pack_real(buffer, position, val%pressure)
352  call pmc_mpi_pack_real(buffer, position, val%longitude)
353  call pmc_mpi_pack_real(buffer, position, val%latitude)
354  call pmc_mpi_pack_real(buffer, position, val%altitude)
355  call pmc_mpi_pack_real(buffer, position, val%start_time)
356  call pmc_mpi_pack_integer(buffer, position, val%start_day)
357  call pmc_mpi_pack_real(buffer, position, val%elapsed_time)
358  call pmc_mpi_pack_real(buffer, position, val%solar_zenith_angle)
359  call pmc_mpi_pack_real(buffer, position, val%height)
360  call assert(464101191, &
361  position - prev_position <= pmc_mpi_pack_size_env_state(val))
362 #endif
363 
364  end subroutine pmc_mpi_pack_env_state
365 
366 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
367 
368  !> Unpacks the given value from the buffer, advancing position.
369  subroutine pmc_mpi_unpack_env_state(buffer, position, val)
371  !> Memory buffer.
372  character, intent(inout) :: buffer(:)
373  !> Current buffer position.
374  integer, intent(inout) :: position
375  !> Value to pack.
376  type(env_state_t), intent(inout) :: val
377 
378 #ifdef PMC_USE_MPI
379  integer :: prev_position
380 
381  prev_position = position
382  call pmc_mpi_unpack_real(buffer, position, val%temp)
383  call pmc_mpi_unpack_real(buffer, position, val%rel_humid)
384  call pmc_mpi_unpack_real(buffer, position, val%pressure)
385  call pmc_mpi_unpack_real(buffer, position, val%longitude)
386  call pmc_mpi_unpack_real(buffer, position, val%latitude)
387  call pmc_mpi_unpack_real(buffer, position, val%altitude)
388  call pmc_mpi_unpack_real(buffer, position, val%start_time)
389  call pmc_mpi_unpack_integer(buffer, position, val%start_day)
390  call pmc_mpi_unpack_real(buffer, position, val%elapsed_time)
391  call pmc_mpi_unpack_real(buffer, position, val%solar_zenith_angle)
392  call pmc_mpi_unpack_real(buffer, position, val%height)
393  call assert(205696745, &
394  position - prev_position <= pmc_mpi_pack_size_env_state(val))
395 #endif
396 
397  end subroutine pmc_mpi_unpack_env_state
398 
399 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
400 
401  !> Computes the average of val across all processes, storing the
402  !> result in val_avg on the root process.
403  subroutine pmc_mpi_reduce_avg_env_state(val, val_avg)
405  !> Value to average.
406  type(env_state_t), intent(in) :: val
407  !> Result.
408  type(env_state_t), intent(inout) :: val_avg
409 
410  val_avg = val
411  call pmc_mpi_reduce_avg_real(val%temp, val_avg%temp)
412  call pmc_mpi_reduce_avg_real(val%rel_humid, val_avg%rel_humid)
413  call pmc_mpi_reduce_avg_real(val%pressure, val_avg%pressure)
414 
415  end subroutine pmc_mpi_reduce_avg_env_state
416 
417 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
418 
419  !> Write full state.
420  subroutine env_state_output_netcdf(env_state, ncid)
421 
422  !> Environment state to write.
423  type(env_state_t), intent(in) :: env_state
424  !> NetCDF file ID, in data mode.
425  integer, intent(in) :: ncid
426 
427  !> \page output_format_env_state Output File Format: Environment State
428  !!
429  !! The environment state NetCDF variables are:
430  !! - \b temperature (unit K): current air temperature
431  !! - \b relative_humidity (dimensionless): current air
432  !! relative humidity (value of 1 means completely saturated)
433  !! - \b pressure (unit Pa): current air pressure
434  !! - \b longitude (unit degrees_east): longitude of simulation location
435  !! - \b latitude (unit degrees_north): latitude of simulation location
436  !! - \b altitude (unit m): altitude of simulation location
437  !! - \b start_time_of_day (unit s): time-of-day of the
438  !! simulation start measured in seconds after midnight UTC
439  !! - \b start_day_of_year: day-in-year number of the simulation start
440  !! (starting from 1 on the first day of the year)
441  !! - \b elapsed_time (unit s): elapsed time since the simulation start
442  !! - \b solar_zenith_angle (unit radians): current angle from
443  !! the zenith to the sun
444  !! - \b height (unit m): current boundary layer mixing height
445  !!
446  !! See also:
447  !! - \ref input_format_env_state and \ref input_format_scenario
448  !! --- the corresponding input formats
449 
450  call pmc_nc_write_real(ncid, env_state%temp, "temperature", unit="K", &
451  standard_name="air_temperature")
452  call pmc_nc_write_real(ncid, env_state%rel_humid, &
453  "relative_humidity", unit="1", standard_name="relative_humidity")
454  call pmc_nc_write_real(ncid, env_state%pressure, "pressure", unit="Pa", &
455  standard_name="air_pressure")
456  call pmc_nc_write_real(ncid, env_state%longitude, "longitude", &
457  unit="degree_east", standard_name="longitude")
458  call pmc_nc_write_real(ncid, env_state%latitude, "latitude", &
459  unit="degree_north", standard_name="latitude")
460  call pmc_nc_write_real(ncid, env_state%altitude, "altitude", unit="m", &
461  standard_name="altitude")
462  call pmc_nc_write_real(ncid, env_state%start_time, &
463  "start_time_of_day", unit="s", description="time-of-day of " &
464  // "simulation start in seconds since midnight")
465  call pmc_nc_write_integer(ncid, env_state%start_day, &
466  "start_day_of_year", &
467  description="day-of-year number of simulation start")
468  call pmc_nc_write_real(ncid, env_state%elapsed_time, "elapsed_time", &
469  unit="s", description="elapsed time since simulation start")
470  call pmc_nc_write_real(ncid, env_state%solar_zenith_angle, &
471  "solar_zenith_angle", unit="radian", &
472  description="current angle from the zenith to the sun")
473  call pmc_nc_write_real(ncid, env_state%height, "height", unit="m", &
474  long_name="boundary layer mixing height")
475 
476  end subroutine env_state_output_netcdf
477 
478 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
479 
480  !> Read full state.
481  subroutine env_state_input_netcdf(env_state, ncid)
483  !> Environment state to read.
484  type(env_state_t), intent(inout) :: env_state
485  !> NetCDF file ID, in data mode.
486  integer, intent(in) :: ncid
487 
488  call pmc_nc_read_real(ncid, env_state%temp, "temperature")
489  call pmc_nc_read_real(ncid, env_state%rel_humid, "relative_humidity")
490  call pmc_nc_read_real(ncid, env_state%pressure, "pressure")
491  call pmc_nc_read_real(ncid, env_state%longitude, "longitude")
492  call pmc_nc_read_real(ncid, env_state%latitude, "latitude")
493  call pmc_nc_read_real(ncid, env_state%altitude, "altitude")
494  call pmc_nc_read_real(ncid, env_state%start_time, &
495  "start_time_of_day")
496  call pmc_nc_read_integer(ncid, env_state%start_day, &
497  "start_day_of_year")
498  call pmc_nc_read_real(ncid, env_state%elapsed_time, "elapsed_time")
499  call pmc_nc_read_real(ncid, env_state%solar_zenith_angle, &
500  "solar_zenith_angle")
501  call pmc_nc_read_real(ncid, env_state%height, "height")
502 
503  end subroutine env_state_input_netcdf
504 
505 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
506 
507 end module pmc_env_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
real(kind=dp) function env_state_a(env_state)
Condensation parameter.
Definition: env_state.F90:177
integer function pmc_mpi_rank()
Returns the rank of the current process.
Definition: mpi.F90:117
subroutine env_state_scale(env_state, alpha)
env_state *= alpha
Definition: env_state.F90:83
subroutine warn_msg(code, warning_msg, already_warned)
Prints a warning message.
Definition: util.F90:37
subroutine env_state_reduce_avg(val)
Average val over all processes, with the result only on the root process.
Definition: env_state.F90:291
The env_state_t structure and associated subroutines.
Definition: env_state.F90:9
integer function pmc_mpi_pack_size_env_state(val)
Determines the number of bytes required to pack the given value.
Definition: env_state.F90:314
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:103
subroutine spec_file_read_real(file, name, var)
Read a real number from a spec file that must have the given name.
Definition: spec_file.F90:562
subroutine pmc_nc_read_integer(ncid, var, name, must_be_present)
Read a single integer from a NetCDF file.
Definition: netcdf.F90:184
subroutine spec_file_read_integer(file, name, var)
Read an integer from a spec file that must have the given name.
Definition: spec_file.F90:540
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:365
Current environment state.
Definition: env_state.F90:26
subroutine pmc_nc_write_integer(ncid, var, name, unit, long_name, standard_name, description)
Write a single integer to a NetCDF file.
Definition: netcdf.F90:459
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:586
subroutine env_state_change_water_volume(env_state, dv)
Adds the given water volume to the water vapor and updates all environment quantities.
Definition: env_state.F90:108
Physical constants.
Definition: constants.F90:9
subroutine pmc_nc_read_real(ncid, var, name, must_be_present)
Read a single real from a NetCDF file.
Definition: netcdf.F90:149
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:561
real(kind=dp) function env_state_air_den(env_state)
Air density (kg m^{-3}).
Definition: env_state.F90:151
Reading formatted text input.
Definition: spec_file.F90:43
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:843
subroutine pmc_mpi_reduce_avg_env_state(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: env_state.F90:404
subroutine pmc_nc_write_real(ncid, var, name, unit, long_name, standard_name, description)
Write a single real to a NetCDF file.
Definition: netcdf.F90:426
subroutine env_state_input_netcdf(env_state, ncid)
Read full state.
Definition: env_state.F90:482
real(kind=dp) function env_state_conc_to_ppb(env_state, conc)
Convert (molecules m^{-3}) to (ppb).
Definition: env_state.F90:205
subroutine pmc_mpi_unpack_env_state(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: env_state.F90:370
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:345
subroutine pmc_mpi_allreduce_average_real(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on all processes...
Definition: mpi.F90:1281
subroutine pmc_mpi_pack_env_state(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: env_state.F90:337
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:818
subroutine env_state_mix(val)
Average val over all processes.
Definition: env_state.F90:269
subroutine pmc_mpi_reduce_avg_real(val, val_avg)
Computes the average of val across all processes, storing the result in val_avg on the root process...
Definition: mpi.F90:1076
type(const_t), save const
Fixed variable for accessing the constant&#39;s values.
Definition: constants.F90:70
Common utility subroutines.
Definition: util.F90:9
Wrapper functions for MPI.
Definition: mpi.F90:13
real(kind=dp) function env_state_ppb_to_conc(env_state, ppb)
Convert (ppb) to (molecules m^{-3}).
Definition: env_state.F90:190
subroutine env_state_add(env_state, env_state_delta)
env_state += env_state_delta
Definition: env_state.F90:57
real(kind=dp) function env_state_sat_vapor_pressure(env_state)
Computes the current saturation vapor pressure (Pa).
Definition: env_state.F90:137