22 #ifndef _aspect_geometry_model_interface_h 23 #define _aspect_geometry_model_interface_h 26 #include <deal.II/base/parameter_handler.h> 27 #include <deal.II/distributed/tria.h> 47 namespace GeometryModel
66 virtual ~Interface() =
default;
79 void create_coarse_mesh (parallel::distributed::Triangulation<dim> &coarse_grid)
const = 0;
95 double length_scale ()
const = 0;
135 double depth(
const Point<dim> &position)
const = 0;
147 double height_above_reference_surface(
const Point<dim> &position)
const = 0;
154 Utilities::NaturalCoordinate<dim>
155 cartesian_to_other_coordinates(
const Point<dim> &position,
172 std::array<double,dim> cartesian_to_natural_coordinates(
const Point<dim> &position)
const;
180 Point<dim> natural_to_cartesian_coordinates(
const std::array<double,dim> &position)
const;
200 Point<dim> representative_point(
const double depth)
const = 0;
211 double maximal_depth()
const = 0;
220 std::set<types::boundary_id>
221 get_used_boundary_indicators ()
const = 0;
259 std::map<std::string,types::boundary_id>
260 get_symbolic_boundary_names_map ()
const;
277 translate_symbolic_boundary_name_to_id (
const std::string &name)
const;
293 std::vector<types::boundary_id>
294 translate_symbolic_boundary_names_to_ids (
const std::vector<std::string> &names)
const;
306 translate_id_to_symbol_name (
const types::boundary_id boundary_id)
const;
315 std::set<std::pair<std::pair<types::boundary_id, types::boundary_id>,
unsigned int>>
316 get_periodic_boundary_pairs ()
const;
341 adjust_positions_for_periodicity (Point<dim> &position,
342 const ArrayView<Point<dim>> &connected_positions = {},
343 const ArrayView<Tensor<1,dim>> &connected_velocities = {})
const;
353 has_curved_elements()
const;
361 point_is_in_domain(
const Point<dim> &p)
const = 0;
381 parse_parameters (ParameterHandler &prm);
391 make_periodicity_constraints(
const DoFHandler<dim> &dof_handler,
392 AffineConstraints<double> &constraints)
const;
415 const std::string &description,
416 void (*declare_parameters_function) (ParameterHandler &),
417 std::unique_ptr<Interface<dim>> (*factory_function) ());
430 std::unique_ptr<Interface<dim>>
464 #define ASPECT_REGISTER_GEOMETRY_MODEL(classname,name,description) \ 465 template class classname<2>; \ 466 template class classname<3>; \ 467 namespace ASPECT_REGISTER_GEOMETRY_MODEL_ ## classname \ 469 aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<2>,classname<2>> \ 470 dummy_ ## classname ## _2d (&aspect::GeometryModel::register_geometry_model<2>, \ 471 name, description); \ 472 aspect::internal::Plugins::RegisterHelper<aspect::GeometryModel::Interface<3>,classname<3>> \ 473 dummy_ ## classname ## _3d (&aspect::GeometryModel::register_geometry_model<3>, \ 474 name, description); \
std::unique_ptr< Interface< dim > > create_geometry_model(ParameterHandler &prm)
void declare_parameters(ParameterHandler &prm)
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)())