Reference documentation for deal.II version GIT relicensing-136-gb80d0be4af 2024-03-18 08:20:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Static Public Member Functions | Private Types | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
Manifold< dim, spacedim > Class Template Referenceabstract

#include <deal.II/grid/manifold.h>

Inheritance diagram for Manifold< dim, spacedim >:
Inheritance graph
[legend]

Public Types

using FaceVertexNormals = std::array< Tensor< 1, spacedim >, GeometryInfo< dim >::vertices_per_face >
 

Public Member Functions

virtual ~Manifold () override=default
 
virtual std::unique_ptr< Manifold< dim, spacedim > > clone () const =0
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 
Computing the location of points.
virtual Point< spacedim > get_intermediate_point (const Point< spacedim > &p1, const Point< spacedim > &p2, const double w) const
 
virtual Point< spacedim > get_new_point (const ArrayView< const Point< spacedim > > &surrounding_points, const ArrayView< const double > &weights) const
 
virtual void get_new_points (const ArrayView< const Point< spacedim > > &surrounding_points, const Table< 2, double > &weights, ArrayView< Point< spacedim > > new_points) const
 
virtual Point< spacedim > project_to_manifold (const ArrayView< const Point< spacedim > > &surrounding_points, const Point< spacedim > &candidate) const
 
virtual Point< spacedim > get_new_point_on_line (const typename Triangulation< dim, spacedim >::line_iterator &line) const
 
virtual Point< spacedim > get_new_point_on_quad (const typename Triangulation< dim, spacedim >::quad_iterator &quad) const
 
virtual Point< spacedim > get_new_point_on_hex (const typename Triangulation< dim, spacedim >::hex_iterator &hex) const
 
Point< spacedim > get_new_point_on_face (const typename Triangulation< dim, spacedim >::face_iterator &face) const
 
Point< spacedim > get_new_point_on_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
 
Computing tangent vectors
virtual Tensor< 1, spacedim > get_tangent_vector (const Point< spacedim > &x1, const Point< spacedim > &x2) const
 
Computing normal vectors
virtual Tensor< 1, spacedim > normal_vector (const typename Triangulation< dim, spacedim >::face_iterator &face, const Point< spacedim > &p) const
 
virtual void get_normals_at_vertices (const typename Triangulation< dim, spacedim >::face_iterator &face, FaceVertexNormals &face_vertex_normals) const
 
Subscriptor functionality

Classes derived from Subscriptor provide a facility to subscribe to this object. This is mostly used by the SmartPointer class.

void subscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
void unsubscribe (std::atomic< bool > *const validity, const std::string &identifier="") const
 
unsigned int n_subscriptions () const
 
template<typename StreamType >
void list_subscribers (StreamType &stream) const
 
void list_subscribers () const
 

Static Public Member Functions

static ::ExceptionBaseExcInUse (int arg1, std::string arg2, std::string arg3)
 
static ::ExceptionBaseExcNoSubscriber (std::string arg1, std::string arg2)
 

Private Types

using map_value_type = decltype(counter_map)::value_type
 
using map_iterator = decltype(counter_map)::iterator
 

Private Member Functions

void check_no_subscribers () const noexcept
 

Private Attributes

std::atomic< unsigned intcounter
 
std::map< std::string, unsigned intcounter_map
 
std::vector< std::atomic< bool > * > validity_pointers
 
const std::type_info * object_info
 

Static Private Attributes

static std::mutex mutex
 

Detailed Description

template<int dim, int spacedim = dim>
class Manifold< dim, spacedim >

Manifolds are used to describe the geometry of boundaries of domains as well as the geometry of the interior. Manifold objects are therefore associated with cells, faces, and/or edges, either by direct user action or, if a user program does not do this explicitly, a default manifold object is used.

Manifolds are best understood by using the language of differential geometry, but their common uses can be easily described simply through examples. An exhaustive discussion of how, where, and why this class is used is provided in the geometry paper.

Common use case: Creating a new vertex

In the most essential use of manifolds, manifold descriptions are used to create a "point between other points". For example, when a triangulation creates a new vertex on a cell, face, or edge, it determines the new vertex' coordinates through the following function call:

...
Point<spacedim> new_vertex = manifold.get_new_point (points,weights);
...

Here, points is a collection of points in spacedim dimension, and a collection of corresponding weights. The points in this context will then be the vertices of the cell, face, or edge, and the weights are typically one over the number of points when a new midpoint of the cell, face, or edge is needed. Derived classes then will implement the Manifold::get_new_point() function in a way that computes the location of this new point. In the simplest case, for example in the FlatManifold class, the function simply computes the arithmetic average (with given weights) of the given points. However, other classes do something differently; for example, the SphericalManifold class, which is used to describe domains that form (part of) the sphere, will ensure that, given the two vertices of an edge at the boundary, the new returned point will lie on the grand circle that connects the two points, rather than choosing a point that is half-way between the two points in \({\mathbb R}^d\).

Note
Unlike almost all other cases in the library, we here interpret the points to be in real space, not on the reference cell.

Manifold::get_new_point() has a default implementation that can simplify this process somewhat: Internally, the function calls the Manifold::get_intermediate_point() to compute pair-wise intermediate points. Internally the Manifold::get_intermediate_point() calls the Manifold::project_to_manifold() function after computing the convex combination of the given points. This allows derived classes to only overload Manifold::project_to_manifold() for simple situations. This is often useful when describing manifolds that are embedded in higher dimensional space, e.g., the surface of a sphere. In those cases, the desired new point may be computed simply by the (weighted) average of the provided points, projected back out onto the sphere.

Common use case: Computing tangent vectors

The second use of this class is in computing directions on domains and boundaries. For example, we may need to compute the normal vector to a face in order to impose the no-flow boundary condition \(\mathbf u \cdot \mathbf n = 0\) (see the VectorTools::compute_no_normal_flux_constraints() as an example). Similarly, we may need normal vectors in the computation of the normal component of the gradient of the numerical solution in order to compute the jump in the gradient of the solution in error estimators (see, for example, the KellyErrorEstimator class).

To make this possible, the Manifold class provides a member function (to be implemented by derived classes) that computes a "vector tangent to the manifold at one point, in direction of another point" via the Manifold::get_tangent_vector() function. For example, in 2d, one would use this function with the two vertices of an edge at the boundary to compute a "tangential" vector along the edge, and then get the normal vector by rotation by 90 degrees. In 3d, one would compute the two vectors "tangential" to the two edges of a boundary face adjacent to a boundary vertex, and then take the cross product of these two to obtain a vector normal to the boundary.

For reasons that are more difficult to understand, these direction vectors are normalized in a very specific way, rather than to have unit norm. See the documentation of Manifold::get_tangent_vector(), as well as below, for more information.

In the simplest case (namely, the FlatManifold class), these tangent vectors are just the difference vector between the two given points. However, in more complicated (and more interesting) cases, the direction may be different. For example, for the SphericalManifold case, if the two given points lie on a common grand circle around the origin, then the tangent vector will be tangential to the grand circle, rather than pointing straight from one point to the other.

A unified description

The "real" way to understand what this class does is to see it in the framework of differential geometry. More specifically, differential geometry is fundamentally based on the assumption that two sufficiently close points are connected via a line of "shortest distance". This line is called a "geodesic", and it is selected from all other lines that connect the two points by the property that it is shortest if distances are measured in terms of the "metric" that describes a manifold. To give examples, recall that the geodesics of a flat manifold (implemented in the FlatManifold class) are simply the straight lines connecting two points, whereas for spherical manifolds (see the SphericalManifold class) geodesics between two points of same distance are the grand circles, and are in general curved lines when connecting two lines of different distance from the origin.

In the following discussion, and for the purposes of implementing the current class, the concept of "metrics" that is so fundamental to differential geometry is no longer of great importance to us. Rather, everything can simply be described by postulating the existence of geodesics connecting points on a manifold.

Given geodesics, the operations discussed in the previous two sections can be described in a more formal way. In essence, they rely on the fact that we can assume that a geodesic is parameterized by a "time" like variable \(t\) so that \(\mathbf s(t)\) describes the curve and so that \(\mathbf s(0)\) is the location of the first and \(\mathbf s(1)\) the location of the second point. Furthermore, \(\mathbf s(t)\) traces out the geodesic at constant speed, covering equal distance in equal time (as measured by the metric). Note that this parameterization uses time, not arc length to denote progress along the geodesic.

In this picture, computing a mid-point between points \(\mathbf x_1\) and \(\mathbf x_2\), with weights \(w_1\) and \(w_2=1-w_1\), simply requires computing the point \(\mathbf s(w_1)\). Computing a new point as a weighted average of more than two points can be done by considering pairwise geodesics, finding suitable points on the geodetic between the first two points, then on the geodetic between this new point and the third given point, etc.

Likewise, the "tangential" vector described above is simply the velocity vector, \(\mathbf s'(t)\), evaluated at one of the end points of a geodesic (i.e., at \(t=0\) or \(t=1\)). In the case of a flat manifold, the geodesic is simply the straight line connecting two points, and the velocity vector is just the connecting vector in that case. On the other hand, for two points on a spherical manifold, the geodesic is a grand circle, and the velocity vector is tangent to the spherical surface.

Note that if we wanted to, we could use this to compute the length of the geodesic that connects two points \(\mathbf x_1\) and \(\mathbf x_2\) by computing \(\int_0^1 \|\mathbf s'(t)\| dt\) along the geodesic that connects them, but this operation will not be of use to us in practice. One could also conceive computing the direction vector using the "new point" operation above, using the formula \(\mathbf s'(0)=\lim_{w\rightarrow 0} \frac{\mathbf s(w)-\mathbf s(0)}{w}\) where all we need to do is compute the new point \(\mathbf s(w)\) with weights \(w\) and \(1-w\) along the geodesic connecting \(\mathbf x_1\) and \(\mathbf x_2\). The default implementation of the function does this, by evaluating the quotient for a small but finite weight \(w\). In practice, however, it is almost always possible to explicitly compute the direction vector, i.e., without the need to numerically approximate the limit process, and derived classes should do so.

Definition at line 285 of file manifold.h.

Member Typedef Documentation

◆ FaceVertexNormals

template<int dim, int spacedim = dim>
using Manifold< dim, spacedim >::FaceVertexNormals = std::array<Tensor<1, spacedim>, GeometryInfo<dim>::vertices_per_face>

Type keeping information about the normals at the vertices of a face of a cell. Thus, there are GeometryInfo<dim>::vertices_per_face normal vectors, that define the tangent spaces of the boundary at the vertices. Note that the vectors stored in this object are not required to be normalized, nor to actually point outward, as one often will only want to check for orthogonality to define the tangent plane; if a function requires the normals to be normalized, then it must do so itself.

For obvious reasons, this type is not useful in 1d.

Definition at line 305 of file manifold.h.

◆ map_value_type

using Subscriptor::map_value_type = decltype(counter_map)::value_type
privateinherited

The data type used in counter_map.

Definition at line 229 of file subscriptor.h.

◆ map_iterator

using Subscriptor::map_iterator = decltype(counter_map)::iterator
privateinherited

The iterator type used in counter_map.

Definition at line 234 of file subscriptor.h.

Constructor & Destructor Documentation

◆ ~Manifold()

template<int dim, int spacedim = dim>
virtual Manifold< dim, spacedim >::~Manifold ( )
overridevirtualdefault

Destructor. Does nothing here, but needs to be declared virtual to make class hierarchies derived from this class possible.

Member Function Documentation

◆ clone()

template<int dim, int spacedim = dim>
virtual std::unique_ptr< Manifold< dim, spacedim > > Manifold< dim, spacedim >::clone ( ) const
pure virtual

◆ get_intermediate_point()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::get_intermediate_point ( const Point< spacedim > &  p1,
const Point< spacedim > &  p2,
const double  w 
) const
virtual

Return an intermediate point between two given points. Overloading this function allows the default pair-wise reduction implementation of the method get_new_point() that takes a Quadrature object as input to work properly.

An implementation of this function should returns a parametric curve on the manifold, joining the points p1 and p2, with parameter w in the interval [0,1]. In particular get_intermediate_point(p1, p2, 0.0) should return p1 and get_intermediate_point(p1, p2, 1.0) should return p2.

In its default implementation, this function calls the project_to_manifold() method with the convex combination of p1 and p2. User classes can get away by simply implementing the project_to_manifold() method.

Reimplemented in ChartManifold< dim, spacedim, chartdim >, ChartManifold< dim, 3, 3 >, ChartManifold< dim, dim, 3 >, ChartManifold< dim, dim, dim >, ChartManifold< dim, spacedim, 1 >, ChartManifold< dim, spacedim, 2 >, ChartManifold< dim, spacedim_A+spacedim_B, chartdim_A+chartdim_B >, ChartManifold< dim_A, spacedim_A, chartdim_A >, ChartManifold< dim_B, spacedim_B, chartdim_B >, and SphericalManifold< dim, spacedim >.

◆ get_new_point()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::get_new_point ( const ArrayView< const Point< spacedim > > &  surrounding_points,
const ArrayView< const double > &  weights 
) const
virtual

Return the point which shall become the new vertex surrounded by the given points surrounding_points. weights contains appropriate weights for the surrounding points according to which the manifold determines the new point's position.

In its default implementation it uses a pair-wise reduction of the points by calling the function get_intermediate_point() on the first two points, then on the resulting point and the next, until all points in the vector have been taken into account. User classes can get away by simply implementing the get_intermediate_point() function. Notice that by default the get_intermediate_point() function calls the project_to_manifold() function with the convex combination of its arguments. For simple situations you may get away by implementing only the project_to_manifold() function.

Reimplemented in FlatManifold< dim, spacedim >, FlatManifold< chartdim, chartdim >, FlatManifold< dim, dim >, FlatManifold< dim, spacedim >, ChartManifold< dim, spacedim, chartdim >, ChartManifold< dim, 3, 3 >, ChartManifold< dim, dim, 3 >, ChartManifold< dim, dim, dim >, ChartManifold< dim, spacedim, 1 >, ChartManifold< dim, spacedim, 2 >, ChartManifold< dim, spacedim_A+spacedim_B, chartdim_A+chartdim_B >, ChartManifold< dim_A, spacedim_A, chartdim_A >, ChartManifold< dim_B, spacedim_B, chartdim_B >, CylindricalManifold< dim, spacedim >, TransfiniteInterpolationManifold< dim, spacedim >, and SphericalManifold< dim, spacedim >.

◆ get_new_points()

template<int dim, int spacedim = dim>
virtual void Manifold< dim, spacedim >::get_new_points ( const ArrayView< const Point< spacedim > > &  surrounding_points,
const Table< 2, double > &  weights,
ArrayView< Point< spacedim > >  new_points 
) const
virtual

Compute a new set of points that interpolate between the given points surrounding_points. weights is a table with as many columns as surrounding_points.size(). The number of rows in weights must match the length of new_points.

In its default implementation, this function simply calls get_new_point() on each row of weights and writes those points into the output array new_points. However, this function is more efficient if multiple new points need to be generated like in MappingQ and the manifold does expensive transformations between a chart space and the physical space, such as ChartManifold. For this function, the surrounding points need to be transformed back to the chart sparse only once, rather than for every call to get_new_point(). If efficiency is not important, you may get away by implementing only the get_new_point() function.

The implementation does not allow for surrounding_points and new_points to point to the same array, so make sure to pass different objects into the function.

Reimplemented in FlatManifold< dim, spacedim >, FlatManifold< chartdim, chartdim >, FlatManifold< dim, dim >, FlatManifold< dim, spacedim >, ChartManifold< dim, spacedim, chartdim >, ChartManifold< dim, 3, 3 >, ChartManifold< dim, dim, 3 >, ChartManifold< dim, dim, dim >, ChartManifold< dim, spacedim, 1 >, ChartManifold< dim, spacedim, 2 >, ChartManifold< dim, spacedim_A+spacedim_B, chartdim_A+chartdim_B >, ChartManifold< dim_A, spacedim_A, chartdim_A >, ChartManifold< dim_B, spacedim_B, chartdim_B >, SphericalManifold< dim, spacedim >, and TransfiniteInterpolationManifold< dim, spacedim >.

◆ project_to_manifold()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::project_to_manifold ( const ArrayView< const Point< spacedim > > &  surrounding_points,
const Point< spacedim > &  candidate 
) const
virtual

Given a point which lies close to the given manifold, it modifies it and projects it to manifold itself.

This class is used by the default implementation of the function get_new_point() and should be implemented by derived classes. The default implementation simply throws an exception if called.

If your manifold is simple, you could implement this function only, and the default behavior should work out of the box.

Reimplemented in FlatManifold< dim, spacedim >, FlatManifold< chartdim, chartdim >, FlatManifold< dim, dim >, FlatManifold< dim, spacedim >, OpenCASCADE::NormalProjectionManifold< dim, spacedim >, OpenCASCADE::DirectionalProjectionManifold< dim, spacedim >, and OpenCASCADE::NormalToMeshProjectionManifold< dim, spacedim >.

◆ get_new_point_on_line()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::get_new_point_on_line ( const typename Triangulation< dim, spacedim >::line_iterator &  line) const
virtual

Backward compatibility interface. Return the point which shall become the new middle vertex of the two children of a regular line. In 2d, this line is a line at the boundary, while in 3d, it is bounding a face at the boundary (the lines therefore is also on the boundary).

The default implementation of this function passes its argument to the Manifolds::get_default_points_and_weights() function, and then calls the Manifold<dim,spacedim>::get_new_point() function. User derived classes can overload Manifold<dim,spacedim>::get_new_point() or Manifold<dim,spacedim>::project_to_manifold(), which is called by the default implementation of Manifold<dim,spacedim>::get_new_point().

◆ get_new_point_on_quad()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::get_new_point_on_quad ( const typename Triangulation< dim, spacedim >::quad_iterator &  quad) const
virtual

Backward compatibility interface. Return the point which shall become the common point of the four children of a quad at the boundary in three or more spatial dimensions. This function therefore is only useful in at least three dimensions and should not be called for lower dimensions.

This function is called after the four lines bounding the given quad are refined, so you may want to use the information provided by quad->line(i)->child(j), i=0...3, j=0,1.

The default implementation of this function passes its argument to the Manifolds::get_default_points_and_weights() function, and then calls the Manifold<dim,spacedim>::get_new_point() function. User derived classes can overload Manifold<dim,spacedim>::get_new_point() or Manifold<dim,spacedim>::project_to_manifold(), which is called by the default implementation of Manifold<dim,spacedim>::get_new_point().

◆ get_new_point_on_hex()

template<int dim, int spacedim = dim>
virtual Point< spacedim > Manifold< dim, spacedim >::get_new_point_on_hex ( const typename Triangulation< dim, spacedim >::hex_iterator &  hex) const
virtual

Backward compatibility interface. Return the point which shall become the common point of the eight children of a hex in three or spatial dimensions. This function therefore is only useful in at least three dimensions and should not be called for lower dimensions.

This function is called after the all the bounding objects of the given hex are refined, so you may want to use the information provided by hex->quad(i)->line(j)->child(k), i=0...5, j=0...3, k=0,1.

The default implementation of this function passes its argument to the Manifolds::get_default_points_and_weights() function, and then calls the Manifold<dim,spacedim>::get_new_point() function. User derived classes can overload Manifold<dim,spacedim>::get_new_point() or Manifold<dim,spacedim>::project_to_manifold(), which is called by the default implementation of Manifold<dim,spacedim>::get_new_point().

◆ get_new_point_on_face()

template<int dim, int spacedim = dim>
Point< spacedim > Manifold< dim, spacedim >::get_new_point_on_face ( const typename Triangulation< dim, spacedim >::face_iterator &  face) const

Backward compatibility interface. Depending on dim=2 or dim=3 this function calls the get_new_point_on_line or the get_new_point_on_quad function. It throws an exception for dim=1. This wrapper allows dimension independent programming.

◆ get_new_point_on_cell()

template<int dim, int spacedim = dim>
Point< spacedim > Manifold< dim, spacedim >::get_new_point_on_cell ( const typename Triangulation< dim, spacedim >::cell_iterator &  cell) const

Backward compatibility interface. Depending on dim=1, dim=2 or dim=3 this function calls the get_new_point_on_line, get_new_point_on_quad or the get_new_point_on_hex function. This wrapper allows dimension independent programming.

◆ get_tangent_vector()

template<int dim, int spacedim = dim>
virtual Tensor< 1, spacedim > Manifold< dim, spacedim >::get_tangent_vector ( const Point< spacedim > &  x1,
const Point< spacedim > &  x2 
) const
virtual

Return a vector that, at \(\mathbf x_1\), is tangential to the geodesic that connects two points \(\mathbf x_1,\mathbf x_2\). The geodesic is the shortest line between these two points, where "shortest" is defined via a metric specific to a particular implementation of this class in a derived class. For example, in the case of a FlatManifold, the shortest line between two points is just the straight line, and in this case the tangent vector is just the difference \(\mathbf d=\mathbf x_2-\mathbf x_1\). On the other hand, for a manifold that describes a surface embedded in a higher dimensional space (e.g., the surface of a sphere), then the tangent vector is tangential to the surface, and consequently may point in a different direction than the straight line that connects the two points.

While tangent vectors are often normalized to unit length, the vectors returned by this function are normalized as described in the introduction of this class. Specifically, if \(\mathbf s(t)\) traces out the geodesic between the two points where \(\mathbf x_1 = \mathbf s(0)\) and \(\mathbf x_2 = \mathbf s(1)\), then the returned vector must equal \(\mathbf s'(0)\). In other words, the norm of the returned vector also encodes, in some sense, the length of the geodesic because a curve \(\mathbf s(t)\) must move "faster" if the two points it connects between arguments \(t=0\) and \(t=1\) are farther apart.

The default implementation of this function approximates \(\mathbf s'(0) \approx \frac{\mathbf s(\epsilon)-\mathbf x_1}{\epsilon}\) for a small value of \(\epsilon\), and the evaluation of \(\mathbf s(\epsilon)\) is done by calling get_new_point(). If possible, derived classes should override this function by an implementation of the exact derivative.

Parameters
x1The first point that describes the geodesic, and the one at which the "direction" is to be evaluated.
x2The second point that describes the geodesic.
Returns
A "direction" vector tangential to the geodesic.

Reimplemented in FlatManifold< dim, spacedim >, FlatManifold< chartdim, chartdim >, FlatManifold< dim, dim >, FlatManifold< dim, spacedim >, ChartManifold< dim, spacedim, chartdim >, ChartManifold< dim, 3, 3 >, ChartManifold< dim, dim, 3 >, ChartManifold< dim, dim, dim >, ChartManifold< dim, spacedim, 1 >, ChartManifold< dim, spacedim, 2 >, ChartManifold< dim, spacedim_A+spacedim_B, chartdim_A+chartdim_B >, ChartManifold< dim_A, spacedim_A, chartdim_A >, ChartManifold< dim_B, spacedim_B, chartdim_B >, and SphericalManifold< dim, spacedim >.

◆ normal_vector()

template<int dim, int spacedim = dim>
virtual Tensor< 1, spacedim > Manifold< dim, spacedim >::normal_vector ( const typename Triangulation< dim, spacedim >::face_iterator &  face,
const Point< spacedim > &  p 
) const
virtual

Return the normal vector to a face embedded in this manifold, at the point p. It is not required that the normals actually point outward from the domain even if the face iterator given points to a face on the boundary of the domain. If p is not in fact on the surface, but only close-by, try to return something reasonable, for example the normal vector at the surface point closest to p. (The point p will in fact not normally lie on the actual surface, but rather be a quadrature point mapped by some polynomial mapping; the mapped surface, however, will not usually coincide with the actual surface.)

This function only makes sense if dim==spacedim because otherwise there is no unique normal vector but in fact a (spacedim-dim+1)-dimensional tangent space of vectors that are all both normal to the face and normal to the dim-dimensional surface that lives in spacedim-dimensional space. For example, think of a two-dimensional mesh that covers a two-dimensional surface in three-dimensional space. In that case, each face (edge) is one-dimensional, and there are two linearly independent vectors that are both normal to the edge: one is normal to the edge and tangent to the surface (intuitively, that would be the one that points from the current cell to the neighboring one, if the surface was locally flat), and the other one is rooted in the edge but points perpendicular to the surface (which is also perpendicular to the edge that lives within the surface). Thus, because there are no obviously correct semantics for this function if spacedim is greater than dim, the function will simply throw an error in that situation.

The face iterator gives an indication which face this function is supposed to compute the normal vector for. This is useful if the boundary of the domain is composed of different nondifferential pieces (for example when using the FlatManifold class to approximate a geometry that is completely described by the coarse mesh, with piecewise (bi-)linear components between the vertices, but where the boundary may have a kink at the vertices itself).

Note
In 2d, the default implementation of this function computes the normal vector by taking the tangent direction from p to the further one of the two vertices that make up an edge, and then rotates it outward (with respect to the coordinate system of the edge) by 90 degrees. In 3d, the default implementation is more complicated, aiming at avoiding problems with numerical round-off for points close to one of the vertices, and avoiding tangent directions that are linearly dependent.

Reimplemented in FlatManifold< dim, spacedim >, FlatManifold< chartdim, chartdim >, FlatManifold< dim, dim >, FlatManifold< dim, spacedim >, PolarManifold< dim, spacedim >, and SphericalManifold< dim, spacedim >.

◆ get_normals_at_vertices()

template<int dim, int spacedim = dim>
virtual void Manifold< dim, spacedim >::get_normals_at_vertices ( const typename Triangulation< dim, spacedim >::face_iterator &  face,
FaceVertexNormals face_vertex_normals 
) const
virtual

Compute the normal vectors to the boundary at each vertex of the given face embedded in the Manifold. It is not required that the normal vectors be normed somehow. Neither is it required that the normals actually point outward.

This function is needed to compute data for C1 mappings. The default implementation calls normal_vector() on each vertex.

Note that when computing normal vectors at a vertex where the boundary is not differentiable, you have to make sure that you compute the one-sided limits, i.e. limit with respect to points inside the given face.

◆ subscribe()

void Subscriptor::subscribe ( std::atomic< bool > *const  validity,
const std::string &  identifier = "" 
) const
inherited

Subscribes a user of the object by storing the pointer validity. The subscriber may be identified by text supplied as identifier.

Definition at line 135 of file subscriptor.cc.

◆ unsubscribe()

void Subscriptor::unsubscribe ( std::atomic< bool > *const  validity,
const std::string &  identifier = "" 
) const
inherited

Unsubscribes a user from the object.

Note
The identifier and the validity pointer must be the same as the one supplied to subscribe().

Definition at line 155 of file subscriptor.cc.

◆ n_subscriptions()

unsigned int Subscriptor::n_subscriptions ( ) const
inlineinherited

Return the present number of subscriptions to this object. This allows to use this class for reference counted lifetime determination where the last one to unsubscribe also deletes the object.

Definition at line 300 of file subscriptor.h.

◆ list_subscribers() [1/2]

template<typename StreamType >
void Subscriptor::list_subscribers ( StreamType &  stream) const
inlineinherited

List the subscribers to the input stream.

Definition at line 317 of file subscriptor.h.

◆ list_subscribers() [2/2]

void Subscriptor::list_subscribers ( ) const
inherited

List the subscribers to deallog.

Definition at line 203 of file subscriptor.cc.

◆ serialize()

template<class Archive >
void Subscriptor::serialize ( Archive &  ar,
const unsigned int  version 
)
inlineinherited

Read or write the data of this object to or from a stream for the purpose of serialization using the BOOST serialization library.

This function does not actually serialize any of the member variables of this class. The reason is that what this class stores is only who subscribes to this object, but who does so at the time of storing the contents of this object does not necessarily have anything to do with who subscribes to the object when it is restored. Consequently, we do not want to overwrite the subscribers at the time of restoring, and then there is no reason to write the subscribers out in the first place.

Definition at line 309 of file subscriptor.h.

◆ check_no_subscribers()

void Subscriptor::check_no_subscribers ( ) const
privatenoexceptinherited

Check that there are no objects subscribing to this object. If this check passes then it is safe to destroy the current object. It this check fails then this function will either abort or print an error message to deallog (by using the AssertNothrow mechanism), but will not throw an exception.

Note
Since this function is just a consistency check it does nothing in release mode.
If this function is called when there is an uncaught exception then, rather than aborting, this function prints an error message to the standard error stream and returns.

Definition at line 52 of file subscriptor.cc.

Member Data Documentation

◆ counter

std::atomic<unsigned int> Subscriptor::counter
mutableprivateinherited

Store the number of objects which subscribed to this object. Initially, this number is zero, and upon destruction it shall be zero again (i.e. all objects which subscribed should have unsubscribed again).

The creator (and owner) of an object is counted in the map below if HE manages to supply identification.

We use the mutable keyword in order to allow subscription to constant objects also.

This counter may be read from and written to concurrently in multithreaded code: hence we use the std::atomic class template.

Definition at line 218 of file subscriptor.h.

◆ counter_map

std::map<std::string, unsigned int> Subscriptor::counter_map
mutableprivateinherited

In this map, we count subscriptions for each different identification string supplied to subscribe().

Definition at line 224 of file subscriptor.h.

◆ validity_pointers

std::vector<std::atomic<bool> *> Subscriptor::validity_pointers
mutableprivateinherited

In this vector, we store pointers to the validity bool in the SmartPointer objects that subscribe to this class.

Definition at line 240 of file subscriptor.h.

◆ object_info

const std::type_info* Subscriptor::object_info
mutableprivateinherited

Pointer to the typeinfo object of this object, from which we can later deduce the class name. Since this information on the derived class is neither available in the destructor, nor in the constructor, we obtain it in between and store it here.

Definition at line 248 of file subscriptor.h.

◆ mutex

std::mutex Subscriptor::mutex
staticprivateinherited

A mutex used to ensure data consistency when accessing the mutable members of this class. This lock is used in the subscribe() and unsubscribe() functions, as well as in list_subscribers().

Definition at line 271 of file subscriptor.h.


The documentation for this class was generated from the following file: