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