PartMC  2.2.0
run_part.F90
Go to the documentation of this file.
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