Defines a 4-node quadrilateral finite element using Q1 isoparametric interpolation. More...

#include <Quad4.h>

Inheritance diagram for Quad4:
FEShape

Public Member Functions

 Quad4 ()
 Default Constructor.
 
 Quad4 (const Element *element)
 Constructor when data of Element el are given.
 
 Quad4 (const Side *side)
 Constructor when data of Side sd are given.
 
 ~Quad4 ()
 Destructor.
 
void set (const Element *el)
 Choose element by giving its pointer.
 
void set (const Side *sd)
 Choose side by giving its pointer.
 
void setLocal (const Point< real_t > &s)
 Initialize local point coordinates in element.
 
void atGauss (int n, std::vector< real_t > &sh, std::vector< Point< real_t > > &dsh, std::vector< real_t > &w)
 Calculate shape functions and their partial derivatives and integration weights.
 
Point< real_t > Grad (const LocalVect< real_t, 4 > &u, const Point< real_t > &s)
 Return gradient of a function defined at element nodes.
 
real_t getMaxEdgeLength () const
 Return maximal edge length of quadrilateral.
 
real_t getMinEdgeLength () const
 Return minimal edge length of quadrilateral.
 
- Public Member Functions inherited from FEShape
 FEShape ()
 Default Constructor.
 
 FEShape (const Element *el)
 Constructor for an element.
 
 FEShape (const Side *sd)
 Constructor for a side.
 
virtual ~FEShape ()
 Destructor.
 
real_t Sh (size_t i) const
 Return shape function of node i at given point.
 
real_t Sh (size_t i, Point< real_t > s) const
 Calculate shape function of node i at a given point s.
 
real_t getDet () const
 Return determinant of jacobian.
 
Point< real_t > getCenter () const
 Return coordinates of center of element.
 
Point< real_t > getLocalPoint () const
 Localize a point in the element.
 
Point< real_t > getLocalPoint (const Point< real_t > &s) const
 Localize a point in the element.
 

Detailed Description

Defines a 4-node quadrilateral finite element using Q1 isoparametric interpolation.

The reference element is the square [-1,1]x[-1,1]. The user must take care to the fact that determinant of jacobian and other quantities depend on the point in the reference element where they are calculated. For this, before any utilization of shape functions or jacobian, function setLocal() must be invoked.

Author
Rachid Touzani

Constructor & Destructor Documentation

◆ Quad4() [1/2]

Quad4 ( const Element element)

Constructor when data of Element el are given.

Parameters
[in]elementPointer to Element

◆ Quad4() [2/2]

Quad4 ( const Side side)

Constructor when data of Side sd are given.

Parameters
[in]sidePointer to Side

Member Function Documentation

◆ atGauss()

void atGauss ( int  n,
std::vector< real_t > &  sh,
std::vector< Point< real_t > > &  dsh,
std::vector< real_t > &  w 
)

Calculate shape functions and their partial derivatives and integration weights.

Parameters
[in]nNumber of Gauss-Legendre integration points in each direction
[in]shVector of shape functions at Gauss points
[in]dshVector of shape function derivatives at Gauss points
[in]wWeights of integration formula at Gauss points

◆ Grad()

Point< real_t > Grad ( const LocalVect< real_t, 4 > &  u,
const Point< real_t > &  s 
)

Return gradient of a function defined at element nodes.

Parameters
[in]uVector of values at nodes
[in]sLocal coordinates (in [-1,1]*[-1,1]) of point where the gradient is evaluated
Returns
Value of gradient
Note
If the derivatives of shape functions were not computed before calling this function (by calling setLocal), this function will compute them

◆ setLocal()

void setLocal ( const Point< real_t > &  s)

Initialize local point coordinates in element.

Parameters
[in]sPoint in the reference element This function computes jacobian, shape functions and their partial derivatives at s. Other member functions only return these values.