25 #ifdef PMC_USE_SUNDIALS
33 integer,
parameter :: PARALLEL_COAG_TYPE_INVALID = 0
35 integer,
parameter :: PARALLEL_COAG_TYPE_LOCAL = 1
37 integer,
parameter :: PARALLEL_COAG_TYPE_DIST = 2
42 real(kind=dp) :: t_max
44 real(kind=dp) :: t_output
46 real(kind=dp) :: t_progress
48 real(kind=dp) :: del_t
50 character(len=300) :: output_prefix
52 integer :: coag_kernel_type
54 integer :: nucleate_type
56 integer :: nucleate_source
58 logical :: do_coagulation
62 integer :: loss_function_type
66 real(kind=dp) :: loss_alg_threshold
68 logical :: do_nucleation
70 logical :: allow_doubling
72 logical :: allow_halving
74 logical :: do_condensation
84 real(kind=dp) :: t_wall_start
86 logical :: record_removals
88 logical :: do_parallel
90 integer :: output_type
92 real(kind=dp) :: mix_timescale
94 logical :: gas_average
96 logical :: env_average
98 integer :: parallel_coag_type
100 character(len=PMC_UUID_LEN) :: uuid
108 subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, &
109 gas_state, run_part_opt)
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
140 integer :: n_time, i_time, i_time_start, pre_i_time
141 integer :: i_state, i_state_netcdf, i_output
154 progress_n_dil_in = 0
155 progress_n_dil_out = 0
159 "del_t", run_part_opt%del_t)
161 "del_t", run_part_opt%del_t)
163 "del_t", run_part_opt%del_t)
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)
172 if (run_part_opt%t_output > 0d0)
then
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)
182 run_part_opt%allow_halving, initial_state_warning=.true.)
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
193 if (run_part_opt%i_repeat == 1)
then
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 &
205 global_n_part, 0, 0, 0, 0, 0, t_wall_elapsed, t_wall_remain)
208 do i_time = i_time_start,n_time
210 time =
real(i_time, kind=dp) * run_part_opt%del_t
215 if (run_part_opt%do_nucleation)
then
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)
223 progress_n_nuc = progress_n_nuc + n_nuc
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)
236 call
die_msg(323011762,
"unknown parallel coagulation type: " &
239 progress_n_samp = progress_n_samp + n_samp
240 progress_n_coag = progress_n_coag + n_coag
243 if (run_part_opt%do_loss)
then
245 run_part_opt%del_t, aero_data, aero_state, env_state)
248 #ifdef PMC_USE_SUNDIALS
249 if (run_part_opt%do_condensation)
then
251 env_state, run_part_opt%del_t)
256 env_state, old_env_state, gas_data, gas_state)
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
265 if (run_part_opt%do_mosaic)
then
267 gas_state, run_part_opt%do_optical)
270 if (run_part_opt%mix_timescale > 0d0)
then
272 run_part_opt%mix_timescale, aero_data)
274 if (run_part_opt%gas_average)
then
277 if (run_part_opt%gas_average)
then
282 run_part_opt%allow_halving, initial_state_warning=.false.)
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)
292 i_output = i_output + 1
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)
302 if (.not. run_part_opt%record_removals)
then
309 if (run_part_opt%t_progress > 0d0)
then
311 run_part_opt%t_progress, last_progress_time, do_progress)
312 if (do_progress)
then
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 &
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)
341 progress_n_dil_in = 0
342 progress_n_dil_out = 0
349 if (run_part_opt%do_mosaic)
then
361 n_emit, n_dil_in, n_dil_out, n_nuc, t_wall_elapsed, t_wall_remain)
364 integer,
intent(in) :: i_repeat
366 real(kind=dp),
intent(in) :: t_sim_elapsed
368 integer,
intent(in) :: n_part
370 integer,
intent(in) :: n_coag
372 integer,
intent(in) :: n_emit
374 integer,
intent(in) :: n_dil_in
376 integer,
intent(in) :: n_dil_out
378 integer,
intent(in) :: n_nuc
380 real(kind=dp),
intent(in) :: t_wall_elapsed
382 real(kind=dp),
intent(in) :: t_wall_remain
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)') &
449 character,
intent(inout) :: buffer(:)
451 integer,
intent(inout) :: position
456 integer :: prev_position
458 prev_position = position
500 character,
intent(inout) :: buffer(:)
502 integer,
intent(inout) :: position
507 integer :: prev_position
509 prev_position = position
548 subroutine spec_file_read_parallel_coag_type(file, parallel_coag_type)
553 integer,
intent(out) :: parallel_coag_type
555 character(len=SPEC_LINE_MAX_VAR_LEN) :: parallel_coag_type_name
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
577 "Unknown parallel coagulation type: " &
578 // trim(parallel_coag_type_name))
581 end subroutine spec_file_read_parallel_coag_type
The scenario_t structure and associated subroutines.
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.
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
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.
subroutine pmc_mpi_pack_logical(buffer, position, val)
Packs the given value into the buffer, advancing position.
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.
The aero_data_t structure and associated subroutines.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, gas_state, run_part_opt)
Do a particle-resolved Monte Carlo simulation.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine pmc_mpi_unpack_run_part_opt(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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() ...
Water condensation onto aerosol particles.
The env_state_t structure and associated subroutines.
Aerosol particle coagulation.
subroutine pmc_mpi_pack_string(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine scenario_update_env_state(scenario, env_state, time)
Update time-dependent contents of the environment. scenario_init_env_state() should have been called ...
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.
integer function pmc_mpi_pack_size_string(val)
Determines the number of bytes required to pack the given value.
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.
integer function pmc_mpi_pack_size_run_part_opt(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
The gas_data_t structure and associated subroutines.
integer function pmc_mpi_size()
Returns the total number of processes.
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...
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...
Interface to the MOSAIC aerosol and gas phase chemistry code.
Common utility subroutines.
subroutine env_state_mix(val)
Average val over all processes.
Current environment state.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Options controlling the execution of run_part().
subroutine env_state_deallocate(env_state)
Free all storage.
The aero_state_t structure and assocated subroutines.
integer function pmc_mpi_rank()
Returns the rank of the current process.
Wrapper functions for MPI.
subroutine mosaic_timestep(env_state, aero_data, aero_state, gas_data, gas_state, do_optical)
Do one timestep with MOSAIC.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
subroutine pmc_mpi_unpack_logical(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
The gas_state_t structure and associated subroutines.
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.
integer function aero_state_total_particles(aero_state, i_group, i_class)
Returns the total number of particles in an aerosol distribution.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
The current collection of aerosol particles.
Aerosol nucleation functions.
subroutine scenario_update_gas_state(scenario, delta_t, env_state, old_env_state, gas_data, gas_state)
Do gas emissions and background dilution.
subroutine env_state_copy(env_from, env_to)
env_to = env_from
subroutine pmc_mpi_unpack_string(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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...
subroutine pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine mosaic_init(env_state, aero_data, del_t, do_optical)
Initialize all MOSAIC data-structures.
subroutine check_event(time, timestep, interval, last_time, do_event)
Computes whether an event is scheduled to take place.
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
subroutine gas_state_mix(val)
Average val over all processes.
Current state of the gas mixing ratios in the system.
subroutine mosaic_cleanup()
Clean-up after running MOSAIC, deallocating memory.
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.
subroutine env_state_allocate(env_state)
Allocate an empty environment.
subroutine pmc_mpi_pack_run_part_opt(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer function aero_state_total_particles_all_procs(aero_state, i_group, i_class)
Returns the total number of particles across all processes.
Aerosol material properties and associated data.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
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.
Write data in NetCDF format.
Generic coagulation kernel.
integer function pmc_mpi_pack_size_logical(val)
Determines the number of bytes required to pack the given value.
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...