ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2018 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>
28 #include <array>
29 #include <aspect/utilities.h>
31 
32 #include <set>
33 
34 
35 namespace aspect
36 {
47  namespace GeometryModel
48  {
49  using namespace dealii;
50 
58  template <int dim>
59  class Interface
60  {
61  public:
66  virtual ~Interface();
67 
73  virtual void initialize ();
74 
78  virtual
79  void create_coarse_mesh (parallel::distributed::Triangulation<dim> &coarse_grid) const = 0;
80 
94  virtual
95  double length_scale () const = 0;
96 
132  virtual
133  double depth(const Point<dim> &position) const = 0;
134 
144  virtual
145  double height_above_reference_surface(const Point<dim> &position) const = 0;
146 
152  Utilities::NaturalCoordinate<dim>
153  cartesian_to_other_coordinates(const Point<dim> &position,
154  const Utilities::Coordinates::CoordinateSystem &coordinate_system) const;
155 
159  virtual
160  aspect::Utilities::Coordinates::CoordinateSystem natural_coordinate_system() const = 0;
161 
169  virtual
170  std::array<double,dim> cartesian_to_natural_coordinates(const Point<dim> &position) const;
171 
177  virtual
178  Point<dim> natural_to_cartesian_coordinates(const std::array<double,dim> &position) const;
179 
197  virtual
198  Point<dim> representative_point(const double depth) const = 0;
199 
203  virtual
204  double maximal_depth() const = 0;
205 
206 
212  virtual
213  std::set<types::boundary_id>
214  get_used_boundary_indicators () const = 0;
215 
251  virtual
252  std::map<std::string,types::boundary_id>
253  get_symbolic_boundary_names_map () const;
254 
270  translate_symbolic_boundary_name_to_id (const std::string &name) const;
271 
286  std::vector<types::boundary_id>
287  translate_symbolic_boundary_names_to_ids (const std::vector<std::string> &names) const;
288 
298  std::string
299  translate_id_to_symbol_name (const types::boundary_id boundary_id) const;
300 
307  virtual
308  std::set< std::pair< std::pair<types::boundary_id, types::boundary_id>, unsigned int> >
309  get_periodic_boundary_pairs () const;
310 
317  virtual
318  bool
319  has_curved_elements() const;
320 
325  virtual
326  bool
327  point_is_in_domain(const Point<dim> &p) const = 0;
328 
335  static
336  void
338 
345  virtual
346  void
347  parse_parameters (ParameterHandler &prm);
348 
349  };
350 
351 
352 
368  template <int dim>
369  void
370  register_geometry_model (const std::string &name,
371  const std::string &description,
372  void (*declare_parameters_function) (ParameterHandler &),
373  Interface<dim> *(*factory_function) ());
374 
385  template <int dim>
386  Interface<dim> *
388 
394  template <int dim>
395  void
397 
407  template <int dim>
408  void
409  write_plugin_graph (std::ostream &output_stream);
410 
411 
412 
420 #define ASPECT_REGISTER_GEOMETRY_MODEL(classname,name,description) \
421  template class classname<2>; \
422  template class classname<3>; \
423  namespace ASPECT_REGISTER_GEOMETRY_MODEL_ ## classname \
424  { \
425  aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<2>,classname<2> > \
426  dummy_ ## classname ## _2d (&aspect::GeometryModel::register_geometry_model<2>, \
427  name, description); \
428  aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<3>,classname<3> > \
429  dummy_ ## classname ## _3d (&aspect::GeometryModel::register_geometry_model<3>, \
430  name, description); \
431  }
432  }
433 }
434 
435 
436 #endif
void register_geometry_model(const std::string &name, const std::string &description, void(*declare_parameters_function)(ParameterHandler &), Interface< dim > *(*factory_function)())
void declare_parameters(ParameterHandler &prm)
Definition: compat.h:38
void write_plugin_graph(std::ostream &output_stream)
Interface< dim > * create_geometry_model(ParameterHandler &prm)
unsigned int boundary_id