25 #ifdef PMC_USE_SUNDIALS 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
60 logical :: do_nucleation
62 logical :: allow_doubling
64 logical :: allow_halving
66 logical :: do_condensation
72 logical :: do_select_weighting
74 integer :: weighting_type
76 real(kind=dp) :: weighting_exponent
82 real(kind=dp) :: t_wall_start
84 logical :: record_removals
86 logical :: do_parallel
88 integer :: output_type
90 real(kind=dp) :: mix_timescale
92 logical :: gas_average
94 logical :: env_average
96 integer :: parallel_coag_type
98 character(len=PMC_UUID_LEN) :: uuid
106 subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, &
107 gas_state, run_part_opt)
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
138 integer :: n_time, i_time, i_time_start, pre_i_time
139 integer :: i_state, i_state_netcdf, i_output
152 progress_n_dil_in = 0
153 progress_n_dil_out = 0
157 "del_t", run_part_opt%del_t)
159 "del_t", run_part_opt%del_t)
161 "del_t", run_part_opt%del_t)
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 168 aero_state, gas_data, gas_state)
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)
178 call aero_info_array_zero(aero_state%aero_info_array)
182 run_part_opt%allow_doubling, &
183 run_part_opt%allow_halving, initial_state_warning=.true.)
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
194 if (run_part_opt%i_repeat == 1)
then 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 &
206 global_n_part, 0, 0, 0, 0, 0, t_wall_elapsed, t_wall_remain)
209 do i_time = i_time_start,n_time
211 time =
real(i_time, kind=dp) * run_part_opt%del_t
213 old_env_state = env_state
216 if (run_part_opt%do_nucleation)
then 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)
224 progress_n_nuc = progress_n_nuc + n_nuc
227 if (run_part_opt%do_coagulation)
then 228 if (run_part_opt%parallel_coag_type &
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 &
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)
237 call die_msg(323011762,
"unknown parallel coagulation type: " &
240 progress_n_samp = progress_n_samp + n_samp
241 progress_n_coag = progress_n_coag + n_coag
244 #ifdef PMC_USE_SUNDIALS 245 if (run_part_opt%do_condensation)
then 247 env_state, run_part_opt%del_t)
252 env_state, old_env_state, gas_data, gas_state)
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
261 if (run_part_opt%do_mosaic)
then 263 gas_state, run_part_opt%do_optical)
266 if (run_part_opt%mix_timescale > 0d0)
then 268 run_part_opt%mix_timescale, aero_data)
270 if (run_part_opt%gas_average)
then 273 if (run_part_opt%gas_average)
then 278 run_part_opt%allow_doubling, &
279 run_part_opt%allow_halving, initial_state_warning=.false.)
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)
289 i_output = i_output + 1
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)
299 if (.not. run_part_opt%record_removals)
then 303 call aero_info_array_zero(aero_state%aero_info_array)
306 if (run_part_opt%t_progress > 0d0)
then 308 run_part_opt%t_progress, last_progress_time, do_progress)
309 if (do_progress)
then 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 &
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)
338 progress_n_dil_in = 0
339 progress_n_dil_out = 0
346 if (run_part_opt%do_mosaic)
then 356 n_emit, n_dil_in, n_dil_out, n_nuc, t_wall_elapsed, t_wall_remain)
359 integer,
intent(in) :: i_repeat
361 real(kind=dp),
intent(in) :: t_sim_elapsed
363 integer,
intent(in) :: n_part
365 integer,
intent(in) :: n_coag
367 integer,
intent(in) :: n_emit
369 integer,
intent(in) :: n_dil_in
371 integer,
intent(in) :: n_dil_out
373 integer,
intent(in) :: n_nuc
375 real(kind=dp),
intent(in) :: t_wall_elapsed
377 real(kind=dp),
intent(in) :: t_wall_remain
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)') &
444 character,
intent(inout) :: buffer(:)
446 integer,
intent(inout) :: position
451 integer :: prev_position
453 prev_position = position
495 character,
intent(inout) :: buffer(:)
497 integer,
intent(inout) :: position
502 integer :: prev_position
504 prev_position = position
543 subroutine spec_file_read_parallel_coag_type(file, parallel_coag_type)
548 integer,
intent(out) :: parallel_coag_type
550 character(len=SPEC_LINE_MAX_VAR_LEN) :: parallel_coag_type_name
565 parallel_coag_type_name)
566 if (trim(parallel_coag_type_name) ==
'local')
then 568 elseif (trim(parallel_coag_type_name) ==
'dist')
then 572 "Unknown parallel coagulation type: " &
573 // trim(parallel_coag_type_name))
576 end subroutine spec_file_read_parallel_coag_type
integer function aero_state_total_particles_all_procs(aero_state, i_group, i_class)
Returns the total number of particles across all processes.
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 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 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.
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.
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 pmc_mpi_rank()
Returns the rank of the current process.
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.
subroutine pmc_mpi_pack_run_part_opt(buffer, position, val)
Packs the given value into 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.
The env_state_t structure and associated subroutines.
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.
subroutine scenario_update_gas_state(scenario, delta_t, env_state, old_env_state, gas_data, gas_state)
Do gas emissions and background dilution.
Water condensation onto aerosol particles.
Generic coagulation kernel.
The gas_data_t structure and associated subroutines.
Interface to the MOSAIC aerosol and gas phase chemistry code.
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.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine mosaic_cleanup()
Clean-up after running MOSAIC, deallocating memory.
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.
subroutine scenario_update_env_state(scenario, env_state, time)
Update time-dependent contents of the environment. scenario_init_env_state() should have been called ...
The aero_state_t structure and assocated subroutines.
Aerosol particle coagulation.
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...
integer function pmc_mpi_pack_size_real(val)
Determines the number of bytes required to pack the given value.
Parallel aerosol particle coagulation with MPI.
Current environment state.
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.
Options controlling the execution of run_part().
subroutine pmc_mpi_pack_logical(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine check_event(time, timestep, interval, last_time, do_event)
Computes whether an event is scheduled to take place.
subroutine die_msg(code, error_msg)
Error immediately.
integer function aero_state_total_particles(aero_state, i_group, i_class)
Returns the total number of particles in an aerosol distribution.
Write data in NetCDF format.
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 pmc_mpi_pack_real(buffer, position, val)
Packs the given value into the buffer, advancing position.
integer, parameter parallel_coag_type_dist
Type code for distributed parallel coagulation.
The current collection of aerosol particles.
subroutine spec_file_die_msg(code, file, msg)
Exit with an error message containing filename and line number.
The aero_data_t structure and associated subroutines.
integer, parameter parallel_coag_type_local
Type code for local parallel coagulation.
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...
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
The gas_state_t structure and associated subroutines.
integer function pmc_mpi_pack_size_logical(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_mpi_pack_string(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine pmc_mpi_unpack_real(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine gas_state_mix(val)
Average val over all processes.
subroutine mosaic_timestep(env_state, aero_data, aero_state, gas_data, gas_state, do_optical)
Do one timestep with MOSAIC.
integer, parameter dp
Kind of a double precision real number.
subroutine pmc_mpi_unpack_logical(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Current state of the gas mixing ratios in the system.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_run_part_opt(val)
Determines the number of bytes required to pack the given value.
subroutine spec_file_read_string(file, name, var)
Read a string from a spec file that must have a given name.
subroutine env_state_mix(val)
Average val over all processes.
integer, parameter parallel_coag_type_invalid
Type code for undefined or invalid parallel coagulation method.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
integer function pmc_mpi_size()
Returns the total number of processes.
Aerosol material properties and associated data.
Common utility subroutines.
Wrapper functions for MPI.
subroutine pmc_mpi_unpack_run_part_opt(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine run_part(scenario, env_state, aero_data, aero_state, gas_data, gas_state, run_part_opt)
Do a particle-resolved Monte Carlo simulation.
Aerosol nucleation functions.
The scenario_t structure and associated subroutines.
subroutine mosaic_init(env_state, aero_data, del_t, do_optical)
Initialize all MOSAIC data-structures.