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