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