ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2015 - 2022 by the authors of the ASPECT code.
3 
4  This file is part of ASPECT.
5 
6  ASPECT is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2, or (at your option)
9  any later version.
10 
11  ASPECT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with ASPECT; see the file LICENSE. If not see
18  <http://www.gnu.org/licenses/>.
19  */
20 
21 #ifndef _aspect_particle_property_interface_h
22 #define _aspect_particle_property_interface_h
23 
24 #include <aspect/global.h>
25 #include <aspect/plugins.h>
26 
29 
30 #include <deal.II/particles/particle.h>
31 #include <deal.II/particles/particle_handler.h>
32 #include <deal.II/particles/property_pool.h>
33 #include <deal.II/fe/fe_update_flags.h>
34 
35 #include <memory>
36 
37 namespace aspect
38 {
39  namespace Particle
40  {
41  namespace Property
42  {
43  using namespace ::Particles;
44 
80  {
81  public:
86 
99  ParticlePropertyInformation(const std::vector<std::vector<std::pair<std::string,unsigned int>>> &property_information);
100 
105  bool
106  fieldname_exists(const std::string &name) const;
107 
111  unsigned int
112  get_field_index_by_name(const std::string &name) const;
113 
117  std::string
118  get_field_name_by_index(const unsigned int field_index) const;
119 
124  unsigned int
125  get_position_by_field_name(const std::string &name) const;
126 
131  unsigned int
132  get_components_by_field_name(const std::string &name) const;
133 
139  unsigned int
140  get_position_by_field_index(const unsigned int field_index) const;
141 
146  unsigned int
147  get_components_by_field_index(const unsigned int field_index) const;
148 
154  unsigned int
155  get_position_by_plugin_index(const unsigned int plugin_index) const;
156 
161  unsigned int
162  get_components_by_plugin_index(const unsigned int plugin_index) const;
163 
168  unsigned int
169  get_fields_by_plugin_index(const unsigned int plugin_index) const;
170 
171 
175  unsigned int
176  n_plugins() const;
177 
181  unsigned int
182  n_fields() const;
183 
187  unsigned int
188  n_components() const;
189 
190  private:
194  std::vector<std::string> field_names;
195 
199  std::vector<unsigned int> components_per_field;
200 
205  std::vector<unsigned int> position_per_field;
206 
211  std::vector<unsigned int> fields_per_plugin;
212 
216  std::vector<unsigned int> components_per_plugin;
217 
222  std::vector<unsigned int> position_per_plugin;
223 
228  unsigned int number_of_components;
229 
233  unsigned int number_of_fields;
234 
238  unsigned int number_of_plugins;
239  };
240 
242  {
263  };
264 
271  {
297  };
298 
307  template <int dim>
308  class Interface
309  {
310  public:
315  virtual ~Interface () = default;
316 
321  virtual
322  void
323  initialize ();
324 
336  virtual
337  void
338  initialize_one_particle_property (const Point<dim> &position,
339  std::vector<double> &particle_properties) const;
340 
367  virtual
368  void
369  update_particle_property (const unsigned int data_position,
370  const Vector<double> &solution,
371  const std::vector<Tensor<1,dim>> &gradients,
372  typename ParticleHandler<dim>::particle_iterator &particle) const;
373 
403  virtual
404  void
405  update_one_particle_property (const unsigned int data_position,
406  const Point<dim> &position,
407  const Vector<double> &solution,
408  const std::vector<Tensor<1,dim>> &gradients,
409  const ArrayView<double> &particle_properties) const;
410 
423  virtual
425  need_update () const;
426 
435  virtual
436  UpdateFlags
437  get_needed_update_flags () const;
438 
453  virtual
455  late_initialization_mode () const;
456 
470  virtual
471  std::vector<std::pair<std::string, unsigned int>>
472  get_property_information() const = 0;
473 
474 
488  static
489  void
490  declare_parameters (ParameterHandler &prm);
491 
498  virtual
499  void
500  parse_parameters (ParameterHandler &prm);
501  };
502 
512  template <int dim>
513  class IntegratorProperties : public Interface<dim>
514  {
515  public:
521  void
522  initialize_one_particle_property (const Point<dim> &position,
523  std::vector<double> &particle_properties) const override;
524 
532  std::vector<std::pair<std::string, unsigned int>>
533  get_property_information() const override;
534 
539  void
540  parse_parameters (ParameterHandler &prm) override;
541 
542  private:
548  };
549 
550 
551 
557  template <int dim>
558  class Manager : public SimulatorAccess<dim>
559  {
560  public:
564  Manager ();
565 
569  ~Manager () override;
570 
575  void
576  initialize ();
577 
583  void
584  initialize_one_particle (typename ParticleHandler<dim>::particle_iterator &particle) const;
585 
592  std::vector<double>
593  initialize_late_particle (const Point<dim> &particle_location,
594  const ParticleHandler<dim> &particle_handler,
595  const Interpolator::Interface<dim> &interpolator,
596  const typename parallel::distributed::Triangulation<dim>::active_cell_iterator &cell = typename parallel::distributed::Triangulation<dim>::active_cell_iterator()) const;
597 
602  void
603  update_one_particle (typename ParticleHandler<dim>::particle_iterator &particle,
604  const Vector<double> &solution,
605  const std::vector<Tensor<1,dim>> &gradients) const;
606 
619  need_update () const;
620 
627  UpdateFlags
628  get_needed_update_flags () const;
629 
634  bool
635  plugin_name_exists(const std::string &name) const;
636 
646  bool
647  check_plugin_order(const std::string &first, const std::string &second) const;
648 
652  unsigned int get_plugin_index_by_name(const std::string &name) const;
653 
663  template <typename ParticlePropertyType,
664  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,ParticlePropertyType>::value>>
665  bool
666  has_matching_property () const;
667 
679  template <typename ParticlePropertyType,
680  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,ParticlePropertyType>::value>>
681  const ParticlePropertyType &
682  get_matching_property () const;
683 
691  unsigned int
692  get_n_property_components () const;
693 
702  std::size_t
703  get_particle_size () const;
704 
712  get_data_info() const;
713 
722  unsigned int
723  get_property_component_by_name(const std::string &name) const;
724 
743  static
744  void
745  register_particle_property (const std::string &name,
746  const std::string &description,
747  void (*declare_parameters_function) (ParameterHandler &),
748  std::unique_ptr<Property::Interface<dim>> (*factory_function) ());
749 
750 
754  static
755  void
756  declare_parameters (ParameterHandler &prm);
757 
761  void
762  parse_parameters (ParameterHandler &prm);
763 
764 
774  static
775  void
776  write_plugin_graph (std::ostream &output_stream);
777 
778  private:
783  std::vector<std::string> plugin_names;
784 
789  std::list<std::unique_ptr<Interface<dim>>> property_list;
790 
796  };
797 
798  /* -------------------------- inline and template functions ---------------------- */
799 
800 
801  template <int dim>
802  template <typename ParticlePropertyType, typename>
803  inline
804  bool
806  {
807  for (const auto &p : property_list)
808  if (Plugins::plugin_type_matches<ParticlePropertyType>(*p))
809  return true;
810 
811  return false;
812  }
813 
814 
815  template <int dim>
816  template <typename ParticlePropertyType, typename>
817  inline
818  const ParticlePropertyType &
820  {
821  AssertThrow(has_matching_property<ParticlePropertyType> (),
822  ExcMessage("You asked Particle::Property::Manager::has_matching_property() for a "
823  "particle property of type <" + boost::core::demangle(typeid(ParticlePropertyType).name()) + "> "
824  "that could not be found in the current model. Activate this "
825  "particle property in the input file."));
826 
827  typename std::vector<std::unique_ptr<Interface<dim>>>::const_iterator property;
828  for (const auto &p : property_list)
829  if (Plugins::plugin_type_matches<ParticlePropertyType>(*p))
830  return Plugins::get_plugin_as_type<ParticlePropertyType>(*p);
831 
832  // We will never get here, because we had the Assert above. Just to avoid warnings.
833  return Plugins::get_plugin_as_type<ParticlePropertyType>(*(*property));
834  }
835 
836 
843 #define ASPECT_REGISTER_PARTICLE_PROPERTY(classname,name,description) \
844  template class classname<2>; \
845  template class classname<3>; \
846  namespace ASPECT_REGISTER_PARTICLE_PROPERTY_ ## classname \
847  { \
848  aspect::internal::Plugins::RegisterHelper<aspect::Particle::Property::Interface<2>,classname<2>> \
849  dummy_ ## classname ## _2d (&aspect::Particle::Property::Manager<2>::register_particle_property, \
850  name, description); \
851  aspect::internal::Plugins::RegisterHelper<aspect::Particle::Property::Interface<3>,classname<3>> \
852  dummy_ ## classname ## _3d (&aspect::Particle::Property::Manager<3>::register_particle_property, \
853  name, description); \
854  }
855 
856  }
857  }
858 }
859 
860 #endif
std::string get_field_name_by_index(const unsigned int field_index) const
StructuredDataLookup< dim > DEAL_II_DEPRECATED
ParticlePropertyInformation property_information
Definition: interface.h:795
unsigned int get_components_by_field_index(const unsigned int field_index) const
void write_plugin_graph(std::ostream &output_stream)
unsigned int get_components_by_field_name(const std::string &name) const
unsigned int get_position_by_field_name(const std::string &name) const
bool fieldname_exists(const std::string &name) const
unsigned int get_field_index_by_name(const std::string &name) const
std::list< std::unique_ptr< Interface< dim > > > property_list
Definition: interface.h:789
const ParticlePropertyType & get_matching_property() const
Definition: interface.h:819
unsigned int get_position_by_plugin_index(const unsigned int plugin_index) const
unsigned int get_components_by_plugin_index(const unsigned int plugin_index) const
unsigned int get_fields_by_plugin_index(const unsigned int plugin_index) const
void declare_parameters(ParameterHandler &prm)
std::vector< std::string > plugin_names
Definition: interface.h:783
unsigned int get_position_by_field_index(const unsigned int field_index) const