21 #ifndef _aspect_particle_world_h 22 #define _aspect_particle_world_h 26 #include <deal.II/particles/particle.h> 27 #include <deal.II/particles/particle_accessor.h> 28 #include <deal.II/particles/particle_iterator.h> 29 #include <deal.II/particles/particle_handler.h> 30 #include <deal.II/particles/property_pool.h> 32 #include <deal.II/matrix_free/fe_point_evaluation.h> 42 #include <deal.II/base/timer.h> 43 #include <deal.II/base/array_view.h> 45 #include <boost/serialization/unique_ptr.hpp> 50 struct SimulatorSignals;
55 using namespace ::Particles;
93 reinit(
const typename DoFHandler<dim>::active_cell_iterator &cell,
94 const ArrayView<Point<dim>> &positions,
95 const ArrayView<double> &solution_values,
96 const UpdateFlags update_flags) = 0;
104 virtual void get_solution(
const unsigned int evaluation_point,
105 Vector<double> &solution) = 0;
113 virtual void get_gradients(
const unsigned int evaluation_point,
114 std::vector<Tensor<1, dim>> &gradients) = 0;
120 virtual FEPointEvaluation<dim, dim> &
121 get_velocity_or_fluid_velocity_evaluator(
const bool use_fluid_velocity) = 0;
126 virtual NonMatching::MappingInfo<dim> &
127 get_mapping_info() = 0;
134 std::unique_ptr<internal::SolutionEvaluators<dim>>
136 const UpdateFlags update_flags);
175 get_property_manager()
const;
182 const Particles::ParticleHandler<dim> &
183 get_particle_handler()
const;
196 Particles::ParticleHandler<dim> &
197 get_particle_handler();
217 void copy_particle_handler (
const Particles::ParticleHandler<dim> &from_particle_handler,
218 Particles::ParticleHandler<dim> &to_particle_handler)
const;
225 void backup_particles ();
233 void restore_particles ();
239 void setup_initial_state ();
247 get_interpolator()
const;
252 void generate_particles();
256 void initialize_particles();
265 void advance_timestep();
278 types::particle_index n_global_particles()
const;
298 #if DEAL_II_VERSION_GTE(9,6,0) 300 cell_weight(
const typename parallel::distributed::Triangulation<dim>::cell_iterator &cell,
301 const CellStatus status);
304 cell_weight(
const typename parallel::distributed::Triangulation<dim>::cell_iterator &cell,
305 const typename parallel::distributed::Triangulation<dim>::CellStatus status);
311 void update_particles();
316 template <
class Archive>
317 void serialize (Archive &ar,
const unsigned int version);
324 save (std::ostringstream &os)
const;
331 load (std::istringstream &is);
345 parse_parameters (ParameterHandler &prm);
353 remove_particles = 0x1,
356 remove_and_add_particles = remove_particles | add_particles
453 std::map<types::subdomain_id, unsigned int>
454 get_subdomain_id_to_neighbor_map()
const;
462 apply_particle_per_cell_bounds();
468 void advect_particles();
474 local_initialize_particles(
const typename ParticleHandler<dim>::particle_iterator &begin_particle,
475 const typename ParticleHandler<dim>::particle_iterator &end_particle);
481 local_update_particles(
const typename DoFHandler<dim>::active_cell_iterator &cell,
482 const typename ParticleHandler<dim>::particle_iterator &begin_particle,
483 const typename ParticleHandler<dim>::particle_iterator &end_particle);
494 local_update_particles(
const typename DoFHandler<dim>::active_cell_iterator &cell,
495 const typename ParticleHandler<dim>::particle_iterator &begin_particle,
496 const typename ParticleHandler<dim>::particle_iterator &end_particle,
508 local_advect_particles(
const typename DoFHandler<dim>::active_cell_iterator &cell,
509 const typename ParticleHandler<dim>::particle_iterator &begin_particle,
510 const typename ParticleHandler<dim>::particle_iterator &end_particle);
526 local_advect_particles(
const typename DoFHandler<dim>::active_cell_iterator &cell,
527 const typename ParticleHandler<dim>::particle_iterator &begin_particle,
528 const typename ParticleHandler<dim>::particle_iterator &end_particle,
537 ParticleHandler<dim> &particle_handler,
538 const bool connect_to_checkpoint_signals =
true)
const;
544 template <
class Archive>
void serialize(Archive &ar, const unsigned int version)
std::unique_ptr< internal::SolutionEvaluators< dim > > construct_solution_evaluators(const SimulatorAccess< dim > &simulator_access, const UpdateFlags update_flags)
unsigned int min_particles_per_cell
std::unique_ptr< Particles::ParticleHandler< dim > > particle_handler
unsigned int max_particles_per_cell
std::unique_ptr< Interpolator::Interface< dim > > interpolator
unsigned int particle_weight
std::unique_ptr< Generator::Interface< dim > > generator
void declare_parameters(ParameterHandler &prm)
std::unique_ptr< Integrator::Interface< dim > > integrator
std::unique_ptr< Property::Manager< dim > > property_manager
ParticleLoadBalancing::Kind particle_load_balancing
bool update_ghost_particles
Particles::ParticleHandler< dim > particle_handler_backup