PartMC
2.1.5
|
00001 ! Copyright (C) 2005-2011 Nicole Riemer and Matthew West 00002 ! Licensed under the GNU General Public License version 2 or (at your 00003 ! option) any later version. See the file COPYING for details. 00004 00005 !> \file 00006 !> The pmc_run_part module. 00007 00008 !> Monte Carlo simulation. 00009 module pmc_run_part 00010 00011 use pmc_util 00012 use pmc_aero_state 00013 use pmc_bin_grid 00014 use pmc_env_data 00015 use pmc_env_state 00016 use pmc_aero_data 00017 use pmc_gas_data 00018 use pmc_gas_state 00019 use pmc_output 00020 use pmc_mosaic 00021 use pmc_coagulation 00022 use pmc_coagulation_dist 00023 use pmc_coag_kernel 00024 use pmc_nucleate 00025 use pmc_mpi 00026 #ifdef PMC_USE_SUNDIALS 00027 use pmc_condense 00028 #endif 00029 #ifdef PMC_USE_MPI 00030 use mpi 00031 #endif 00032 00033 !> Type code for undefined or invalid parallel coagulation method. 00034 integer, parameter :: PARALLEL_COAG_TYPE_INVALID = 0 00035 !> Type code for local parallel coagulation. 00036 integer, parameter :: PARALLEL_COAG_TYPE_LOCAL = 1 00037 !> Type code for distributed parallel coagulation. 00038 integer, parameter :: PARALLEL_COAG_TYPE_DIST = 2 00039 00040 !> Options controlling the execution of run_part(). 00041 type run_part_opt_t 00042 !> Preferred number of particles. 00043 integer :: n_part_ideal 00044 !> Final time (s). 00045 real(kind=dp) :: t_max 00046 !> Output interval (0 disables) (s). 00047 real(kind=dp) :: t_output 00048 !> Progress interval (0 disables) (s). 00049 real(kind=dp) :: t_progress 00050 !> Timestep for coagulation. 00051 real(kind=dp) :: del_t 00052 !> Prefix for output files. 00053 character(len=300) :: output_prefix 00054 !> Type of coagulation kernel. 00055 integer :: coag_kernel_type 00056 !> Type of nucleation. 00057 integer :: nucleate_type 00058 !> Whether to do coagulation. 00059 logical :: do_coagulation 00060 !> Whether to do nucleation. 00061 logical :: do_nucleation 00062 !> Allow doubling if needed. 00063 logical :: allow_doubling 00064 !> Allow halving if needed. 00065 logical :: allow_halving 00066 !> Whether to do condensation. 00067 logical :: do_condensation 00068 !> Whether to do MOSAIC. 00069 logical :: do_mosaic 00070 !> Whether to compute optical properties. 00071 logical :: do_optical 00072 !> Repeat number of run. 00073 integer :: i_repeat 00074 !> Total number of repeats. 00075 integer :: n_repeat 00076 !> Cpu_time() of start. 00077 real(kind=dp) :: t_wall_start 00078 !> Whether to record particle removal information. 00079 logical :: record_removals 00080 !> Whether to run in parallel. 00081 logical :: do_parallel 00082 !> Parallel output type. 00083 integer :: output_type 00084 !> Mixing timescale between processes (s). 00085 real(kind=dp) :: mix_timescale 00086 !> Whether to average gases each timestep. 00087 logical :: gas_average 00088 !> Whether to average environment each timestep. 00089 logical :: env_average 00090 !> Parallel coagulation method type. 00091 integer :: parallel_coag_type 00092 !> UUID for this simulation. 00093 character(len=PMC_UUID_LEN) :: uuid 00094 end type run_part_opt_t 00095 00096 contains 00097 00098 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00099 00100 !> Do a particle-resolved Monte Carlo simulation. 00101 subroutine run_part(bin_grid, env_data, env_state, aero_data, & 00102 aero_weight, aero_state, gas_data, gas_state, run_part_opt) 00103 00104 !> Bin grid. 00105 type(bin_grid_t), intent(in) :: bin_grid 00106 !> Environment state. 00107 type(env_data_t), intent(in) :: env_data 00108 !> Environment state. 00109 type(env_state_t), intent(inout) :: env_state 00110 !> Aerosol data. 00111 type(aero_data_t), intent(in) :: aero_data 00112 !> Aerosol weight. 00113 type(aero_weight_t), intent(in) :: aero_weight 00114 !> Aerosol state. 00115 type(aero_state_t), intent(inout) :: aero_state 00116 !> Gas data. 00117 type(gas_data_t), intent(in) :: gas_data 00118 !> Gas state. 00119 type(gas_state_t), intent(inout) :: gas_state 00120 !> Monte Carlo options. 00121 type(run_part_opt_t), intent(in) :: run_part_opt 00122 00123 real(kind=dp) :: time, pre_time, pre_del_t, prop_done 00124 real(kind=dp) :: last_output_time, last_progress_time 00125 real(kind=dp) :: k_max(bin_grid%n_bin, bin_grid%n_bin) 00126 integer :: rank, n_proc, pre_index, ncid 00127 integer :: pre_i_repeat 00128 integer :: n_samp, n_coag, n_emit, n_dil_in, n_dil_out, n_nuc 00129 integer :: progress_n_samp, progress_n_coag 00130 integer :: progress_n_emit, progress_n_dil_in, progress_n_dil_out 00131 integer :: progress_n_nuc, n_part_before 00132 integer :: global_n_part, global_n_samp, global_n_coag 00133 integer :: global_n_emit, global_n_dil_in, global_n_dil_out 00134 integer :: global_n_nuc 00135 logical :: do_output, do_state, do_state_netcdf, do_progress, did_coag 00136 logical :: update_rel_humid 00137 real(kind=dp) :: t_start, t_wall_now, t_wall_elapsed, t_wall_remain 00138 type(env_state_t) :: old_env_state 00139 integer :: n_time, i_time, i_time_start, pre_i_time 00140 integer :: i_state, i_state_netcdf, i_output 00141 00142 rank = pmc_mpi_rank() 00143 n_proc = pmc_mpi_size() 00144 00145 i_time = 0 00146 i_output = 1 00147 i_state = 1 00148 i_state_netcdf = 1 00149 time = 0d0 00150 progress_n_samp = 0 00151 progress_n_coag = 0 00152 progress_n_emit = 0 00153 progress_n_dil_in = 0 00154 progress_n_dil_out = 0 00155 progress_n_nuc = 0 00156 00157 call env_state_allocate(old_env_state) 00158 00159 if (run_part_opt%do_coagulation) then 00160 call est_k_max_binned(bin_grid, run_part_opt%coag_kernel_type, & 00161 aero_data, aero_weight, env_state, k_max) 00162 end if 00163 00164 if (run_part_opt%do_mosaic) then 00165 call mosaic_init(bin_grid, env_state, run_part_opt%del_t, & 00166 run_part_opt%do_optical) 00167 end if 00168 00169 if (run_part_opt%t_output > 0d0) then 00170 call output_state(run_part_opt%output_prefix, & 00171 run_part_opt%output_type, bin_grid, aero_data, aero_weight, & 00172 aero_state, gas_data, gas_state, env_state, i_state, time, & 00173 run_part_opt%del_t, run_part_opt%i_repeat, & 00174 run_part_opt%record_removals, run_part_opt%do_optical, & 00175 run_part_opt%uuid) 00176 call aero_info_array_zero(aero_state%aero_info_array) 00177 end if 00178 00179 ! Do an initial double/halve test. This shouldn't happen (except 00180 ! for restart with inconsistent number) so issue a warning. 00181 if (run_part_opt%allow_doubling) then 00182 global_n_part = aero_state_total_particles_all_procs(aero_state) 00183 do while ((global_n_part & 00184 < run_part_opt%n_part_ideal * pmc_mpi_size() / 2) & 00185 .and. (global_n_part > 0)) 00186 call warn_msg(716882783, "doubling particles in initial condition") 00187 call aero_state_double(aero_state) 00188 global_n_part = aero_state_total_particles_all_procs(aero_state) 00189 end do 00190 end if 00191 if (run_part_opt%allow_halving) then 00192 do while (aero_state_total_particles_all_procs(aero_state) & 00193 > run_part_opt%n_part_ideal * pmc_mpi_size() * 2) 00194 call warn_msg(661936373, "halving particles in initial condition") 00195 call aero_state_halve(aero_state, bin_grid) 00196 end do 00197 end if 00198 00199 t_start = env_state%elapsed_time 00200 last_output_time = time 00201 last_progress_time = time 00202 n_time = nint(run_part_opt%t_max / run_part_opt%del_t) 00203 i_time_start = nint(time / run_part_opt%del_t) + 1 00204 00205 global_n_part = aero_state_total_particles_all_procs(aero_state) 00206 if (rank == 0) then 00207 ! progress only printed from root process 00208 if (run_part_opt%i_repeat == 1) then 00209 t_wall_elapsed = 0d0 00210 t_wall_remain = 0d0 00211 else 00212 call cpu_time(t_wall_now) 00213 prop_done = real(run_part_opt%i_repeat - 1, kind=dp) 00214 / real(run_part_opt%n_repeat, kind=dp) 00215 t_wall_elapsed = t_wall_now - run_part_opt%t_wall_start 00216 t_wall_remain = (1d0 - prop_done) / prop_done & 00217 * t_wall_elapsed 00218 end if 00219 call print_part_progress(run_part_opt%i_repeat, time, & 00220 global_n_part, 0, 0, 0, 0, 0, t_wall_elapsed, t_wall_remain) 00221 end if 00222 00223 do i_time = i_time_start,n_time 00224 00225 time = real(i_time, kind=dp) * run_part_opt%del_t 00226 00227 call env_state_copy(env_state, old_env_state) 00228 update_rel_humid = .not. run_part_opt%do_condensation 00229 call env_data_update_state(env_data, env_state, time + t_start, & 00230 update_rel_humid) 00231 00232 if (run_part_opt%do_nucleation) then 00233 n_part_before = aero_state_total_particles(aero_state) 00234 call nucleate(run_part_opt%nucleate_type, bin_grid, env_state, & 00235 gas_data, aero_data, aero_weight, aero_state, gas_state, & 00236 run_part_opt%del_t) 00237 n_nuc = aero_state_total_particles(aero_state) & 00238 - n_part_before 00239 progress_n_nuc = progress_n_nuc + n_nuc 00240 end if 00241 00242 if (run_part_opt%do_coagulation) then 00243 if (run_part_opt%parallel_coag_type & 00244 == PARALLEL_COAG_TYPE_LOCAL) then 00245 call mc_coag(run_part_opt%coag_kernel_type, bin_grid, & 00246 env_state, aero_data, aero_weight, aero_state, & 00247 run_part_opt%del_t, k_max, n_samp, n_coag) 00248 elseif (run_part_opt%parallel_coag_type & 00249 == PARALLEL_COAG_TYPE_DIST) then 00250 call mc_coag_dist(run_part_opt%coag_kernel_type, bin_grid, & 00251 env_state, aero_data, aero_weight, aero_state, & 00252 run_part_opt%del_t, k_max, n_samp, n_coag) 00253 else 00254 call die_msg(323011762, "unknown parallel coagulation type: " & 00255 // trim(integer_to_string(run_part_opt%parallel_coag_type))) 00256 end if 00257 progress_n_samp = progress_n_samp + n_samp 00258 progress_n_coag = progress_n_coag + n_coag 00259 end if 00260 00261 call env_state_update_gas_state(env_state, run_part_opt%del_t, & 00262 old_env_state, gas_data, gas_state) 00263 call env_state_update_aero_state(env_state, run_part_opt%del_t, & 00264 old_env_state, bin_grid, aero_data, aero_weight, aero_state, & 00265 n_emit, n_dil_in, n_dil_out) 00266 progress_n_emit = progress_n_emit + n_emit 00267 progress_n_dil_in = progress_n_dil_in + n_dil_in 00268 progress_n_dil_out = progress_n_dil_out + n_dil_out 00269 00270 #ifdef PMC_USE_SUNDIALS 00271 if (run_part_opt%do_condensation) then 00272 call condense_particles(bin_grid, env_state, env_data, aero_data, & 00273 aero_weight, aero_state, run_part_opt%del_t) 00274 end if 00275 #endif 00276 00277 if (run_part_opt%do_mosaic) then 00278 call mosaic_timestep(bin_grid, env_state, aero_data, & 00279 aero_weight, aero_state, gas_data, gas_state, & 00280 run_part_opt%do_optical) 00281 end if 00282 00283 if (run_part_opt%mix_timescale > 0d0) then 00284 call aero_state_mix(aero_state, run_part_opt%del_t, & 00285 run_part_opt%mix_timescale, aero_data, bin_grid) 00286 end if 00287 if (run_part_opt%gas_average) then 00288 call gas_state_mix(gas_state) 00289 end if 00290 if (run_part_opt%gas_average) then 00291 call env_state_mix(env_state) 00292 end if 00293 00294 ! if we have less than half the maximum number of particles then 00295 ! double until we fill up the array 00296 if (run_part_opt%allow_doubling) then 00297 global_n_part = aero_state_total_particles_all_procs(aero_state) 00298 do while ((global_n_part & 00299 < run_part_opt%n_part_ideal * pmc_mpi_size() / 2) & 00300 .and. (global_n_part > 0)) 00301 call aero_state_double(aero_state) 00302 global_n_part = aero_state_total_particles_all_procs(aero_state) 00303 end do 00304 end if 00305 ! same for halving if we have too many particles 00306 if (run_part_opt%allow_halving) then 00307 do while (aero_state_total_particles_all_procs(aero_state) & 00308 > run_part_opt%n_part_ideal * pmc_mpi_size() * 2) 00309 call aero_state_halve(aero_state, bin_grid) 00310 end do 00311 end if 00312 00313 ! DEBUG: enable to check array handling 00314 ! call aero_state_check(bin_grid, aero_data, aero_state) 00315 ! DEBUG: end 00316 00317 if (run_part_opt%t_output > 0d0) then 00318 call check_event(time, run_part_opt%del_t, run_part_opt%t_output, & 00319 last_output_time, do_output) 00320 if (do_output) then 00321 i_output = i_output + 1 00322 call output_state(run_part_opt%output_prefix, & 00323 run_part_opt%output_type, bin_grid, aero_data, & 00324 aero_weight, aero_state, gas_data, gas_state, env_state, & 00325 i_output, time, run_part_opt%del_t, & 00326 run_part_opt%i_repeat, run_part_opt%record_removals, & 00327 run_part_opt%do_optical, run_part_opt%uuid) 00328 call aero_info_array_zero(aero_state%aero_info_array) 00329 end if 00330 end if 00331 00332 if (.not. run_part_opt%record_removals) then 00333 ! If we are not recording removals then we can zero them as 00334 ! often as possible to minimize the cost of maintaining 00335 ! them. 00336 call aero_info_array_zero(aero_state%aero_info_array) 00337 end if 00338 00339 if (run_part_opt%t_progress > 0d0) then 00340 call check_event(time, run_part_opt%del_t, & 00341 run_part_opt%t_progress, last_progress_time, do_progress) 00342 if (do_progress) then 00343 global_n_part = aero_state_total_particles_all_procs(aero_state) 00344 call pmc_mpi_reduce_sum_integer(progress_n_samp, global_n_samp) 00345 call pmc_mpi_reduce_sum_integer(progress_n_coag, global_n_coag) 00346 call pmc_mpi_reduce_sum_integer(progress_n_emit, global_n_emit) 00347 call pmc_mpi_reduce_sum_integer(progress_n_dil_in, & 00348 global_n_dil_in) 00349 call pmc_mpi_reduce_sum_integer(progress_n_dil_out, & 00350 global_n_dil_out) 00351 call pmc_mpi_reduce_sum_integer(progress_n_nuc, global_n_nuc) 00352 if (rank == 0) then 00353 ! progress only printed from root process 00354 call cpu_time(t_wall_now) 00355 prop_done = (real(run_part_opt%i_repeat - 1, kind=dp) 00356 + time / run_part_opt%t_max) 00357 / real(run_part_opt%n_repeat, kind=dp) 00358 t_wall_elapsed = t_wall_now - run_part_opt%t_wall_start 00359 t_wall_remain = (1d0 - prop_done) / prop_done & 00360 * t_wall_elapsed 00361 call print_part_progress(run_part_opt%i_repeat, time, & 00362 global_n_part, global_n_coag, global_n_emit, & 00363 global_n_dil_in, global_n_dil_out, global_n_nuc, & 00364 t_wall_elapsed, t_wall_remain) 00365 end if 00366 ! reset counters so they show information since last 00367 ! progress display 00368 progress_n_samp = 0 00369 progress_n_coag = 0 00370 progress_n_emit = 0 00371 progress_n_dil_in = 0 00372 progress_n_dil_out = 0 00373 progress_n_nuc = 0 00374 end if 00375 end if 00376 00377 end do 00378 00379 if (run_part_opt%do_mosaic) then 00380 call mosaic_cleanup() 00381 end if 00382 00383 call env_state_deallocate(old_env_state) 00384 00385 end subroutine run_part 00386 00387 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00388 00389 !> Print the current simulation progress to the screen. 00390 subroutine print_part_progress(i_repeat, t_sim_elapsed, n_part, n_coag, & 00391 n_emit, n_dil_in, n_dil_out, n_nuc, t_wall_elapsed, t_wall_remain) 00392 00393 !> Repeat number of simulation. 00394 integer, intent(in) :: i_repeat 00395 !> Elapsed simulation time (s). 00396 real(kind=dp), intent(in) :: t_sim_elapsed 00397 !> Number of particles. 00398 integer, intent(in) :: n_part 00399 !> Number of coagulated particles since last progress printing. 00400 integer, intent(in) :: n_coag 00401 !> Number of emitted particles since last progress printing. 00402 integer, intent(in) :: n_emit 00403 !> Number of diluted-in particles since last progress printing. 00404 integer, intent(in) :: n_dil_in 00405 !> Number of diluted-out particles since last progress printing. 00406 integer, intent(in) :: n_dil_out 00407 !> Number of nucleated particles since last progress printing. 00408 integer, intent(in) :: n_nuc 00409 !> Elapsed wall time (s). 00410 real(kind=dp), intent(in) :: t_wall_elapsed 00411 !> Estimated remaining wall time (s). 00412 real(kind=dp), intent(in) :: t_wall_remain 00413 00414 write(*,'(a6,a1,a6,a1,a7,a1,a7,a1,a7,a1,a8,a1,a9,a1,a7,a1,a6,a1,a6)') & 00415 "repeat", " ", "t_sim", " ", "n_part", " ", "n_coag", " ", & 00416 "n_emit", " ", "n_dil_in", " ", "n_dil_out", " ", "n_nuc", " ", & 00417 "t_wall", " ", "t_rem" 00418 write(*,'(a6,a1,a6,a1,a7,a1,a7,a1,a7,a1,a8,a1,a9,a1,a7,a1,a6,a1,a6)') & 00419 trim(integer_to_string_max_len(i_repeat, 6)), " ", & 00420 trim(time_to_string_max_len(t_sim_elapsed, 6)), " ", & 00421 trim(integer_to_string_max_len(n_part, 7)), " ", & 00422 trim(integer_to_string_max_len(n_coag, 7)), " ", & 00423 trim(integer_to_string_max_len(n_emit, 7)), " ", & 00424 trim(integer_to_string_max_len(n_dil_in, 7)), " ", & 00425 trim(integer_to_string_max_len(n_dil_out, 7)), " ", & 00426 trim(integer_to_string_max_len(n_nuc, 7)), " ", & 00427 trim(time_to_string_max_len(t_wall_elapsed, 6)), " ", & 00428 trim(time_to_string_max_len(t_wall_remain, 6)) 00429 00430 end subroutine print_part_progress 00431 00432 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00433 00434 !> Determines the number of bytes required to pack the given value. 00435 integer function pmc_mpi_pack_size_run_part_opt(val) 00436 00437 !> Value to pack. 00438 type(run_part_opt_t), intent(in) :: val 00439 00440 pmc_mpi_pack_size_run_part_opt = & 00441 pmc_mpi_pack_size_integer(val%n_part_ideal) & 00442 + pmc_mpi_pack_size_real(val%t_max) & 00443 + pmc_mpi_pack_size_real(val%t_output) & 00444 + pmc_mpi_pack_size_real(val%t_progress) & 00445 + pmc_mpi_pack_size_real(val%del_t) & 00446 + pmc_mpi_pack_size_string(val%output_prefix) & 00447 + pmc_mpi_pack_size_integer(val%coag_kernel_type) & 00448 + pmc_mpi_pack_size_integer(val%nucleate_type) & 00449 + pmc_mpi_pack_size_logical(val%do_coagulation) & 00450 + pmc_mpi_pack_size_logical(val%do_nucleation) & 00451 + pmc_mpi_pack_size_logical(val%allow_doubling) & 00452 + pmc_mpi_pack_size_logical(val%allow_halving) & 00453 + pmc_mpi_pack_size_logical(val%do_condensation) & 00454 + pmc_mpi_pack_size_logical(val%do_mosaic) & 00455 + pmc_mpi_pack_size_logical(val%do_optical) & 00456 + pmc_mpi_pack_size_integer(val%i_repeat) & 00457 + pmc_mpi_pack_size_integer(val%n_repeat) & 00458 + pmc_mpi_pack_size_real(val%t_wall_start) & 00459 + pmc_mpi_pack_size_logical(val%record_removals) & 00460 + pmc_mpi_pack_size_logical(val%do_parallel) & 00461 + pmc_mpi_pack_size_integer(val%output_type) & 00462 + pmc_mpi_pack_size_real(val%mix_timescale) & 00463 + pmc_mpi_pack_size_logical(val%gas_average) & 00464 + pmc_mpi_pack_size_logical(val%env_average) & 00465 + pmc_mpi_pack_size_integer(val%parallel_coag_type) & 00466 + pmc_mpi_pack_size_string(val%uuid) 00467 00468 end function pmc_mpi_pack_size_run_part_opt 00469 00470 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00471 00472 !> Packs the given value into the buffer, advancing position. 00473 subroutine pmc_mpi_pack_run_part_opt(buffer, position, val) 00474 00475 !> Memory buffer. 00476 character, intent(inout) :: buffer(:) 00477 !> Current buffer position. 00478 integer, intent(inout) :: position 00479 !> Value to pack. 00480 type(run_part_opt_t), intent(in) :: val 00481 00482 #ifdef PMC_USE_MPI 00483 integer :: prev_position 00484 00485 prev_position = position 00486 call pmc_mpi_pack_integer(buffer, position, val%n_part_ideal) 00487 call pmc_mpi_pack_real(buffer, position, val%t_max) 00488 call pmc_mpi_pack_real(buffer, position, val%t_output) 00489 call pmc_mpi_pack_real(buffer, position, val%t_progress) 00490 call pmc_mpi_pack_real(buffer, position, val%del_t) 00491 call pmc_mpi_pack_string(buffer, position, val%output_prefix) 00492 call pmc_mpi_pack_integer(buffer, position, val%coag_kernel_type) 00493 call pmc_mpi_pack_integer(buffer, position, val%nucleate_type) 00494 call pmc_mpi_pack_logical(buffer, position, val%do_coagulation) 00495 call pmc_mpi_pack_logical(buffer, position, val%do_nucleation) 00496 call pmc_mpi_pack_logical(buffer, position, val%allow_doubling) 00497 call pmc_mpi_pack_logical(buffer, position, val%allow_halving) 00498 call pmc_mpi_pack_logical(buffer, position, val%do_condensation) 00499 call pmc_mpi_pack_logical(buffer, position, val%do_mosaic) 00500 call pmc_mpi_pack_logical(buffer, position, val%do_optical) 00501 call pmc_mpi_pack_integer(buffer, position, val%i_repeat) 00502 call pmc_mpi_pack_integer(buffer, position, val%n_repeat) 00503 call pmc_mpi_pack_real(buffer, position, val%t_wall_start) 00504 call pmc_mpi_pack_logical(buffer, position, val%record_removals) 00505 call pmc_mpi_pack_logical(buffer, position, val%do_parallel) 00506 call pmc_mpi_pack_integer(buffer, position, val%output_type) 00507 call pmc_mpi_pack_real(buffer, position, val%mix_timescale) 00508 call pmc_mpi_pack_logical(buffer, position, val%gas_average) 00509 call pmc_mpi_pack_logical(buffer, position, val%env_average) 00510 call pmc_mpi_pack_integer(buffer, position, val%parallel_coag_type) 00511 call pmc_mpi_pack_string(buffer, position, val%uuid) 00512 call assert(946070052, & 00513 position - prev_position <= pmc_mpi_pack_size_run_part_opt(val)) 00514 #endif 00515 00516 end subroutine pmc_mpi_pack_run_part_opt 00517 00518 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00519 00520 !> Unpacks the given value from the buffer, advancing position. 00521 subroutine pmc_mpi_unpack_run_part_opt(buffer, position, val) 00522 00523 !> Memory buffer. 00524 character, intent(inout) :: buffer(:) 00525 !> Current buffer position. 00526 integer, intent(inout) :: position 00527 !> Value to pack. 00528 type(run_part_opt_t), intent(inout) :: val 00529 00530 #ifdef PMC_USE_MPI 00531 integer :: prev_position 00532 00533 prev_position = position 00534 call pmc_mpi_unpack_integer(buffer, position, val%n_part_ideal) 00535 call pmc_mpi_unpack_real(buffer, position, val%t_max) 00536 call pmc_mpi_unpack_real(buffer, position, val%t_output) 00537 call pmc_mpi_unpack_real(buffer, position, val%t_progress) 00538 call pmc_mpi_unpack_real(buffer, position, val%del_t) 00539 call pmc_mpi_unpack_string(buffer, position, val%output_prefix) 00540 call pmc_mpi_unpack_integer(buffer, position, val%coag_kernel_type) 00541 call pmc_mpi_unpack_integer(buffer, position, val%nucleate_type) 00542 call pmc_mpi_unpack_logical(buffer, position, val%do_coagulation) 00543 call pmc_mpi_unpack_logical(buffer, position, val%do_nucleation) 00544 call pmc_mpi_unpack_logical(buffer, position, val%allow_doubling) 00545 call pmc_mpi_unpack_logical(buffer, position, val%allow_halving) 00546 call pmc_mpi_unpack_logical(buffer, position, val%do_condensation) 00547 call pmc_mpi_unpack_logical(buffer, position, val%do_mosaic) 00548 call pmc_mpi_unpack_logical(buffer, position, val%do_optical) 00549 call pmc_mpi_unpack_integer(buffer, position, val%i_repeat) 00550 call pmc_mpi_unpack_integer(buffer, position, val%n_repeat) 00551 call pmc_mpi_unpack_real(buffer, position, val%t_wall_start) 00552 call pmc_mpi_unpack_logical(buffer, position, val%record_removals) 00553 call pmc_mpi_unpack_logical(buffer, position, val%do_parallel) 00554 call pmc_mpi_unpack_integer(buffer, position, val%output_type) 00555 call pmc_mpi_unpack_real(buffer, position, val%mix_timescale) 00556 call pmc_mpi_unpack_logical(buffer, position, val%gas_average) 00557 call pmc_mpi_unpack_logical(buffer, position, val%env_average) 00558 call pmc_mpi_unpack_integer(buffer, position, val%parallel_coag_type) 00559 call pmc_mpi_unpack_string(buffer, position, val%uuid) 00560 call assert(480118362, & 00561 position - prev_position <= pmc_mpi_pack_size_run_part_opt(val)) 00562 #endif 00563 00564 end subroutine pmc_mpi_unpack_run_part_opt 00565 00566 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00567 00568 !> Read the specification for a parallel coagulation type from a spec file. 00569 subroutine spec_file_read_parallel_coag_type(file, parallel_coag_type) 00570 00571 !> Spec file. 00572 type(spec_file_t), intent(inout) :: file 00573 !> Kernel type. 00574 integer, intent(out) :: parallel_coag_type 00575 00576 character(len=SPEC_LINE_MAX_VAR_LEN) :: parallel_coag_type_name 00577 00578 !> \page input_format_parallel_coag Input File Format: Parallel Coagulation Type 00579 !! 00580 !! The output type is specified by the parameter: 00581 !! - \b parallel_coag (string): type of parallel coagulation --- 00582 !! must be one of: \c local for only within-process 00583 !! coagulation or \c dist to have all processes perform 00584 !! coagulation globally, requesting particles from other 00585 !! processes as needed 00586 !! 00587 !! See also: 00588 !! - \ref spec_file_format --- the input file text format 00589 00590 call spec_file_read_string(file, 'parallel_coag', & 00591 parallel_coag_type_name) 00592 if (trim(parallel_coag_type_name) == 'local') then 00593 parallel_coag_type = PARALLEL_COAG_TYPE_LOCAL 00594 elseif (trim(parallel_coag_type_name) == 'dist') then 00595 parallel_coag_type = PARALLEL_COAG_TYPE_DIST 00596 else 00597 call spec_file_die_msg(494684716, file, & 00598 "Unknown parallel coagulation type: " & 00599 // trim(parallel_coag_type_name)) 00600 end if 00601 00602 end subroutine spec_file_read_parallel_coag_type 00603 00604 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00605 00606 end module pmc_run_part