PartMC  2.4.0
run_part.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2017 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_run_part module.
7 
8 !> Monte Carlo simulation.
10 
11  use pmc_util
12  use pmc_aero_state
13  use pmc_scenario
14  use pmc_env_state
15  use pmc_aero_data
16  use pmc_gas_data
17  use pmc_gas_state
18  use pmc_output
19  use pmc_mosaic
20  use pmc_coagulation
22  use pmc_coag_kernel
23  use pmc_nucleate
24  use pmc_mpi
25 #ifdef PMC_USE_SUNDIALS
26  use pmc_condense
27 #endif
28 #ifdef PMC_USE_MPI
29  use mpi
30 #endif
31 
32  !> Type code for undefined or invalid parallel coagulation method.
33  integer, parameter :: parallel_coag_type_invalid = 0
34  !> Type code for local parallel coagulation.
35  integer, parameter :: parallel_coag_type_local = 1
36  !> Type code for distributed parallel coagulation.
37  integer, parameter :: parallel_coag_type_dist = 2
38 
39  !> Options controlling the execution of run_part().
41  !> Final time (s).
42  real(kind=dp) :: t_max
43  !> Output interval (0 disables) (s).
44  real(kind=dp) :: t_output
45  !> Progress interval (0 disables) (s).
46  real(kind=dp) :: t_progress
47  !> Timestep for coagulation.
48  real(kind=dp) :: del_t
49  !> Prefix for output files.
50  character(len=300) :: output_prefix
51  !> Type of coagulation kernel.
52  integer :: coag_kernel_type
53  !> Type of nucleation.
54  integer :: nucleate_type
55  !> Source of nucleation.
56  integer :: nucleate_source
57  !> Whether to do coagulation.
58  logical :: do_coagulation
59  !> Whether to do nucleation.
60  logical :: do_nucleation
61  !> Allow doubling if needed.
62  logical :: allow_doubling
63  !> Allow halving if needed.
64  logical :: allow_halving
65  !> Whether to do condensation.
66  logical :: do_condensation
67  !> Whether to do MOSAIC.
68  logical :: do_mosaic
69  !> Whether to compute optical properties.
70  logical :: do_optical
71  !> Whether to have explicitly selected weighting.
72  logical :: do_select_weighting
73  !> Type of particle weighting scheme.
74  integer :: weighting_type
75  !> Weighting exponent for power weighting scheme.
76  real(kind=dp) :: weighting_exponent
77  !> Repeat number of run.
78  integer :: i_repeat
79  !> Total number of repeats.
80  integer :: n_repeat
81  !> Cpu_time() of start.
82  real(kind=dp) :: t_wall_start
83  !> Whether to record particle removal information.
84  logical :: record_removals
85  !> Whether to run in parallel.
86  logical :: do_parallel
87  !> Parallel output type.
88  integer :: output_type
89  !> Mixing timescale between processes (s).
90  real(kind=dp) :: mix_timescale
91  !> Whether to average gases each timestep.
92  logical :: gas_average
93  !> Whether to average environment each timestep.
94  logical :: env_average
95  !> Parallel coagulation method type.
96  integer :: parallel_coag_type
97  !> UUID for this simulation.
98  character(len=PMC_UUID_LEN) :: uuid
99  end type run_part_opt_t
100 
101 contains
102 
103 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
104 
105  !> Do a particle-resolved Monte Carlo simulation.
106  subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, &
107  gas_state, run_part_opt)
109  !> Environment state.
110  type(scenario_t), intent(in) :: scenario
111  !> Environment state.
112  type(env_state_t), intent(inout) :: env_state
113  !> Aerosol data.
114  type(aero_data_t), intent(in) :: aero_data
115  !> Aerosol state.
116  type(aero_state_t), intent(inout) :: aero_state
117  !> Gas data.
118  type(gas_data_t), intent(in) :: gas_data
119  !> Gas state.
120  type(gas_state_t), intent(inout) :: gas_state
121  !> Monte Carlo options.
122  type(run_part_opt_t), intent(in) :: run_part_opt
123 
124  real(kind=dp) :: time, pre_time, pre_del_t, prop_done
125  real(kind=dp) :: last_output_time, last_progress_time
126  integer :: rank, n_proc, pre_index, ncid
127  integer :: pre_i_repeat
128  integer :: n_samp, n_coag, n_emit, n_dil_in, n_dil_out, n_nuc
129  integer :: progress_n_samp, progress_n_coag
130  integer :: progress_n_emit, progress_n_dil_in, progress_n_dil_out
131  integer :: progress_n_nuc, n_part_before
132  integer :: global_n_part, global_n_samp, global_n_coag
133  integer :: global_n_emit, global_n_dil_in, global_n_dil_out
134  integer :: global_n_nuc
135  logical :: do_output, do_state, do_state_netcdf, do_progress, did_coag
136  real(kind=dp) :: t_start, t_wall_now, t_wall_elapsed, t_wall_remain
137  type(env_state_t) :: old_env_state
138  integer :: n_time, i_time, i_time_start, pre_i_time
139  integer :: i_state, i_state_netcdf, i_output
140 
141  rank = pmc_mpi_rank()
142  n_proc = pmc_mpi_size()
143 
144  i_time = 0
145  i_output = 1
146  i_state = 1
147  i_state_netcdf = 1
148  time = 0d0
149  progress_n_samp = 0
150  progress_n_coag = 0
151  progress_n_emit = 0
152  progress_n_dil_in = 0
153  progress_n_dil_out = 0
154  progress_n_nuc = 0
155 
156  call check_time_multiple("t_max", run_part_opt%t_max, &
157  "del_t", run_part_opt%del_t)
158  call check_time_multiple("t_output", run_part_opt%t_output, &
159  "del_t", run_part_opt%del_t)
160  call check_time_multiple("t_progress", run_part_opt%t_progress, &
161  "del_t", run_part_opt%del_t)
162 
163  if (run_part_opt%do_mosaic) then
164  call mosaic_init(env_state, aero_data, run_part_opt%del_t, &
165  run_part_opt%do_optical)
166  if (run_part_opt%do_optical) then
167  call mosaic_aero_optical_init(env_state, aero_data, &
168  aero_state, gas_data, gas_state)
169  end if
170  end if
171 
172  if (run_part_opt%t_output > 0d0) then
173  call output_state(run_part_opt%output_prefix, &
174  run_part_opt%output_type, aero_data, aero_state, gas_data, &
175  gas_state, env_state, i_state, time, run_part_opt%del_t, &
176  run_part_opt%i_repeat, run_part_opt%record_removals, &
177  run_part_opt%do_optical, run_part_opt%uuid)
178  call aero_info_array_zero(aero_state%aero_info_array)
179  end if
180 
181  call aero_state_rebalance(aero_state, aero_data, &
182  run_part_opt%allow_doubling, &
183  run_part_opt%allow_halving, initial_state_warning=.true.)
184 
185  t_start = env_state%elapsed_time
186  last_output_time = time
187  last_progress_time = time
188  n_time = nint(run_part_opt%t_max / run_part_opt%del_t)
189  i_time_start = nint(time / run_part_opt%del_t) + 1
190 
191  global_n_part = aero_state_total_particles_all_procs(aero_state)
192  if (rank == 0) then
193  ! progress only printed from root process
194  if (run_part_opt%i_repeat == 1) then
195  t_wall_elapsed = 0d0
196  t_wall_remain = 0d0
197  else
198  call cpu_time(t_wall_now)
199  prop_done = real(run_part_opt%i_repeat - 1, kind=dp) &
200  / real(run_part_opt%n_repeat, kind=dp)
201  t_wall_elapsed = t_wall_now - run_part_opt%t_wall_start
202  t_wall_remain = (1d0 - prop_done) / prop_done &
203  * t_wall_elapsed
204  end if
205  call print_part_progress(run_part_opt%i_repeat, time, &
206  global_n_part, 0, 0, 0, 0, 0, t_wall_elapsed, t_wall_remain)
207  end if
208 
209  do i_time = i_time_start,n_time
210 
211  time = real(i_time, kind=dp) * run_part_opt%del_t
212 
213  old_env_state = env_state
214  call scenario_update_env_state(scenario, env_state, time + t_start)
215 
216  if (run_part_opt%do_nucleation) then
217  n_part_before = aero_state_total_particles(aero_state)
218  call nucleate(run_part_opt%nucleate_type, &
219  run_part_opt%nucleate_source, env_state, gas_data, aero_data, &
220  aero_state, gas_state, run_part_opt%del_t, &
221  run_part_opt%allow_doubling, run_part_opt%allow_halving)
222  n_nuc = aero_state_total_particles(aero_state) &
223  - n_part_before
224  progress_n_nuc = progress_n_nuc + n_nuc
225  end if
226 
227  if (run_part_opt%do_coagulation) then
228  if (run_part_opt%parallel_coag_type &
229  == parallel_coag_type_local) then
230  call mc_coag(run_part_opt%coag_kernel_type, env_state, &
231  aero_data, aero_state, run_part_opt%del_t, n_samp, n_coag)
232  elseif (run_part_opt%parallel_coag_type &
233  == parallel_coag_type_dist) then
234  call mc_coag_dist(run_part_opt%coag_kernel_type, env_state, &
235  aero_data, aero_state, run_part_opt%del_t, n_samp, n_coag)
236  else
237  call die_msg(323011762, "unknown parallel coagulation type: " &
238  // trim(integer_to_string(run_part_opt%parallel_coag_type)))
239  end if
240  progress_n_samp = progress_n_samp + n_samp
241  progress_n_coag = progress_n_coag + n_coag
242  end if
243 
244 #ifdef PMC_USE_SUNDIALS
245  if (run_part_opt%do_condensation) then
246  call condense_particles(aero_state, aero_data, old_env_state, &
247  env_state, run_part_opt%del_t)
248  end if
249 #endif
250 
251  call scenario_update_gas_state(scenario, run_part_opt%del_t, &
252  env_state, old_env_state, gas_data, gas_state)
253  call scenario_update_aero_state(scenario, run_part_opt%del_t, &
254  env_state, old_env_state, aero_data, aero_state, n_emit, &
255  n_dil_in, n_dil_out, run_part_opt%allow_doubling, &
256  run_part_opt%allow_halving)
257  progress_n_emit = progress_n_emit + n_emit
258  progress_n_dil_in = progress_n_dil_in + n_dil_in
259  progress_n_dil_out = progress_n_dil_out + n_dil_out
260 
261  if (run_part_opt%do_mosaic) then
262  call mosaic_timestep(env_state, aero_data, aero_state, gas_data, &
263  gas_state, run_part_opt%do_optical)
264  end if
265 
266  if (run_part_opt%mix_timescale > 0d0) then
267  call aero_state_mix(aero_state, run_part_opt%del_t, &
268  run_part_opt%mix_timescale, aero_data)
269  end if
270  if (run_part_opt%gas_average) then
271  call gas_state_mix(gas_state)
272  end if
273  if (run_part_opt%gas_average) then
274  call env_state_mix(env_state)
275  end if
276 
277  call aero_state_rebalance(aero_state, aero_data, &
278  run_part_opt%allow_doubling, &
279  run_part_opt%allow_halving, initial_state_warning=.false.)
280 
281  ! DEBUG: enable to check consistency
282  ! call aero_state_check(aero_state, aero_data)
283  ! DEBUG: end
284 
285  if (run_part_opt%t_output > 0d0) then
286  call check_event(time, run_part_opt%del_t, run_part_opt%t_output, &
287  last_output_time, do_output)
288  if (do_output) then
289  i_output = i_output + 1
290  call output_state(run_part_opt%output_prefix, &
291  run_part_opt%output_type, aero_data, aero_state, gas_data, &
292  gas_state, env_state, i_output, time, run_part_opt%del_t, &
293  run_part_opt%i_repeat, run_part_opt%record_removals, &
294  run_part_opt%do_optical, run_part_opt%uuid)
295  call aero_info_array_zero(aero_state%aero_info_array)
296  end if
297  end if
298 
299  if (.not. run_part_opt%record_removals) then
300  ! If we are not recording removals then we can zero them as
301  ! often as possible to minimize the cost of maintaining
302  ! them.
303  call aero_info_array_zero(aero_state%aero_info_array)
304  end if
305 
306  if (run_part_opt%t_progress > 0d0) then
307  call check_event(time, run_part_opt%del_t, &
308  run_part_opt%t_progress, last_progress_time, do_progress)
309  if (do_progress) then
310  global_n_part = aero_state_total_particles_all_procs(aero_state)
311  call pmc_mpi_reduce_sum_integer(progress_n_samp, global_n_samp)
312  call pmc_mpi_reduce_sum_integer(progress_n_coag, global_n_coag)
313  call pmc_mpi_reduce_sum_integer(progress_n_emit, global_n_emit)
314  call pmc_mpi_reduce_sum_integer(progress_n_dil_in, &
315  global_n_dil_in)
316  call pmc_mpi_reduce_sum_integer(progress_n_dil_out, &
317  global_n_dil_out)
318  call pmc_mpi_reduce_sum_integer(progress_n_nuc, global_n_nuc)
319  if (rank == 0) then
320  ! progress only printed from root process
321  call cpu_time(t_wall_now)
322  prop_done = (real(run_part_opt%i_repeat - 1, kind=dp) &
323  + time / run_part_opt%t_max) &
324  / real(run_part_opt%n_repeat, kind=dp)
325  t_wall_elapsed = t_wall_now - run_part_opt%t_wall_start
326  t_wall_remain = (1d0 - prop_done) / prop_done &
327  * t_wall_elapsed
328  call print_part_progress(run_part_opt%i_repeat, time, &
329  global_n_part, global_n_coag, global_n_emit, &
330  global_n_dil_in, global_n_dil_out, global_n_nuc, &
331  t_wall_elapsed, t_wall_remain)
332  end if
333  ! reset counters so they show information since last
334  ! progress display
335  progress_n_samp = 0
336  progress_n_coag = 0
337  progress_n_emit = 0
338  progress_n_dil_in = 0
339  progress_n_dil_out = 0
340  progress_n_nuc = 0
341  end if
342  end if
343 
344  end do
345 
346  if (run_part_opt%do_mosaic) then
347  call mosaic_cleanup()
348  end if
349 
350  end subroutine run_part
351 
352 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
353 
354  !> Print the current simulation progress to the screen.
355  subroutine print_part_progress(i_repeat, t_sim_elapsed, n_part, n_coag, &
356  n_emit, n_dil_in, n_dil_out, n_nuc, t_wall_elapsed, t_wall_remain)
358  !> Repeat number of simulation.
359  integer, intent(in) :: i_repeat
360  !> Elapsed simulation time (s).
361  real(kind=dp), intent(in) :: t_sim_elapsed
362  !> Number of particles.
363  integer, intent(in) :: n_part
364  !> Number of coagulated particles since last progress printing.
365  integer, intent(in) :: n_coag
366  !> Number of emitted particles since last progress printing.
367  integer, intent(in) :: n_emit
368  !> Number of diluted-in particles since last progress printing.
369  integer, intent(in) :: n_dil_in
370  !> Number of diluted-out particles since last progress printing.
371  integer, intent(in) :: n_dil_out
372  !> Number of nucleated particles since last progress printing.
373  integer, intent(in) :: n_nuc
374  !> Elapsed wall time (s).
375  real(kind=dp), intent(in) :: t_wall_elapsed
376  !> Estimated remaining wall time (s).
377  real(kind=dp), intent(in) :: t_wall_remain
378 
379  write(*,'(a6,a1,a6,a1,a7,a1,a7,a1,a7,a1,a8,a1,a9,a1,a7,a1,a6,a1,a6)') &
380  "repeat", " ", "t_sim", " ", "n_part", " ", "n_coag", " ", &
381  "n_emit", " ", "n_dil_in", " ", "n_dil_out", " ", "n_nuc", " ", &
382  "t_wall", " ", "t_rem"
383  write(*,'(a6,a1,a6,a1,a7,a1,a7,a1,a7,a1,a8,a1,a9,a1,a7,a1,a6,a1,a6)') &
384  trim(integer_to_string_max_len(i_repeat, 6)), " ", &
385  trim(time_to_string_max_len(t_sim_elapsed, 6)), " ", &
386  trim(integer_to_string_max_len(n_part, 7)), " ", &
387  trim(integer_to_string_max_len(n_coag, 7)), " ", &
388  trim(integer_to_string_max_len(n_emit, 7)), " ", &
389  trim(integer_to_string_max_len(n_dil_in, 7)), " ", &
390  trim(integer_to_string_max_len(n_dil_out, 7)), " ", &
391  trim(integer_to_string_max_len(n_nuc, 7)), " ", &
392  trim(time_to_string_max_len(t_wall_elapsed, 6)), " ", &
393  trim(time_to_string_max_len(t_wall_remain, 6))
394 
395  end subroutine print_part_progress
396 
397 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
398 
399  !> Determines the number of bytes required to pack the given value.
400  integer function pmc_mpi_pack_size_run_part_opt(val)
402  !> Value to pack.
403  type(run_part_opt_t), intent(in) :: val
404 
406  pmc_mpi_pack_size_real(val%t_max) &
407  + pmc_mpi_pack_size_real(val%t_output) &
408  + pmc_mpi_pack_size_real(val%t_progress) &
409  + pmc_mpi_pack_size_real(val%del_t) &
410  + pmc_mpi_pack_size_string(val%output_prefix) &
411  + pmc_mpi_pack_size_integer(val%coag_kernel_type) &
412  + pmc_mpi_pack_size_integer(val%nucleate_type) &
413  + pmc_mpi_pack_size_integer(val%nucleate_source) &
414  + pmc_mpi_pack_size_logical(val%do_coagulation) &
415  + pmc_mpi_pack_size_logical(val%do_nucleation) &
416  + pmc_mpi_pack_size_logical(val%allow_doubling) &
417  + pmc_mpi_pack_size_logical(val%allow_halving) &
418  + pmc_mpi_pack_size_logical(val%do_condensation) &
419  + pmc_mpi_pack_size_logical(val%do_mosaic) &
420  + pmc_mpi_pack_size_logical(val%do_optical) &
421  + pmc_mpi_pack_size_logical(val%do_select_weighting) &
422  + pmc_mpi_pack_size_integer(val%weighting_type) &
423  + pmc_mpi_pack_size_real(val%weighting_exponent) &
424  + pmc_mpi_pack_size_integer(val%i_repeat) &
425  + pmc_mpi_pack_size_integer(val%n_repeat) &
426  + pmc_mpi_pack_size_real(val%t_wall_start) &
427  + pmc_mpi_pack_size_logical(val%record_removals) &
428  + pmc_mpi_pack_size_logical(val%do_parallel) &
429  + pmc_mpi_pack_size_integer(val%output_type) &
430  + pmc_mpi_pack_size_real(val%mix_timescale) &
431  + pmc_mpi_pack_size_logical(val%gas_average) &
432  + pmc_mpi_pack_size_logical(val%env_average) &
433  + pmc_mpi_pack_size_integer(val%parallel_coag_type) &
434  + pmc_mpi_pack_size_string(val%uuid)
435 
436  end function pmc_mpi_pack_size_run_part_opt
437 
438 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
439 
440  !> Packs the given value into the buffer, advancing position.
441  subroutine pmc_mpi_pack_run_part_opt(buffer, position, val)
443  !> Memory buffer.
444  character, intent(inout) :: buffer(:)
445  !> Current buffer position.
446  integer, intent(inout) :: position
447  !> Value to pack.
448  type(run_part_opt_t), intent(in) :: val
449 
450 #ifdef PMC_USE_MPI
451  integer :: prev_position
452 
453  prev_position = position
454  call pmc_mpi_pack_real(buffer, position, val%t_max)
455  call pmc_mpi_pack_real(buffer, position, val%t_output)
456  call pmc_mpi_pack_real(buffer, position, val%t_progress)
457  call pmc_mpi_pack_real(buffer, position, val%del_t)
458  call pmc_mpi_pack_string(buffer, position, val%output_prefix)
459  call pmc_mpi_pack_integer(buffer, position, val%coag_kernel_type)
460  call pmc_mpi_pack_integer(buffer, position, val%nucleate_type)
461  call pmc_mpi_pack_integer(buffer, position, val%nucleate_source)
462  call pmc_mpi_pack_logical(buffer, position, val%do_coagulation)
463  call pmc_mpi_pack_logical(buffer, position, val%do_nucleation)
464  call pmc_mpi_pack_logical(buffer, position, val%allow_doubling)
465  call pmc_mpi_pack_logical(buffer, position, val%allow_halving)
466  call pmc_mpi_pack_logical(buffer, position, val%do_condensation)
467  call pmc_mpi_pack_logical(buffer, position, val%do_mosaic)
468  call pmc_mpi_pack_logical(buffer, position, val%do_optical)
469  call pmc_mpi_pack_logical(buffer, position, val%do_select_weighting)
470  call pmc_mpi_pack_integer(buffer, position, val%weighting_type)
471  call pmc_mpi_pack_real(buffer, position, val%weighting_exponent)
472  call pmc_mpi_pack_integer(buffer, position, val%i_repeat)
473  call pmc_mpi_pack_integer(buffer, position, val%n_repeat)
474  call pmc_mpi_pack_real(buffer, position, val%t_wall_start)
475  call pmc_mpi_pack_logical(buffer, position, val%record_removals)
476  call pmc_mpi_pack_logical(buffer, position, val%do_parallel)
477  call pmc_mpi_pack_integer(buffer, position, val%output_type)
478  call pmc_mpi_pack_real(buffer, position, val%mix_timescale)
479  call pmc_mpi_pack_logical(buffer, position, val%gas_average)
480  call pmc_mpi_pack_logical(buffer, position, val%env_average)
481  call pmc_mpi_pack_integer(buffer, position, val%parallel_coag_type)
482  call pmc_mpi_pack_string(buffer, position, val%uuid)
483  call assert(946070052, &
484  position - prev_position <= pmc_mpi_pack_size_run_part_opt(val))
485 #endif
486 
487  end subroutine pmc_mpi_pack_run_part_opt
488 
489 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
490 
491  !> Unpacks the given value from the buffer, advancing position.
492  subroutine pmc_mpi_unpack_run_part_opt(buffer, position, val)
494  !> Memory buffer.
495  character, intent(inout) :: buffer(:)
496  !> Current buffer position.
497  integer, intent(inout) :: position
498  !> Value to pack.
499  type(run_part_opt_t), intent(inout) :: val
500 
501 #ifdef PMC_USE_MPI
502  integer :: prev_position
503 
504  prev_position = position
505  call pmc_mpi_unpack_real(buffer, position, val%t_max)
506  call pmc_mpi_unpack_real(buffer, position, val%t_output)
507  call pmc_mpi_unpack_real(buffer, position, val%t_progress)
508  call pmc_mpi_unpack_real(buffer, position, val%del_t)
509  call pmc_mpi_unpack_string(buffer, position, val%output_prefix)
510  call pmc_mpi_unpack_integer(buffer, position, val%coag_kernel_type)
511  call pmc_mpi_unpack_integer(buffer, position, val%nucleate_type)
512  call pmc_mpi_unpack_integer(buffer, position, val%nucleate_source)
513  call pmc_mpi_unpack_logical(buffer, position, val%do_coagulation)
514  call pmc_mpi_unpack_logical(buffer, position, val%do_nucleation)
515  call pmc_mpi_unpack_logical(buffer, position, val%allow_doubling)
516  call pmc_mpi_unpack_logical(buffer, position, val%allow_halving)
517  call pmc_mpi_unpack_logical(buffer, position, val%do_condensation)
518  call pmc_mpi_unpack_logical(buffer, position, val%do_mosaic)
519  call pmc_mpi_unpack_logical(buffer, position, val%do_optical)
520  call pmc_mpi_unpack_logical(buffer, position, val%do_select_weighting)
521  call pmc_mpi_unpack_integer(buffer, position, val%weighting_type)
522  call pmc_mpi_unpack_real(buffer, position, val%weighting_exponent)
523  call pmc_mpi_unpack_integer(buffer, position, val%i_repeat)
524  call pmc_mpi_unpack_integer(buffer, position, val%n_repeat)
525  call pmc_mpi_unpack_real(buffer, position, val%t_wall_start)
526  call pmc_mpi_unpack_logical(buffer, position, val%record_removals)
527  call pmc_mpi_unpack_logical(buffer, position, val%do_parallel)
528  call pmc_mpi_unpack_integer(buffer, position, val%output_type)
529  call pmc_mpi_unpack_real(buffer, position, val%mix_timescale)
530  call pmc_mpi_unpack_logical(buffer, position, val%gas_average)
531  call pmc_mpi_unpack_logical(buffer, position, val%env_average)
532  call pmc_mpi_unpack_integer(buffer, position, val%parallel_coag_type)
533  call pmc_mpi_unpack_string(buffer, position, val%uuid)
534  call assert(480118362, &
535  position - prev_position <= pmc_mpi_pack_size_run_part_opt(val))
536 #endif
537 
538  end subroutine pmc_mpi_unpack_run_part_opt
539 
540 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
541 
542  !> Read the specification for a parallel coagulation type from a spec file.
543  subroutine spec_file_read_parallel_coag_type(file, parallel_coag_type)
544 
545  !> Spec file.
546  type(spec_file_t), intent(inout) :: file
547  !> Kernel type.
548  integer, intent(out) :: parallel_coag_type
549 
550  character(len=SPEC_LINE_MAX_VAR_LEN) :: parallel_coag_type_name
551 
552  !> \page input_format_parallel_coag Input File Format: Parallel Coagulation Type
553  !!
554  !! The output type is specified by the parameter:
555  !! - \b parallel_coag (string): type of parallel coagulation ---
556  !! must be one of: \c local for only within-process
557  !! coagulation or \c dist to have all processes perform
558  !! coagulation globally, requesting particles from other
559  !! processes as needed
560  !!
561  !! See also:
562  !! - \ref spec_file_format --- the input file text format
563 
564  call spec_file_read_string(file, 'parallel_coag', &
565  parallel_coag_type_name)
566  if (trim(parallel_coag_type_name) == 'local') then
567  parallel_coag_type = parallel_coag_type_local
568  elseif (trim(parallel_coag_type_name) == 'dist') then
569  parallel_coag_type = parallel_coag_type_dist
570  else
571  call spec_file_die_msg(494684716, file, &
572  "Unknown parallel coagulation type: " &
573  // trim(parallel_coag_type_name))
574  end if
575 
576  end subroutine spec_file_read_parallel_coag_type
577 
578 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
579 
580 end module pmc_run_part
integer function aero_state_total_particles_all_procs(aero_state, i_group, i_class)
Returns the total number of particles across all processes.
Definition: aero_state.F90:299
subroutine check_time_multiple(first_name, first_time, second_name, second_time)
Check that the first time interval is close to an integer multiple of the second, and warn if it is n...
Definition: util.F90:376
subroutine aero_state_mix(aero_state, del_t, mix_timescale, aero_data, specify_prob_transfer)
Mix the aero_states between all processes. Currently uses a simple all-to-all diffusion.
An input file with extra data for printing messages.
Definition: spec_file.F90:59
subroutine output_state(prefix, output_type, aero_data, aero_state, gas_data, gas_state, env_state, index, time, del_t, i_repeat, record_removals, record_optical, uuid)
Write the current state.
Definition: output.F90:112
subroutine print_part_progress(i_repeat, t_sim_elapsed, n_part, n_coag, n_emit, n_dil_in, n_dil_out, n_nuc, t_wall_elapsed, t_wall_remain)
Print the current simulation progress to the screen.
Definition: run_part.F90:357
integer function pmc_mpi_pack_size_string(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:385
subroutine mc_coag(coag_kernel_type, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag)
Do coagulation for time del_t.
Definition: coagulation.F90:45
integer function pmc_mpi_rank()
Returns the rank of the current process.
Definition: mpi.F90:117
subroutine aero_state_rebalance(aero_state, aero_data, allow_doubling, allow_halving, initial_state_warning)
Double or halve the particle population in each weight group to maintain close to n_part_ideal partic...
subroutine pmc_mpi_unpack_string(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:868
subroutine pmc_mpi_pack_run_part_opt(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: run_part.F90:442
subroutine mc_coag_dist(coag_kernel_type, env_state, aero_data, aero_state, del_t, tot_n_samp, tot_n_coag)
Do coagulation for time del_t.
The env_state_t structure and associated subroutines.
Definition: env_state.F90:9
character(len=pmc_util_convert_string_len) function time_to_string_max_len(time, max_len)
Convert a time to a string format of maximum length.
Definition: util.F90:912
subroutine scenario_update_gas_state(scenario, delta_t, env_state, old_env_state, gas_data, gas_state)
Do gas emissions and background dilution.
Definition: scenario.F90:209
Water condensation onto aerosol particles.
Definition: condense.F90:29
Generic coagulation kernel.
Definition: coag_kernel.F90:9
The gas_data_t structure and associated subroutines.
Definition: gas_data.F90:9
Interface to the MOSAIC aerosol and gas phase chemistry code.
Definition: mosaic.F90:9
character(len=pmc_util_convert_string_len) function integer_to_string_max_len(val, max_len)
Convert an integer to a string format of maximum length.
Definition: util.F90:836
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:103
subroutine mosaic_cleanup()
Clean-up after running MOSAIC, deallocating memory.
Definition: mosaic.F90:130
subroutine scenario_update_aero_state(scenario, delta_t, env_state, old_env_state, aero_data, aero_state, n_emit, n_dil_in, n_dil_out, allow_doubling, allow_halving)
Do emissions and background dilution for a particle aerosol distribution.
Definition: scenario.F90:259
subroutine scenario_update_env_state(scenario, env_state, time)
Update time-dependent contents of the environment. scenario_init_env_state() should have been called ...
Definition: scenario.F90:134
The aero_state_t structure and assocated subroutines.
Definition: aero_state.F90:9
Aerosol particle coagulation.
Definition: coagulation.F90:11
subroutine condense_particles(aero_state, aero_data, env_state_initial, env_state_final, del_t)
Do condensation to all the particles for a given time interval, including updating the environment to...
Definition: condense.F90:149
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:365
Parallel aerosol particle coagulation with MPI.
Current environment state.
Definition: env_state.F90:26
subroutine mosaic_aero_optical_init(env_state, aero_data, aero_state, gas_data, gas_state)
Compute the optical properties of each aerosol particle for initial timestep.
Definition: mosaic.F90:490
Options controlling the execution of run_part().
Definition: run_part.F90:40
subroutine pmc_mpi_pack_logical(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:638
subroutine check_event(time, timestep, interval, last_time, do_event)
Computes whether an event is scheduled to take place.
Definition: util.F90:407
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:134
integer function aero_state_total_particles(aero_state, i_group, i_class)
Returns the total number of particles in an aerosol distribution.
Definition: aero_state.F90:260
Write data in NetCDF format.
Definition: output.F90:68
subroutine nucleate(nucleate_type, nucleate_source, env_state, gas_data, aero_data, aero_state, gas_state, del_t, allow_doubling, allow_halving)
Do nucleation of the type given by the first argument.
Definition: nucleate.F90:33
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:586
integer, parameter parallel_coag_type_dist
Type code for distributed parallel coagulation.
Definition: run_part.F90:37
The current collection of aerosol particles.
Definition: aero_state.F90:63
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
Definition: spec_file.F90:74
Scenario data.
Definition: scenario.F90:54
The aero_data_t structure and associated subroutines.
Definition: aero_data.F90:9
integer, parameter parallel_coag_type_local
Type code for local parallel coagulation.
Definition: run_part.F90:35
subroutine pmc_mpi_reduce_sum_integer(val, val_sum)
Computes the sum of val across all processes, storing the result in val_sum on the root process...
Definition: mpi.F90:1180
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:561
The gas_state_t structure and associated subroutines.
Definition: gas_state.F90:9
integer function pmc_mpi_pack_size_logical(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:407
subroutine pmc_mpi_pack_string(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:611
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:843
subroutine gas_state_mix(val)
Average val over all processes.
Definition: gas_state.F90:427
subroutine mosaic_timestep(env_state, aero_data, aero_state, gas_data, gas_state, do_optical)
Do one timestep with MOSAIC.
Definition: mosaic.F90:368
integer, parameter dp
Kind of a double precision real number.
Definition: constants.F90:12
subroutine pmc_mpi_unpack_logical(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:896
Current state of the gas mixing ratios in the system.
Definition: gas_state.F90:30
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:345
Constant gas data.
Definition: gas_data.F90:29
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:818
integer function pmc_mpi_pack_size_run_part_opt(val)
Determines the number of bytes required to pack the given value.
Definition: run_part.F90:401
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
Definition: spec_file.F90:604
subroutine env_state_mix(val)
Average val over all processes.
Definition: env_state.F90:269
integer, parameter parallel_coag_type_invalid
Type code for undefined or invalid parallel coagulation method.
Definition: run_part.F90:33
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
Definition: util.F90:766
integer function pmc_mpi_size()
Returns the total number of processes.
Definition: mpi.F90:134
Aerosol material properties and associated data.
Definition: aero_data.F90:41
Common utility subroutines.
Definition: util.F90:9
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine pmc_mpi_unpack_run_part_opt(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: run_part.F90:493
subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, gas_state, run_part_opt)
Do a particle-resolved Monte Carlo simulation.
Definition: run_part.F90:108
Monte Carlo simulation.
Definition: run_part.F90:9
Aerosol nucleation functions.
Definition: nucleate.F90:9
The scenario_t structure and associated subroutines.
Definition: scenario.F90:9
subroutine mosaic_init(env_state, aero_data, del_t, do_optical)
Initialize all MOSAIC data-structures.
Definition: mosaic.F90:38