PeriDEM 0.3.0
PeriDEM -- Peridynamics-based high-fidelity model for granular media
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions
fe::QuadElem Class Reference

A class for mapping and quadrature related operations for bi-linear quadrangle element. More...

#include <quadElem.h>

Inheritance diagram for fe::QuadElem:
Collaboration diagram for fe::QuadElem:

Public Member Functions

 QuadElem (size_t order)
 Constructor for quadrangle element.
 
double elemSize (const std::vector< util::Point > &nodes) override
 Returns the area of element.
 
- Public Member Functions inherited from fe::BaseElem
 BaseElem (size_t order, size_t element_type)
 Constructor.
 
size_t getElemType ()
 Get element type.
 
size_t getQuadOrder ()
 Get order of quadrature approximation.
 
size_t getNumQuadPoints ()
 Get number of quadrature points in the data.
 
virtual std::vector< double > getShapes (const util::Point &p, const std::vector< util::Point > &nodes)
 Returns the values of shape function at point p.
 
virtual std::vector< std::vector< double > > getDerShapes (const util::Point &p, const std::vector< util::Point > &nodes)
 Returns the values of derivative of shape function at point p.
 
std::vector< fe::QuadData > getQuadDatas (const std::vector< util::Point > &nodes)
 Get quadrature data mapped to the physical element (N, dN/dx, x, w, J). Shared isoparametric map; types only fill reference N, dN/dξ, and w.
 
std::vector< fe::QuadData > getQuadPoints (const std::vector< util::Point > &nodes)
 Same as getQuadDatas without transforming dN/dξ to dN/dx.
 

Private Member Functions

std::vector< double > getShapes (const util::Point &p) override
 Returns the values of shape function at point p on reference element.
 
std::vector< std::vector< double > > getDerShapes (const util::Point &p) override
 Returns the values of derivative of shape function at point p on reference element.
 
double getJacobian (const util::Point &p, const std::vector< util::Point > &nodes, std::vector< std::vector< double > > *J) override
 Computes the Jacobian of map \( \Phi: T^0 \to T \).
 
void init () override
 Compute the quadrature points for quadrangle element.
 

Additional Inherited Members

- Protected Member Functions inherited from fe::BaseElem
virtual util::Point mapPointToRefElem (const util::Point &p, const std::vector< util::Point > &nodes)
 Maps point p in a given element to the reference element and returns the mapped point.
 
- Protected Attributes inherited from fe::BaseElem
size_t d_quadOrder
 Order of quadrature point integration approximation.
 
size_t d_elemType
 Element type.
 
std::vector< fe::QuadData > d_quads
 Quadrature data collection.
 

Detailed Description

A class for mapping and quadrature related operations for bi-linear quadrangle element.

The reference quadrangle element \( T^0 \) is given by vertices \( (-1, -1), \, (1,-1), \, (1,1), \, (-1,1) \).

  1. The shape functions at point \( (\xi, \eta ) \in T^0 \) are

    \[N^0_1(\xi, \eta) = \frac{(1- \xi)(1 - \eta)}{4}, \quad N^0_2(\xi, \eta) = \frac{(1+ \xi)(1 - \eta)}{4}, \]

    \[ N^0_3(\xi, \eta) = \frac{(1+ \xi)(1 + \eta)}{4},\quad N^0_4(\xi, \eta) = \frac{(1- \xi)(1 + \eta)}{4}. \]

  2. Derivative of shape functions at point \( (\xi, \eta ) \in T^0 \) are as follows

    \[\frac{d N^0_1(\xi, \eta)}{d\xi} = \frac{-(1 - \eta)}{4}, \quad \frac{d N^0_1(\xi, \eta)}{d\eta} = \frac{-(1 - \xi)}{4}, \]

    \[\frac{d N^0_2(\xi, \eta)}{d\xi} = \frac{(1 - \eta)}{4}, \quad \frac{d N^0_2 (\xi, \eta)}{d\eta} = \frac{-(1 + \xi)}{4}, \]

    \[\frac{d N^0_3(\xi, \eta)}{d\xi} = \frac{(1 + \eta)}{4}, \quad \frac{d N^0_3 (\xi, \eta)}{d\eta} = \frac{(1 + \xi)}{4}, \]

    \[\frac{d N^0_4(\xi, \eta)}{d\xi} = \frac{-(1 + \eta)}{4}, \quad \frac{d N^0_4 (\xi, \eta)}{d\eta} = \frac{(1 - \xi)}{4}. \]

  3. Map \( \Phi: T^0 \to T\) is given by

    \[ x(\xi, \eta) = \sum_{i=1}^4 N^0_i(\xi, \eta) v^i_x, \quad y(\xi, \eta) = \sum_{i=1}^4 N^0_i(\xi, \eta) v^i_y \]

    where \( v^1, v^2, v^3, v^4\) are vertices of element \( T \).
  4. Jacobian of the map \( \Phi: T^0 \to T\) is given by

    \[ J = \left[ { \begin{array}{cc} \frac{dx}{d\xi} &\frac{dy}{d\xi} \\ \frac{dx}{d\eta} & \frac{dy}{d\eta} \\ \end{array} } \right] \]

    and determinant of Jacobian is

    \[ det(J) = \frac{dx}{d\xi} \times \frac{dy}{d\eta} - \frac{dy}{d\xi}\times \frac{dx}{d\eta}. \]

Definition at line 64 of file quadElem.h.

Constructor & Destructor Documentation

◆ QuadElem()

fe::QuadElem::QuadElem ( size_t  order)
explicit

Constructor for quadrangle element.

Parameters
orderOrder of quadrature point approximation

Definition at line 14 of file quadElem.cpp.

15 : fe::BaseElem(order, util::vtk_type_quad) {
16
17 // compute quad data
18 this->init();
19}
fe::BaseElem
A base class which provides methods to map points to/from reference element and to compute quadrature...
Definition baseElem.h:84
fe::QuadElem::init
void init() override
Compute the quadrature points for quadrangle element.
Definition quadElem.cpp:100
util::vtk_type_quad
static const int vtk_type_quad
Integer flag for quad element.
Definition feElementDefs.h:67

References init().

Here is the call graph for this function:

Member Function Documentation

◆ elemSize()

double fe::QuadElem::elemSize ( const std::vector< util::Point > &  nodes)
overridevirtual

Returns the area of element.

If quadrangle \( T \) is given by points \( v^1, v^2, v^3, v^4 \) then the area is

\[ area(T) = \frac{(-v^1_1 + v^2_1 + v^3_1 - v^4_1) (-v^1_2 - v^2_2 + v^3_2 + v^4_2) - (-v^1_1 - v^2_1 + v^3_1 + v^4_1) (-v^1_2 + v^2_2 + v^3_2 - v^4_2)}{4}, \]

where \( v^i_1, v^i_2 \) are the x and y component of point \( v^i \).

Note that area and Jacobian of map \( \Phi: T^0 \to T \) are related as

\[ area(T) = area(T^0) \times det(J(\xi = 0, \eta = 0)), \]

where \( area(T^0) = 4 \) and \( J(\xi = 0, \eta = 0) \) is the Jacobian of map at point \( (0,0) \in T^0 \).

Parameters
nodesVertices of element
Returns
vector Vector of shape functions at point p

Implements fe::BaseElem.

Definition at line 21 of file quadElem.cpp.

21 {
22 return 0.25 *
23 ((-nodes[0].d_x + nodes[1].d_x + nodes[2].d_x - nodes[3].d_x) *
24 (-nodes[0].d_y - nodes[1].d_y + nodes[2].d_y + nodes[3].d_y) -
25 (-nodes[0].d_x - nodes[1].d_x + nodes[2].d_x + nodes[3].d_x) *
26 (-nodes[0].d_y + nodes[1].d_y + nodes[2].d_y - nodes[3].d_y));
27}

◆ getDerShapes()

std::vector< std::vector< double > > fe::QuadElem::getDerShapes ( const util::Point &  p)
overrideprivatevirtual

Returns the values of derivative of shape function at point p on reference element.

Parameters
pLocation of point
Returns
vector Vector of derivative of shape functions

Implements fe::BaseElem.

Definition at line 43 of file quadElem.cpp.

43 {
44
45 // N1 = (1 - xi)(1 - eta)/4
46 // --> d N1/d xi = -(1 - eta)/4, d N1/d eta = -(1 - xi)/4
47 //
48 // N2 = (1 + xi)(1 - eta)/4
49 // --> d N2/d xi = (1 - eta)/4, d N2/d eta = -(1 + xi)/4
50 //
51 // N3 = (1 + xi)(1 + eta)/4
52 // --> d N3/d xi = (1 + eta)/4, d N3/d eta = (1 + xi)/4
53 //
54 // N4 = (1 - xi)(1 + eta)/4
55 // --> d N4/d xi = -(1 + eta)/4, d N4/d eta = (1 - xi)/4
56 std::vector<std::vector<double>> r;
57 r.push_back(std::vector<double>{-0.25 * (1. - p.d_y), -0.25 * (1. - p.d_x)});
58 r.push_back(std::vector<double>{0.25 * (1. - p.d_y), -0.25 * (1. + p.d_x)});
59 r.push_back(std::vector<double>{0.25 * (1. + p.d_y), 0.25 * (1. + p.d_x)});
60 r.push_back(std::vector<double>{-0.25 * (1. + p.d_y), 0.25 * (1. - p.d_x)});
61
62 return r;
63}
util::Point::d_y
double d_y
the y coordinate
Definition point.h:36
util::Point::d_x
double d_x
the x coordinate
Definition point.h:33

References util::Point::d_x, and util::Point::d_y.

◆ getJacobian()

double fe::QuadElem::getJacobian ( const util::Point &  p,
const std::vector< util::Point > &  nodes,
std::vector< std::vector< double > > *  J 
)
overrideprivatevirtual

Computes the Jacobian of map \( \Phi: T^0 \to T \).

Parameters
pLocation of point in reference element
nodesVertices of element
JMatrix to store the Jacobian (if not nullptr)
Returns
det(J) Determinant of the Jacobain

Implements fe::BaseElem.

Definition at line 65 of file quadElem.cpp.

67 {
68
69 auto der_shapes = getDerShapes(p);
70 if (J != nullptr) {
71 J->resize(2);
72 (*J)[0] = std::vector<double>{
73 der_shapes[0][0] * nodes[0].d_x + der_shapes[1][0] * nodes[1].d_x +
74 der_shapes[2][0] * nodes[2].d_x + der_shapes[3][0] * nodes[3].d_x,
75 der_shapes[0][0] * nodes[0].d_y + der_shapes[1][0] * nodes[1].d_y +
76 der_shapes[2][0] * nodes[2].d_y + der_shapes[3][0] * nodes[3].d_y};
77 (*J)[1] = std::vector<double>{
78 der_shapes[0][1] * nodes[0].d_x + der_shapes[1][1] * nodes[1].d_x +
79 der_shapes[2][1] * nodes[2].d_x + der_shapes[3][1] * nodes[3].d_x,
80 der_shapes[0][1] * nodes[0].d_y + der_shapes[1][1] * nodes[1].d_y +
81 der_shapes[2][1] * nodes[2].d_y + der_shapes[3][1] * nodes[3].d_y};
82
83 return (*J)[0][0] * (*J)[1][1] - (*J)[0][1] * (*J)[1][0];
84 }
85
86 return (der_shapes[0][0] * nodes[0].d_x + der_shapes[1][0] * nodes[1].d_x +
87 der_shapes[2][0] * nodes[2].d_x + der_shapes[3][0] * nodes[3].d_x) *
88 (der_shapes[0][1] * nodes[0].d_y +
89 der_shapes[1][1] * nodes[1].d_y +
90 der_shapes[2][1] * nodes[2].d_y +
91 der_shapes[3][1] * nodes[3].d_y) -
92 (der_shapes[0][0] * nodes[0].d_y + der_shapes[1][0] * nodes[1].d_y +
93 der_shapes[2][0] * nodes[2].d_y + der_shapes[3][0] * nodes[3].d_y) *
94 (der_shapes[0][1] * nodes[0].d_x +
95 der_shapes[1][1] * nodes[1].d_x +
96 der_shapes[2][1] * nodes[2].d_x +
97 der_shapes[3][1] * nodes[3].d_x);
98}
fe::QuadElem::getDerShapes
std::vector< std::vector< double > > getDerShapes(const util::Point &p) override
Returns the values of derivative of shape function at point p on reference element.
Definition quadElem.cpp:43

◆ getShapes()

std::vector< double > fe::QuadElem::getShapes ( const util::Point &  p)
overrideprivatevirtual

Returns the values of shape function at point p on reference element.

Parameters
pLocation of point
Returns
vector Vector of shape functions at point p

Implements fe::BaseElem.

Definition at line 31 of file quadElem.cpp.

31 {
32
33 // N1 = (1 - xi)(1 - eta)/4
34 // N2 = (1 + xi)(1 - eta)/4
35 // N3 = (1 + xi)(1 + eta)/4
36 // N4 = (1 - xi)(1 + eta)/4
37 return std::vector<double>{
38 0.25 * (1. - p.d_x) * (1. - p.d_y), 0.25 * (1. + p.d_x) * (1. - p.d_y),
39 0.25 * (1. + p.d_x) * (1. + p.d_y), 0.25 * (1. - p.d_x) * (1. + p.d_y)};
40}

References util::Point::d_x, and util::Point::d_y.

◆ init()

void fe::QuadElem::init ( )
overrideprivatevirtual

Compute the quadrature points for quadrangle element.

Implements fe::BaseElem.

Definition at line 100 of file quadElem.cpp.

100 {
101
102 //
103 // compute quad data for reference quadrangle with vertex at
104 // p1 = (-1,-1), p2 = (1,-1), p3 = (1,1), p4 = (-1,1)
105 //
106 // Shape functions are
107 // N1 = (1 - xi)(1 - eta)/4
108 // N2 = (1 + xi)(1 - eta)/4
109 // N3 = (1 + xi)(1 + eta)/4
110 // N4 = (1 - xi)(1 + eta)/4
111 //
112 //
113 // Let [-1,1] is the 1-d reference element and {x1, x2, x3,.., xN} are N
114 // quad points for 1-d domain and {w1, w2, w3,..., wN} are respective
115 // weights. Then, the Nth order quad points in Quadrangle [-1,1]x[-1,1] is
116 // given by N^2 points and
117 //
118 // (i,j) point is (xi, xj) and weight is wi \times wj
119 //
120
121 if (!d_quads.empty())
122 return;
123
124 // no point in zeroth order
125 if (d_quadOrder == 0)
126 d_quads.resize(0);
127
128 // 2x2 identity matrix
129 std::vector<std::vector<double>> ident_mat;
130 ident_mat.push_back(std::vector<double>{1., 0.});
131 ident_mat.push_back(std::vector<double>{0., 1.});
132
133 //
134 // first order quad points
135 //
136 if (d_quadOrder == 1) {
137 d_quads.clear();
138 // 1-d points are: {0} and weights are: {2}
139 int npts = 1;
140 std::vector<double> x = std::vector<double>(1, 0.);
141 std::vector<double> w = std::vector<double>(1, 2.);
142 for (size_t i = 0; i < npts; i++)
143 for (size_t j = 0; j < npts; j++) {
144
145 fe::QuadData qd;
146 qd.d_w = w[i] * w[j];
147 qd.d_p = util::Point(x[i], x[j], 0.);
148 qd.d_shapes = getShapes(qd.d_p);
149 qd.d_derShapes = getDerShapes(qd.d_p);
150 qd.d_J = ident_mat;
151 qd.d_detJ = 1.;
152 d_quads.push_back(qd);
153 }
154 }
155
156 //
157 // second order quad points
158 //
159 if (d_quadOrder == 2) {
160 d_quads.clear();
161 // 1-d points are: {-1/sqrt{3], 1/sqrt{3}} and weights are: {1,1}
162 int npts = 2;
163 std::vector<double> x =
164 std::vector<double>{-1. / std::sqrt(3.), 1. / std::sqrt(3.)};
165 std::vector<double> w = std::vector<double>{1., 1.};
166 for (size_t i = 0; i < npts; i++)
167 for (size_t j = 0; j < npts; j++) {
168
169 fe::QuadData qd;
170 qd.d_w = w[i] * w[j];
171 qd.d_p = util::Point(x[i], x[j], 0.);
172 qd.d_shapes = getShapes(qd.d_p);
173 qd.d_derShapes = getDerShapes(qd.d_p);
174 qd.d_J = ident_mat;
175 qd.d_detJ = 1.;
176 d_quads.push_back(qd);
177 }
178 }
179
180 //
181 // third order quad points
182 //
183 if (d_quadOrder == 3) {
184 d_quads.clear();
185 // 1-d points are: {-sqrt{3}/sqrt{5}, 0, sqrt{3}/sqrt{5}}
186 // weights are: {5/9, 8/9, 5/9}
187 int npts = 3;
188 std::vector<double> x = std::vector<double>{
189 -std::sqrt(3.) / std::sqrt(5.), 0., std::sqrt(3.) / std::sqrt(5.)};
190 std::vector<double> w = std::vector<double>{5. / 9., 8. / 9., 5. / 9.};
191 for (size_t i = 0; i < npts; i++)
192 for (size_t j = 0; j < npts; j++) {
193
194 fe::QuadData qd;
195 qd.d_w = w[i] * w[j];
196 qd.d_p = util::Point(x[i], x[j], 0.);
197 qd.d_shapes = getShapes(qd.d_p);
198 qd.d_derShapes = getDerShapes(qd.d_p);
199 qd.d_J = ident_mat;
200 qd.d_detJ = 1.;
201 d_quads.push_back(qd);
202 }
203 }
204
205 //
206 // fourth order quad points
207 //
208 if (d_quadOrder == 4) {
209 d_quads.clear();
210 int npts = 4;
211 std::vector<double> x =
212 std::vector<double>{-0.3399810435848563, 0.3399810435848563,
213 -0.8611363115940526, 0.8611363115940526};
214 std::vector<double> w =
215 std::vector<double>{0.6521451548625461, 0.6521451548625461,
216 0.3478548451374538, 0.3478548451374538};
217 for (size_t i = 0; i < npts; i++)
218 for (size_t j = 0; j < npts; j++) {
219
220 fe::QuadData qd;
221 qd.d_w = w[i] * w[j];
222 qd.d_p = util::Point(x[i], x[j], 0.);
223 qd.d_shapes = getShapes(qd.d_p);
224 qd.d_derShapes = getDerShapes(qd.d_p);
225 qd.d_J = ident_mat;
226 qd.d_detJ = 1.;
227 d_quads.push_back(qd);
228 }
229 }
230
231 //
232 // fifth order quad points
233 //
234 if (d_quadOrder == 5) {
235 d_quads.clear();
236 int npts = 5;
237 std::vector<double> x =
238 std::vector<double>{0., -0.5384693101056831, 0.5384693101056831,
239 -0.9061798459386640, 0.9061798459386640};
240 std::vector<double> w = std::vector<double>{
241 0.5688888888888889, 0.4786286704993665, 0.4786286704993665,
242 0.2369268850561891, 0.2369268850561891};
243 for (size_t i = 0; i < npts; i++)
244 for (size_t j = 0; j < npts; j++) {
245
246 fe::QuadData qd;
247 qd.d_w = w[i] * w[j];
248 qd.d_p = util::Point(x[i], x[j], 0.);
249 qd.d_shapes = getShapes(qd.d_p);
250 qd.d_derShapes = getDerShapes(qd.d_p);
251 qd.d_J = ident_mat;
252 qd.d_detJ = 1.;
253 d_quads.push_back(qd);
254 }
255 }
256}
fe::BaseElem::d_quads
std::vector< fe::QuadData > d_quads
Quadrature data collection.
Definition baseElem.h:224
fe::BaseElem::d_quadOrder
size_t d_quadOrder
Order of quadrature point integration approximation.
Definition baseElem.h:218
fe::QuadElem::getShapes
std::vector< double > getShapes(const util::Point &p) override
Returns the values of shape function at point p on reference element.
Definition quadElem.cpp:31
fe::QuadData
A struct to store the quadrature data. List of data are.
Definition quadData.h:23
fe::QuadData::d_shapes
std::vector< double > d_shapes
Value of shape functions at quad point p.
Definition quadData.h:37
fe::QuadData::d_derShapes
std::vector< std::vector< double > > d_derShapes
Derivatives of shape functions at the quad point.
Definition quadData.h:46
fe::QuadData::d_w
double d_w
Quadrature weight.
Definition quadData.h:26
fe::QuadData::d_p
util::Point d_p
Quadrature point in 1-d, 2-d or 3-d.
Definition quadData.h:29
fe::QuadData::d_J
std::vector< std::vector< double > > d_J
Jacobian of the map from reference element to the element.
Definition quadData.h:54
fe::QuadData::d_detJ
double d_detJ
Determinant of the Jacobian of the map from reference element to the element.
Definition quadData.h:60
util::Point
A structure to represent 3d vectors.
Definition point.h:30

References fe::QuadData::d_derShapes, fe::QuadData::d_detJ, fe::QuadData::d_J, fe::QuadData::d_p, fe::QuadData::d_shapes, and fe::QuadData::d_w.

Referenced by QuadElem().

Here is the caller graph for this function:
The documentation for this class was generated from the following files: