PartMC 2.1.2
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_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