![]() |
Reference documentation for deal.II version GIT 7b2de2f2f9 2023-09-24 11:00:02+00:00
|
#include <deal.II/non_matching/immersed_surface_quadrature.h>
Public Types | |
using | SubQuadrature = Quadrature< dim==0 ? 0 :dim - 1 > |
Public Member Functions | |
ImmersedSurfaceQuadrature ()=default | |
ImmersedSurfaceQuadrature (const std::vector< Point< dim >> &points, const std::vector< double > &weights, const std::vector< Tensor< 1, spacedim >> &normals) | |
void | clear () |
void | push_back (const Point< dim > &point, const double weight, const Tensor< 1, spacedim > &normal) |
const Tensor< 1, spacedim > & | normal_vector (const unsigned int i) const |
const std::vector< Tensor< 1, spacedim > > & | get_normal_vectors () const |
bool | operator== (const Quadrature< dim > &p) const |
void | initialize (const std::vector< Point< dim >> &points, const std::vector< double > &weights) |
unsigned int | size () const |
bool | empty () const |
const Point< dim > & | point (const unsigned int i) const |
const std::vector< Point< dim > > & | get_points () const |
double | weight (const unsigned int i) const |
const std::vector< double > & | get_weights () const |
std::size_t | memory_consumption () const |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
bool | is_tensor_product () const |
const std::array< Quadrature< 1 >, dim > & | get_tensor_basis () 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 ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
Protected Attributes | |
std::vector< Tensor< 1, spacedim > > | normals |
std::vector< Point< dim > > | quadrature_points |
std::vector< double > | weights |
bool | is_tensor_product_flag |
std::unique_ptr< std::array< Quadrature< 1 >, dim > > | tensor_basis |
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 int > | counter |
std::map< std::string, unsigned int > | counter_map |
std::vector< std::atomic< bool > * > | validity_pointers |
const std::type_info * | object_info |
Static Private Attributes | |
static std::mutex | mutex |
This class defines a quadrature formula to integrate over the intersection between an oriented surface, \(\hat{S}\), and a cell or face. The word "immersed" in the class name reflects that the surface may intersect the cell/face in an arbitrary way.
The spacedim template parameter of this class is the dimension that the (spacedim-1)-dimensional surface is embedded in: \(\hat{S} \subset \mathbb{R}^{\text{spacedim}}\). The dim parameter describes the dimension of the "object" that the surface intersects. That is, dim = spacedim corresponds to the surface intersecting a cell and dim = spacedim - 1 corresponds to the surface intersecting a face. The quadrature formula is described by a set of quadrature points, \(\hat{x}_q \in \mathbb{R}^{\text{dim}}\), weights, \(w_q\), and normalized surface normals, \(\hat{n}_q \in \mathbb{R}^{\text{spacedim}}\).
Consider first the case dim = spacedim. We typically want to compute integrals in real space. A surface, \(S\), intersecting a cell, \(K\), in real space can be mapped onto a surface, \(\hat{S}\), intersecting the unit cell, \(\hat{K}\). Thus an integral over \(S\cap K\) in real space can be transformed to an integral over \(\hat{S} \cap \hat{K}\) according to
\[ \int_{S\cap K} f dS = \int_{S\cap K} f |d\bar{S}| = \int_{\hat{S}\cap\hat{K}} f \circ F_{K} \det(J) |\left( J^{-1} \right )^T d\hat{S}|, \]
where \(F_K\) is the mapping from reference to real space and \(J\) is its Jacobian matrix. This transformation is possible since the continuous surface elements are vectors: \(d\bar{S}, d\hat{S} \in \mathbb{R}^{spacedim}\), which are parallel to the normals of \(S\) and \(\hat{S}\). That is, the normal is needed to do the transformation. Thus, in addition to storing points and weights, this quadrature stores also the normalized normal for each quadrature point. This can be viewed as storing a discrete surface element,
\[ \Delta \hat{S}_q \dealcoloneq w_q \hat{n}_q \approx d\hat{S}(\hat{x}_q), \]
for each quadrature point. The surface integral in real space would then be approximated as
\[ \int_{S\cap K} f dS \approx \sum_{q} f \left(F_{K}(\hat{x}_{q}) \right) \det(J_q) |\left( J_q^{-1} \right)^T \hat{n}_q| w_q. \]
When dim = spacedim - 1, this class represents a (spacedim-2)-dimensional integral. That is, if spacedim = 3 we have a line integral immersed in a face. Let \(\hat{r}(t)\), \(t \in [0,T]\) be an arc-length parameterizations of \(\hat{F}\cap \hat{S}\), i.e., the part of the surface that intersects the face in reference space. This means that \(\bar{r}(t) = F_K(\hat{r}(t))\) is a parameterization of \(S\cap F\). The transformation of the line integral now reads
\[ \int_{S\cap F} f dr = \int_{0}^T f(\bar{r}(t)) \left \|\frac{d\bar{r}}{dt} \right \| dt = \int_{0}^T f(F_K(\hat{r}(t))) \left \| J \frac{d\hat{r}}{dt} \right \| dt \approx \sum_{q} f \left(F_{K}(\hat{x}_{q}) \right) \|J(\hat{x}_q) \hat{t}_q \| w_q, \]
where \(\hat{t}_q = \frac{d\hat{r}}{dt}(x_q) \) is the tangent to the curve at \(\hat{x}_q\). This tangent can also be computed as \(t_q = \hat{n}_q \times \hat{n}_F / \| \hat{n}_q \times \hat{n}_F \|\) where \(\hat{n}_F\) is the face normal. It would be possible to compute the tangent by only knowing the normal to the curve in the face plane (i.e. the dim-dimensional normal). However, when these quadratures are used, the weak form typically involves the so-called conormal, which can not be computed without knowing the surface normal in \(\mathbb{R}^{\text{spacedim}}\). The conormal is the unit vector parallel to the projection of the face normal into the surface plane. This is the same as the normalized boundary form.
Definition at line 107 of file immersed_surface_quadrature.h.
|
inherited |
Define an alias for a quadrature that acts on an object of one dimension less. For cells, this would then be a face quadrature. A sub quadrature of a 0-dimensional quadrature is defined as still being 0-dimensional.
Definition at line 130 of file quadrature.h.
|
privateinherited |
The data type used in counter_map.
Definition at line 230 of file subscriptor.h.
|
privateinherited |
The iterator type used in counter_map.
Definition at line 235 of file subscriptor.h.
|
default |
Default constructor to initialize the quadrature with no quadrature points.
NonMatching::ImmersedSurfaceQuadrature< dim, spacedim >::ImmersedSurfaceQuadrature | ( | const std::vector< Point< dim >> & | points, |
const std::vector< double > & | weights, | ||
const std::vector< Tensor< 1, spacedim >> & | normals | ||
) |
Construct a quadrature formula from vectors of points, weights and surface normals. The points, weights and normals should be with respect to reference space, and the normals should be normalized.
Definition at line 22 of file immersed_surface_quadrature.cc.
|
inline |
Clears weights, points and normals vectors.
Definition at line 43 of file immersed_surface_quadrature.cc.
void NonMatching::ImmersedSurfaceQuadrature< dim, spacedim >::push_back | ( | const Point< dim > & | point, |
const double | weight, | ||
const Tensor< 1, spacedim > & | normal | ||
) |
Extend the given formula by an additional quadrature point. The point, weight and normal should be with respect to reference space, and the normal should be normalized.
This function exists since immersed quadrature rules can be rather complicated to construct. Often the construction is done by partitioning the cell into regions and constructing points on each region separately. This can make it cumbersome to create the quadrature from the constructor since all quadrature points have to be known at time of creation of the object.
Definition at line 54 of file immersed_surface_quadrature.cc.
const Tensor< 1, spacedim > & NonMatching::ImmersedSurfaceQuadrature< dim, spacedim >::normal_vector | ( | const unsigned int | i | ) | const |
Return a reference to the i
th surface normal.
Definition at line 70 of file immersed_surface_quadrature.cc.
const std::vector< Tensor< 1, spacedim > > & NonMatching::ImmersedSurfaceQuadrature< dim, spacedim >::get_normal_vectors |
Return a reference to the whole vector of normals.
Definition at line 81 of file immersed_surface_quadrature.cc.
|
inherited |
Test for equality of two quadratures.
Definition at line 317 of file quadrature.cc.
|
inherited |
Set the quadrature points and weights to the values provided in the arguments.
Definition at line 52 of file quadrature.cc.
|
inherited |
Number of quadrature points.
|
inherited |
Return if quadrature is empty.
|
inherited |
Return the i
th quadrature point.
|
inherited |
Return a reference to the whole array of quadrature points.
|
inherited |
Return the weight of the i
th quadrature point.
|
inherited |
Return a reference to the whole array of weights.
|
inherited |
Determine an estimate for the memory consumption (in bytes) of this object.
Definition at line 326 of file quadrature.cc.
|
inherited |
Write or read the data of this object to or from a stream for the purpose of serialization using the BOOST serialization library.
|
inherited |
This function returns true if the quadrature object is a tensor product of one-dimensional formulas and the quadrature points are sorted lexicographically.
|
inherited |
In case the quadrature formula is a tensor product, this function returns the dim
one-dimensional basis objects. Otherwise, calling this function is not allowed.
For dim
equal to one, we can not return the std::array as a const reference and have to return it by value. In this case, the array will always contain a single element (this
).
Definition at line 338 of file quadrature.cc.
|
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 136 of file subscriptor.cc.
|
inherited |
Unsubscribes a user from the object.
identifier
and the validity
pointer must be the same as the one supplied to subscribe(). Definition at line 156 of file subscriptor.cc.
|
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.
|
inlineinherited |
List the subscribers to the input stream
.
Definition at line 317 of file subscriptor.h.
|
inherited |
List the subscribers to deallog
.
Definition at line 204 of file subscriptor.cc.
|
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.
Definition at line 53 of file subscriptor.cc.
|
protected |
Vector of surface normals at each quadrature point.
Definition at line 167 of file immersed_surface_quadrature.h.
|
protectedinherited |
List of quadrature points. To be filled by the constructors of derived classes.
Definition at line 339 of file quadrature.h.
|
protectedinherited |
List of weights of the quadrature points. To be filled by the constructors of derived classes.
Definition at line 345 of file quadrature.h.
|
protectedinherited |
Indicates if this object represents quadrature formula that is a tensor product of one-dimensional formulas. This flag is set if dim==1 or the constructors taking a Quadrature<1> (and possibly a Quadrature<dim-1> object) is called. This implies that the quadrature points are sorted lexicographically.
Definition at line 354 of file quadrature.h.
|
protectedinherited |
Stores the one-dimensional tensor basis objects in case this object can be represented by a tensor product.
Definition at line 360 of file quadrature.h.
|
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 219 of file subscriptor.h.
|
mutableprivateinherited |
In this map, we count subscriptions for each different identification string supplied to subscribe().
Definition at line 225 of file subscriptor.h.
|
mutableprivateinherited |
In this vector, we store pointers to the validity bool in the SmartPointer objects that subscribe to this class.
Definition at line 241 of file subscriptor.h.
|
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 249 of file subscriptor.h.
|
staticprivateinherited |
A mutex used to ensure data consistency when printing out the list of subscribers.
Definition at line 271 of file subscriptor.h.