LehrFEM++ 1.0.0
A simple Finite Element Library for teaching
Loading...
Searching...
No Matches
quad_o2.h
1
9#ifndef QUAD_O2_H
10#define QUAD_O2_H
11
12#include "geometry_interface.h"
13
14namespace lf::geometry {
15
29class QuadO2 : public Geometry {
30 public:
36 explicit QuadO2(Eigen::Matrix<double, Eigen::Dynamic, 8> coords);
37
38 [[nodiscard]] dim_t DimLocal() const override { return 2; }
39 [[nodiscard]] dim_t DimGlobal() const override { return coords_.rows(); }
40 [[nodiscard]] base::RefEl RefEl() const override {
41 return base::RefEl::kQuad();
42 }
43
44 [[nodiscard]] Eigen::MatrixXd Global(
45 const Eigen::MatrixXd &local) const override;
46 [[nodiscard]] Eigen::MatrixXd Jacobian(
47 const Eigen::MatrixXd &local) const override;
48 [[nodiscard]] Eigen::MatrixXd JacobianInverseGramian(
49 const Eigen::MatrixXd &local) const override;
50 [[nodiscard]] Eigen::VectorXd IntegrationElement(
51 const Eigen::MatrixXd &local) const override;
52
54 [[nodiscard]] std::unique_ptr<Geometry> SubGeometry(dim_t codim,
55 dim_t i) const override;
56
63 [[nodiscard]] std::vector<std::unique_ptr<Geometry>> ChildGeometry(
64 const RefinementPattern &ref_pat, lf::base::dim_t codim) const override;
65
66 private:
70 Eigen::Matrix<double, Eigen::Dynamic, 8> coords_;
71
72 /*
73 * QuadO2 is parametrized by:
74 * alpha_ + beta_ * [x1, x2] + gamma_ * [x1^2, x2^2] + delta_ * [x1 * x2]
75 * + epsilon * [x1^2 * x2, x1 * x2^2]
76 */
77 Eigen::Matrix<double, Eigen::Dynamic, 1> alpha_;
78 Eigen::Matrix<double, Eigen::Dynamic, 2> beta_;
79 Eigen::Matrix<double, Eigen::Dynamic, 2> gamma_;
80 Eigen::Matrix<double, Eigen::Dynamic, 1> delta_;
81 Eigen::Matrix<double, Eigen::Dynamic, 2> epsilon_;
82
83 /* Coefficients for efficient evaluation of Jacobian() */
84 Eigen::Matrix<double, Eigen::Dynamic, 2> gamma_x_2_;
85 Eigen::Matrix<double, Eigen::Dynamic, 2> epsilon_x_2_;
86};
87} // namespace lf::geometry
88
89#endif /* QUAD_O2_H */
Represents a reference element with all its properties.
Definition ref_el.h:109
static constexpr RefEl kQuad()
Returns the reference quadrilateral.
Definition ref_el.h:169
Interface class for shape information on a mesh cell in the spirit of parametric finite element metho...
base::RefEl::dim_t dim_t
A second-order quadrilateral in the plane or in 3D space.
Definition quad_o2.h:29
Eigen::Matrix< double, Eigen::Dynamic, 2 > gamma_x_2_
Definition quad_o2.h:84
Eigen::MatrixXd Jacobian(const Eigen::MatrixXd &local) const override
Evaluate the jacobian of the mapping simultaneously at numPoints points.
Definition quad_o2.cc:69
Eigen::Matrix< double, Eigen::Dynamic, 2 > epsilon_x_2_
Definition quad_o2.h:85
base::RefEl RefEl() const override
The Reference element that defines the domain of this mapping.
Definition quad_o2.h:40
Eigen::Matrix< double, Eigen::Dynamic, 2 > beta_
Definition quad_o2.h:78
Eigen::VectorXd IntegrationElement(const Eigen::MatrixXd &local) const override
The integration element (factor appearing in integral transformation formula, see below) at number of...
Definition quad_o2.cc:134
Eigen::Matrix< double, Eigen::Dynamic, 2 > gamma_
Definition quad_o2.h:79
Eigen::MatrixXd Global(const Eigen::MatrixXd &local) const override
Map a number of points in local coordinates into the global coordinate system.
Definition quad_o2.cc:56
dim_t DimLocal() const override
Dimension of the domain of this mapping.
Definition quad_o2.h:38
Eigen::Matrix< double, Eigen::Dynamic, 8 > coords_
Coordinates of the 8 vertices/midpoints, stored in matrix columns.
Definition quad_o2.h:70
Eigen::MatrixXd JacobianInverseGramian(const Eigen::MatrixXd &local) const override
Evaluate the Jacobian * Inverse Gramian ( ) simultaneously at numPoints.
Definition quad_o2.cc:100
std::vector< std::unique_ptr< Geometry > > ChildGeometry(const RefinementPattern &ref_pat, lf::base::dim_t codim) const override
Generate geometry objects for child entities created in the course of refinement.
Definition quad_o2.cc:184
Eigen::Matrix< double, Eigen::Dynamic, 1 > alpha_
Definition quad_o2.h:77
Eigen::Matrix< double, Eigen::Dynamic, 1 > delta_
Definition quad_o2.h:80
Eigen::Matrix< double, Eigen::Dynamic, 2 > epsilon_
Definition quad_o2.h:81
QuadO2(Eigen::Matrix< double, Eigen::Dynamic, 8 > coords)
Constructor building quadrilateral from vertex/midpoint coordinates.
Definition quad_o2.cc:19
dim_t DimGlobal() const override
Dimension of the image of this mapping.
Definition quad_o2.h:39
std::unique_ptr< Geometry > SubGeometry(dim_t codim, dim_t i) const override
Construct a new Geometry() object that describes the geometry of the i-th sub-entity with codimension...
Definition quad_o2.cc:160
Abstract interface class for encoding topological local refinement
unsigned int dim_t
type for dimensions and co-dimensions and numbers derived from them
Definition types.h:32
Defines the Geometry interface and provides a number of classes that implement this interface + addit...
Definition compose.h:13