ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2024 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 
22 #ifndef _aspect_geometry_model_interface_h
23 #define _aspect_geometry_model_interface_h
24 
25 #include <aspect/plugins.h>
26 #include <deal.II/base/parameter_handler.h>
27 #include <deal.II/distributed/tria.h>
28 #include <array>
29 #include <aspect/utilities.h>
31 
32 #include <set>
33 
34 
35 namespace aspect
36 {
47  namespace GeometryModel
48  {
56  template <int dim>
57  class Interface : public Plugins::InterfaceBase
58  {
59  public:
63  virtual
64  void create_coarse_mesh (parallel::distributed::Triangulation<dim> &coarse_grid) const = 0;
65 
79  virtual
80  double length_scale () const = 0;
81 
119  virtual
120  double depth(const Point<dim> &position) const = 0;
121 
131  virtual
132  double height_above_reference_surface(const Point<dim> &position) const = 0;
133 
139  Utilities::NaturalCoordinate<dim>
140  cartesian_to_other_coordinates(const Point<dim> &position,
141  const Utilities::Coordinates::CoordinateSystem &coordinate_system) const;
142 
146  virtual
148 
156  virtual
157  std::array<double,dim> cartesian_to_natural_coordinates(const Point<dim> &position) const;
158 
164  virtual
165  Point<dim> natural_to_cartesian_coordinates(const std::array<double,dim> &position) const;
166 
184  virtual
185  Point<dim> representative_point(const double depth) const = 0;
186 
195  virtual
196  double maximal_depth() const = 0;
197 
198 
204  virtual
205  std::set<types::boundary_id>
206  get_used_boundary_indicators () const = 0;
207 
243  virtual
244  std::map<std::string,types::boundary_id>
246 
261  types::boundary_id
262  translate_symbolic_boundary_name_to_id (const std::string &name) const;
263 
278  std::vector<types::boundary_id>
279  translate_symbolic_boundary_names_to_ids (const std::vector<std::string> &names) const;
280 
290  std::string
291  translate_id_to_symbol_name (const types::boundary_id boundary_id) const;
292 
299  virtual
300  std::set<std::pair<std::pair<types::boundary_id, types::boundary_id>, unsigned int>>
302 
324  virtual
325  void
326  adjust_positions_for_periodicity (Point<dim> &position,
327  const ArrayView<Point<dim>> &connected_positions = {},
328  const ArrayView<Tensor<1,dim>> &connected_velocities = {}) const;
329 
336  virtual
337  bool
338  has_curved_elements() const;
339 
344  virtual
345  bool
346  point_is_in_domain(const Point<dim> &p) const = 0;
347 
354  virtual
355  void
356  make_periodicity_constraints(const DoFHandler<dim> &dof_handler,
357  AffineConstraints<double> &constraints) const;
358  };
359 
360 
361 
377  template <int dim>
378  void
379  register_geometry_model (const std::string &name,
380  const std::string &description,
381  void (*declare_parameters_function) (ParameterHandler &),
382  std::unique_ptr<Interface<dim>> (*factory_function) ());
383 
394  template <int dim>
395  std::unique_ptr<Interface<dim>>
396  create_geometry_model (ParameterHandler &prm);
397 
403  template <int dim>
404  void
405  declare_parameters (ParameterHandler &prm);
406 
416  template <int dim>
417  void
418  write_plugin_graph (std::ostream &output_stream);
419 
420 
421 
429 #define ASPECT_REGISTER_GEOMETRY_MODEL(classname,name,description) \
430  template class classname<2>; \
431  template class classname<3>; \
432  namespace ASPECT_REGISTER_GEOMETRY_MODEL_ ## classname \
433  { \
434  aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<2>,classname<2>> \
435  dummy_ ## classname ## _2d (&aspect::GeometryModel::register_geometry_model<2>, \
436  name, description); \
437  aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<3>,classname<3>> \
438  dummy_ ## classname ## _3d (&aspect::GeometryModel::register_geometry_model<3>, \
439  name, description); \
440  }
441  }
442 }
443 
444 
445 #endif
virtual void make_periodicity_constraints(const DoFHandler< dim > &dof_handler, AffineConstraints< double > &constraints) const
std::string translate_id_to_symbol_name(const types::boundary_id boundary_id) const
virtual Point< dim > natural_to_cartesian_coordinates(const std::array< double, dim > &position) const
virtual std::set< std::pair< std::pair< types::boundary_id, types::boundary_id >, unsigned int > > get_periodic_boundary_pairs() const
std::unique_ptr< Interface< dim > > create_geometry_model(ParameterHandler &prm)
virtual bool point_is_in_domain(const Point< dim > &p) const =0
virtual std::set< types::boundary_id > get_used_boundary_indicators() const =0
virtual double maximal_depth() const =0
virtual void create_coarse_mesh(parallel::distributed::Triangulation< dim > &coarse_grid) const =0
virtual double height_above_reference_surface(const Point< dim > &position) const =0
virtual double length_scale() const =0
virtual aspect::Utilities::Coordinates::CoordinateSystem natural_coordinate_system() const =0
Definition: compat.h:59
Utilities::NaturalCoordinate< dim > cartesian_to_other_coordinates(const Point< dim > &position, const Utilities::Coordinates::CoordinateSystem &coordinate_system) const
virtual std::map< std::string, types::boundary_id > get_symbolic_boundary_names_map() const
virtual std::array< double, dim > cartesian_to_natural_coordinates(const Point< dim > &position) const
void write_plugin_graph(std::ostream &output_stream)
void register_geometry_model(const std::string &name, const std::string &description, void(*declare_parameters_function)(ParameterHandler &), std::unique_ptr< Interface< dim >>(*factory_function)())
std::vector< types::boundary_id > translate_symbolic_boundary_names_to_ids(const std::vector< std::string > &names) const
types::boundary_id translate_symbolic_boundary_name_to_id(const std::string &name) const
virtual void adjust_positions_for_periodicity(Point< dim > &position, const ArrayView< Point< dim >> &connected_positions={}, const ArrayView< Tensor< 1, dim >> &connected_velocities={}) const
virtual double depth(const Point< dim > &position) const =0
virtual Point< dim > representative_point(const double depth) const =0
static void declare_parameters(ParameterHandler &prm)
virtual bool has_curved_elements() const