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