15#ifndef dealii_manifold_lib_h
16#define dealii_manifold_lib_h
82template <
int dim,
int spacedim = dim>
97 virtual std::unique_ptr<Manifold<dim, spacedim>>
98 clone()
const override;
148 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
149 "Access the center with get_center() instead.")
164 static
Tensor<1, spacedim>
261template <
int dim,
int spacedim = dim>
276 virtual std::unique_ptr<Manifold<dim, spacedim>>
277 clone()
const override;
289 const double w)
const override;
354 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
355 "Access the center with get_center() instead.")
374 const
ArrayView<const
double> &distances,
375 const
ArrayView<const
double> &weights) const;
418template <
int dim,
int spacedim = dim>
445 virtual std::unique_ptr<Manifold<dim, spacedim>>
446 clone()
const override;
562template <
int dim,
int spacedim = dim>
581 virtual std::unique_ptr<Manifold<dim, spacedim>>
582 clone()
const override;
667template <
int dim,
int spacedim = dim,
int chartdim = dim>
732 const double h = 1e-8);
742 virtual std::unique_ptr<Manifold<dim, spacedim>>
743 clone()
const override;
877 virtual std::unique_ptr<Manifold<dim, 3>>
878 clone()
const override;
1055template <
int dim,
int spacedim = dim>
1072 virtual std::unique_ptr<Manifold<dim, spacedim>>
1073 clone()
const override;
1145 std::array<unsigned int, 20>
1242 std::vector<internal::MappingQImplementation::
1243 InverseQuadraticApproximation<dim, spacedim>>
1255template <
int dim,
int spacedim>
1264template <
int dim,
int spacedim>
1273template <
int dim,
int spacedim>
1282template <
int dim,
int spacedim>
1291template <
int dim,
int spacedim>
1300template <
int dim,
int spacedim>
1309template <
int dim,
int spacedim>
1318template <
int dim,
int spacedim>
ChartManifold(const Tensor< 1, chartdim > &periodicity=Tensor< 1, chartdim >())
virtual void get_new_points(const ArrayView< const Point< spacedim > > &surrounding_points, const Table< 2, double > &weights, ArrayView< Point< spacedim > > new_points) const override
virtual Tensor< 1, spacedim > get_tangent_vector(const Point< spacedim > &x1, const Point< spacedim > &x2) const override
virtual Point< spacedim > get_intermediate_point(const Point< spacedim > &p1, const Point< spacedim > &p2, const double w) const override
virtual Point< spacedim > get_new_point(const ArrayView< const Point< spacedim > > &surrounding_points, const ArrayView< const double > &weights) const override
const Tensor< 1, spacedim > dxn
const Tensor< 1, spacedim > normal_direction
CylindricalManifold(const unsigned int axis=0, const double tolerance=1e-10)
const Point< spacedim > point_on_axis
const Tensor< 1, spacedim > direction
double get_tolerance() const
const Point< spacedim > & get_point_on_axis() const
const Tensor< 1, spacedim > & get_direction() const
const Tensor< 1, spacedim > & get_major_axis_direction() const
const Point< spacedim > & get_center() const
double get_eccentricity() const
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
const Point< spacedim > center
const double eccentricity
virtual Point< spacedim > push_forward(const Point< spacedim > &chart_point) const override
EllipticalManifold(const Point< spacedim > ¢er, const Tensor< 1, spacedim > &major_axis_direction, const double eccentricity)
virtual DerivativeForm< 1, spacedim, spacedim > push_forward_gradient(const Point< spacedim > &chart_point) const override
virtual Point< spacedim > pull_back(const Point< spacedim > &space_point) const override
static Tensor< 1, spacedim > get_periodicity()
const Tensor< 1, spacedim > direction
SmartPointer< const Function< spacedim >, FunctionManifold< dim, spacedim, chartdim > > pull_back_function
const FunctionParser< spacedim >::ConstMap const_map
const std::string chart_vars
const std::string pull_back_expression
virtual DerivativeForm< 1, chartdim, spacedim > push_forward_gradient(const Point< chartdim > &chart_point) const override
SmartPointer< const Function< chartdim >, FunctionManifold< dim, spacedim, chartdim > > push_forward_function
const std::string space_vars
virtual Point< spacedim > push_forward(const Point< chartdim > &chart_point) const override
const double finite_difference_step
FunctionManifold(const Function< chartdim > &push_forward_function, const Function< spacedim > &pull_back_function, const Tensor< 1, chartdim > &periodicity=Tensor< 1, chartdim >(), const double tolerance=1e-10)
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
const std::string push_forward_expression
virtual Point< chartdim > pull_back(const Point< spacedim > &space_point) const override
std::map< std::string, double > ConstMap
static std::string default_variable_names()
virtual void get_normals_at_vertices(const typename Triangulation< dim, spacedim >::face_iterator &face, FaceVertexNormals &face_vertex_normals) const
std::array< Tensor< 1, spacedim >, GeometryInfo< dim >::vertices_per_face > FaceVertexNormals
PolarManifold(const Point< spacedim > center=Point< spacedim >())
static Tensor< 1, spacedim > get_periodicity()
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
virtual DerivativeForm< 1, spacedim, spacedim > push_forward_gradient(const Point< spacedim > &chart_point) const override
const Point< spacedim > center
const Point< spacedim > p_center
virtual Tensor< 1, spacedim > normal_vector(const typename Triangulation< dim, spacedim >::face_iterator &face, const Point< spacedim > &p) const override
virtual Point< spacedim > push_forward(const Point< spacedim > &chart_point) const override
virtual Point< spacedim > pull_back(const Point< spacedim > &space_point) const override
const Point< spacedim > & get_center() const
const Point< spacedim > & get_center() const
const Point< spacedim > p_center
const Point< spacedim > center
const PolarManifold< spacedim > polar_manifold
std::pair< double, Tensor< 1, spacedim > > guess_new_point(const ArrayView< const Tensor< 1, spacedim > > &directions, const ArrayView< const double > &distances, const ArrayView< const double > &weights) const
SphericalManifold(const Point< spacedim > center=Point< spacedim >())
void do_get_new_points(const ArrayView< const Point< spacedim > > &surrounding_points, const ArrayView< const double > &weights, ArrayView< Point< spacedim > > new_points) const
virtual Point< 3 > pull_back(const Point< 3 > &p) const override
double get_centerline_radius() const
static const int spacedim
virtual DerivativeForm< 1, 3, 3 > push_forward_gradient(const Point< 3 > &chart_point) const override
virtual Point< 3 > push_forward(const Point< 3 > &chart_point) const override
static const int chartdim
TorusManifold(const double centerline_radius, const double inner_radius)
double get_inner_radius() const
virtual std::unique_ptr< Manifold< dim, 3 > > clone() const override
const Triangulation< dim, spacedim > * triangulation
boost::signals2::connection clear_signal
DerivativeForm< 1, dim, spacedim > push_forward_gradient(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &chart_point, const Point< spacedim > &pushed_forward_chart_point) const
std::vector< internal::MappingQImplementation::InverseQuadraticApproximation< dim, spacedim > > quadratic_approximation
void initialize(const Triangulation< dim, spacedim > &triangulation)
Triangulation< dim, spacedim >::cell_iterator compute_chart_points(const ArrayView< const Point< spacedim > > &surrounding_points, ArrayView< Point< dim > > chart_points) const
std::array< unsigned int, 20 > get_possible_cells_around_points(const ArrayView< const Point< spacedim > > &surrounding_points) const
FlatManifold< dim > chart_manifold
Point< dim > pull_back(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_guess) const
virtual void get_new_points(const ArrayView< const Point< spacedim > > &surrounding_points, const Table< 2, double > &weights, ArrayView< Point< spacedim > > new_points) const override
virtual Point< spacedim > get_new_point(const ArrayView< const Point< spacedim > > &surrounding_points, const ArrayView< const double > &weights) const override
TransfiniteInterpolationManifold()
std::vector< bool > coarse_cell_is_flat
Point< spacedim > push_forward(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &chart_point) const
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const override
InverseQuadraticApproximation(const ArrayView< const Point< spacedim > > &real_support_points, const std::vector< Point< dim > > &unit_support_points)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face_iterator
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator