Rodin::Geometry::PolytopeTransformation class

Represents the transformation function of a simplex, taking reference coordinates to physical coordinates.

Let $ \tau $ denote a polytope in a triangulation $ \mathcal{T}_h $ which has an associated reference element $ K $ , This class represents the transformation $ x : K \subset \mathbb{R}^k \rightarrow \tau \subset \mathbb{R}^s $ of a reference point $ r $ into a physical point $ p $ :

\[ p = x(r) \]

Here, $ k $ and $ s $ represent the reference and physical dimensions, $ p \in \tau $ denotes the physical coordinates of the point, while $ x : K \rightarrow \tau $ represents the transformation taking reference coordinates $ r \in K $ , for a reference geometry $ K $ .

Derived classes

class IdentityTransformation final
Polytope Identity transformation.
template<class FE>
class IsoparametricTransformation
Polytope isoparametric transformation.

Public functions

auto transform(const Math::SpatialVector<Real>& rc) const -> Math::SpatialVector<Real>
Computes the physical coordinates of the given reference point.
auto jacobian(const Math::SpatialVector<Real>& rc) const -> Math::SpatialMatrix<Real>
Computes the Jacobian matrix of the transformation.
auto inverse(const Math::SpatialVector<Real>& pc) const -> Math::SpatialVector<Real>
Computes the reference coordinates of the given physical point.

Function documentation

Math::SpatialVector<Real> Rodin::Geometry::PolytopeTransformation::transform(const Math::SpatialVector<Real>& rc) const

Computes the physical coordinates of the given reference point.

Parameters
rc in Reference coordinates of the point.
Returns A vector of size $ s $ where $ s $ represents the physical dimension.

Given $ r \in K $ , computes the point:

\[ p = x(r) \]

in physical coordinates.

Math::SpatialMatrix<Real> Rodin::Geometry::PolytopeTransformation::jacobian(const Math::SpatialVector<Real>& rc) const

Computes the Jacobian matrix of the transformation.

Returns A matrix of dimensions $ s \times k $ where $ k $ represents the reference dimension and $ s $ represents the physical dimension.

Given $ r \in K $ , computes the Jacobian matrix:

\[ \mathbf{J}_x (r) = \begin{bmatrix} \dfrac{\partial x_1}{\partial r_1} & \ldots & \dfrac{\partial x_s}{\partial r_k}\\ \vdots & \ddots & \vdots\\ \dfrac{\partial x_s}{\partial r_1} & \ldots & \dfrac{\partial x_s}{\partial r_k} \end{bmatrix} , \]

for the given transformation $ x : K \rightarrow \tau $ .

Math::SpatialVector<Real> Rodin::Geometry::PolytopeTransformation::inverse(const Math::SpatialVector<Real>& pc) const

Computes the reference coordinates of the given physical point.

Parameters
pc in Physical coordinates of the point.

Given $ p \in \tau $ , computes the point:

\[ r = x^{-1}(p) \]

in reference coordinates.