ASPECT
ellipsoidal_chunk.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2022 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_ellipsoidal_chunk_h
23 #define _aspect_geometry_model_ellipsoidal_chunk_h
24 
28 #include <deal.II/grid/manifold.h>
29 
30 
31 namespace aspect
32 {
33  namespace GeometryModel
34  {
35  using namespace dealii;
36 
37  namespace internal
38  {
42  template <int dim>
43  class EllipsoidalChunkGeometry : public ChartManifold<dim,3,3>
44  {
45  public:
50  const double para_semi_major_axis_a,
51  const double para_eccentricity,
52  const double para_semi_minor_axis_b,
53  const double para_bottom_depth,
54  const std::vector<Point<2>> &para_corners);
55 
59  EllipsoidalChunkGeometry(const EllipsoidalChunkGeometry &other) = default;
60 
67  Point<3>
68  pull_back(const Point<3> &space_point) const override;
69 
74  virtual
75  Point<2>
76  pull_back(const Point<2> &space_point) const;
77 
84  Point<3>
85  push_forward(const Point<3> &chart_point) const override;
86 
91  Point<3> pull_back_ellipsoid (const Point<3> &x, const double semi_major_axis_a, const double eccentricity) const;
92 
97  Point<3> push_forward_ellipsoid (const Point<3> &phi_theta_d, const double semi_major_axis_a, const double eccentricity) const;
98 
102  std::unique_ptr<Manifold<dim,3>>
103  clone() const override;
104 
105  private:
110  Point<3> push_forward_topography (const Point<3> &phi_theta_d_hat) const;
111 
116  Point<3> pull_back_topography (const Point<3> &phi_theta_d) const;
117 
122 
123  const double semi_major_axis_a;
124  const double eccentricity;
125  const double semi_minor_axis_b;
126  const double bottom_depth;
127  const std::vector<Point<2>> corners;
128  };
129  }
130 
139  template <int dim>
140  class EllipsoidalChunk : public Interface<dim>, public SimulatorAccess<dim>
141  {
142  public:
146  void
147  initialize () override;
148 
149 
153  void
154  create_coarse_mesh(parallel::distributed::Triangulation<dim> &coarse_grid) const override;
155 
160  double
161  length_scale() const override;
162 
176  double
177  depth(const Point<dim> &position) const override;
178 
183  double
184  height_above_reference_surface(const Point<dim> &position) const override;
185 
189  Point<dim>
190  representative_point(const double depth) const override;
191 
197  bool
198  point_is_in_domain(const Point<dim> &point) const override;
199 
204  double
205  maximal_depth() const override;
206 
216  std::set<types::boundary_id>
217  get_used_boundary_indicators() const override;
218 
222  std::map<std::string,types::boundary_id>
223  get_symbolic_boundary_names_map() const override;
224 
229  aspect::Utilities::Coordinates::CoordinateSystem natural_coordinate_system() const override;
230 
238  std::array<double,dim> cartesian_to_natural_coordinates(const Point<dim> &position) const override;
239 
245  Point<dim> natural_to_cartesian_coordinates(const std::array<double,dim> &position) const override;
246 
250  static
251  void
252  declare_parameters(ParameterHandler &prm);
253 
258  void
259  parse_parameters(ParameterHandler &prm) override;
260 
264  double
265  get_radius(const Point<dim> &position) const;
266 
270  double
271  get_semi_minor_axis_b() const;
272 
276  double
277  get_semi_major_axis_a() const;
278 
279 
283  double
284  get_eccentricity() const;
285 
286 
294  const std::vector<Point<2>> &
295  get_corners() const;
296 
297 
302  get_manifold() const;
303 
304  private:
305  // Declare variables for reading in coordinates of the region of interest.
306  std::vector<Point<2>> corners;
308  double eccentricity;
314  double bottom_depth;
319  // Declare variables for subdividing
320  unsigned int EW_subdiv;
321  unsigned int NS_subdiv;
322  unsigned int depth_subdiv;
323 
324 
325 
337  std::unique_ptr<const internal::EllipsoidalChunkGeometry<dim>> manifold;
338 
339  void
340  set_boundary_ids(parallel::distributed::Triangulation<dim> &coarse_grid) const;
341  };
342  }
343 }
344 
345 
346 #endif
void declare_parameters(ParameterHandler &prm)
const InitialTopographyModel::Interface< dim > * topography
std::unique_ptr< const internal::EllipsoidalChunkGeometry< dim > > manifold
Definition: compat.h:42