ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2015 - 2020 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 
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  {
288  };
289 
298  template <int dim>
299  class Interface
300  {
301  public:
306  virtual ~Interface ();
307 
312  virtual
313  void
314  initialize ();
315 
327  virtual
328  void
329  initialize_one_particle_property (const Point<dim> &position,
330  std::vector<double> &particle_properties) const;
331 
358  virtual
359  void
360  update_particle_property (const unsigned int data_position,
361  const Vector<double> &solution,
362  const std::vector<Tensor<1,dim> > &gradients,
363  typename ParticleHandler<dim>::particle_iterator &particle) const;
364 
394  virtual
395  void
396  update_one_particle_property (const unsigned int data_position,
397  const Point<dim> &position,
398  const Vector<double> &solution,
399  const std::vector<Tensor<1,dim> > &gradients,
400  const ArrayView<double> &particle_properties) const;
401 
414  virtual
416  need_update () const;
417 
426  virtual
428  get_needed_update_flags () const;
429 
445  virtual
447  late_initialization_mode () const;
448 
462  virtual
463  std::vector<std::pair<std::string, unsigned int> >
464  get_property_information() const = 0;
465 
466 
480  static
481  void
483 
490  virtual
491  void
492  parse_parameters (ParameterHandler &prm);
493  };
494 
495 
496 
502  template <int dim>
503  class Manager : public SimulatorAccess<dim>
504  {
505  public:
509  Manager ();
510 
514  ~Manager () override;
515 
520  void
521  initialize ();
522 
528  void
529  initialize_one_particle (typename ParticleHandler<dim>::particle_iterator &particle) const;
530 
537  std::vector<double>
538  initialize_late_particle (const Point<dim> &particle_location,
539  const ParticleHandler<dim> &particle_handler,
540  const Interpolator::Interface<dim> &interpolator,
542 
547  void
548  update_one_particle (typename ParticleHandler<dim>::particle_iterator &particle,
549  const Vector<double> &solution,
550  const std::vector<Tensor<1,dim> > &gradients) const;
551 
564  need_update () const;
565 
573  get_needed_update_flags () const;
574 
579  bool
580  plugin_name_exists(const std::string &name) const;
581 
591  bool
592  check_plugin_order(const std::string &first, const std::string &second) const;
593 
597  unsigned int get_plugin_index_by_name(const std::string &name) const;
598 
606  unsigned int
607  get_n_property_components () const;
608 
617  std::size_t
618  get_particle_size () const;
619 
627  get_data_info() const;
628 
637  unsigned int
638  get_property_component_by_name(const std::string &name) const;
639 
658  static
659  void
660  register_particle_property (const std::string &name,
661  const std::string &description,
662  void (*declare_parameters_function) (ParameterHandler &),
663  Property::Interface<dim> *(*factory_function) ());
664 
665 
669  static
670  void
672 
676  void
677  parse_parameters (ParameterHandler &prm);
678 
679 
689  static
690  void
691  write_plugin_graph (std::ostream &output_stream);
692 
693  private:
698  std::vector<std::string> plugin_names;
699 
704  std::list<std::unique_ptr<Interface<dim> > > property_list;
705 
711  };
712 
713 
720 #define ASPECT_REGISTER_PARTICLE_PROPERTY(classname,name,description) \
721  template class classname<2>; \
722  template class classname<3>; \
723  namespace ASPECT_REGISTER_PARTICLE_PROPERTY_ ## classname \
724  { \
725  aspect::internal::Plugins::RegisterHelper<aspect::Particle::Property::Interface<2>,classname<2> > \
726  dummy_ ## classname ## _2d (&aspect::Particle::Property::Manager<2>::register_particle_property, \
727  name, description); \
728  aspect::internal::Plugins::RegisterHelper<aspect::Particle::Property::Interface<3>,classname<3> > \
729  dummy_ ## classname ## _3d (&aspect::Particle::Property::Manager<3>::register_particle_property, \
730  name, description); \
731  }
732 
733  }
734  }
735 }
736 
737 #endif
std::string get_field_name_by_index(const unsigned int field_index) const
ParticlePropertyInformation property_information
Definition: interface.h:710
unsigned int get_components_by_field_index(const unsigned int field_index) const
Point< 2 > second
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
gradients
UpdateFlags
unsigned int get_field_index_by_name(const std::string &name) const
std::list< std::unique_ptr< Interface< dim > > > property_list
Definition: interface.h:704
unsigned int get_position_by_plugin_index(const unsigned int plugin_index) const
Point< 2 > first
Property
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)
#define DEAL_II_DEPRECATED
std::vector< std::string > plugin_names
Definition: interface.h:698
unsigned int get_position_by_field_index(const unsigned int field_index) const