ASPECT
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
box.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_box_h
23 #define _aspect_geometry_model_box_h
24 
27 
28 
29 namespace aspect
30 {
31  namespace GeometryModel
32  {
37  template <int dim>
38  class Box : public Interface<dim>, public SimulatorAccess<dim>
39  {
40  public:
46  void initialize () override;
47 
51  void create_coarse_mesh (parallel::distributed::Triangulation<dim> &coarse_grid) const override;
52 
57  virtual
58  Point<dim> get_extents () const;
59 
64  const std::array<unsigned int, dim> &
65  get_repetitions () const;
66 
71  virtual
72  Point<dim> get_origin () const;
73 
80  double length_scale () const override;
81 
98  double depth(const Point<dim> &position) const override;
99 
104  double height_above_reference_surface(const Point<dim> &position) const override;
105 
109  Point<dim> representative_point(const double depth) const override;
110 
114  double maximal_depth() const override;
115 
125  std::set<types::boundary_id>
126  get_used_boundary_indicators () const override;
127 
140  std::map<std::string,types::boundary_id>
141  get_symbolic_boundary_names_map () const override;
142 
147  std::set<std::pair<std::pair<types::boundary_id, types::boundary_id>, unsigned int>>
148  get_periodic_boundary_pairs () const override;
149 
156  void
157  adjust_positions_for_periodicity (Point<dim> &position,
158  const ArrayView<Point<dim>> &connected_positions = {},
159  const ArrayView<Tensor<1, dim>> &connected_velocities = {}) const override;
160 
166  bool
167  has_curved_elements() const override;
168 
174  bool
175  point_is_in_domain(const Point<dim> &point) const override;
176 
182 
188  std::array<double,dim> cartesian_to_natural_coordinates(const Point<dim> &position) const override;
189 
195  Point<dim> natural_to_cartesian_coordinates(const std::array<double,dim> &position) const override;
196 
200  static
201  void
202  declare_parameters (ParameterHandler &prm);
203 
207  void
208  parse_parameters (ParameterHandler &prm) override;
209 
210  private:
211 
215  void add_topography_to_mesh (typename parallel::distributed::Triangulation<dim> &grid) const;
216 
222  Point<dim> add_topography_to_point (const Point<dim> &x_y_z) const;
223 
227  Point<dim> extents;
228 
232  Point<dim> box_origin;
233 
237  std::array<bool, dim> periodic;
238 
242  std::array<unsigned int, dim> repetitions;
243  };
244  }
245 }
246 
247 
248 #endif
Point< dim > extents
Definition: box.h:227
std::array< double, dim > cartesian_to_natural_coordinates(const Point< dim > &position) const override
void adjust_positions_for_periodicity(Point< dim > &position, const ArrayView< Point< dim >> &connected_positions={}, const ArrayView< Tensor< 1, dim >> &connected_velocities={}) const override
double depth(const Point< dim > &position) const override
bool point_is_in_domain(const Point< dim > &point) const override
Point< dim > box_origin
Definition: box.h:232
std::map< std::string, types::boundary_id > get_symbolic_boundary_names_map() const override
void add_topography_to_mesh(typename parallel::distributed::Triangulation< dim > &grid) const
std::set< std::pair< std::pair< types::boundary_id, types::boundary_id >, unsigned int > > get_periodic_boundary_pairs() const override
bool has_curved_elements() const override
std::array< unsigned int, dim > repetitions
Definition: box.h:242
double length_scale() const override
aspect::Utilities::Coordinates::CoordinateSystem natural_coordinate_system() const override
static void declare_parameters(ParameterHandler &prm)
std::set< types::boundary_id > get_used_boundary_indicators() const override
double height_above_reference_surface(const Point< dim > &position) const override
const std::array< unsigned int, dim > & get_repetitions() const
Point< dim > representative_point(const double depth) const override
void create_coarse_mesh(parallel::distributed::Triangulation< dim > &coarse_grid) const override
Point< dim > natural_to_cartesian_coordinates(const std::array< double, dim > &position) const override
std::array< bool, dim > periodic
Definition: box.h:237
virtual Point< dim > get_origin() const
void parse_parameters(ParameterHandler &prm) override
void initialize() override
virtual Point< dim > get_extents() const
double maximal_depth() const override
Point< dim > add_topography_to_point(const Point< dim > &x_y_z) const