Element on an edge representing the normal derivatives of neighbouring elements, especially their mean or jump. More...

#include <neumannTraceElement.hh>

Inheritance diagram for hp2D::NeumannTraceElement< F >:
hp1D::BaseElement< Real > concepts::ElementWithCell< Real > hp1D::IntegrableElm concepts::IntegrationCell

Classes

struct  mapPart
 
class  NTShapeFunction
 

Public Types

typedef Real FieldT
 
enum  intFormType { ZERO, ONE, TWO, THREE }
 Integration form, which determines terms coming from integration over reference element. More...
 
typedef Real type
 
typedef concepts::ElementAndFacette< hp2D::Element< F > > UnderlyingElement
 

Public Member Functions

void addElement (const hp2D::Quad< Real > &quad, uint k, Real weight=1.0)
 Adds the contribution to the Neumann trace from the element on one of the (at most two) sides. More...
 
void addIrrElement (const hp2D::Quad< Real > &coarseQuad, uint k_coarse, Real weight)
 Basically the same as addElement, but with the difference in the meaning. More...
 
void buildT ()
 After adding the Elements this routine builds the TMatrix and polynomial informations on the element. More...
 
virtual const concepts::EdgeNdcell () const
 Returns the cell on which the element is built. More...
 
const concepts::Real3d chi (const Real x) const
 Computes the element map. More...
 
void computeShpfctVals (const concepts::Array< Real > &xP, concepts::Array< Real > &val) const
 Evaluation method to get basisfunctions evaluated on requested points this should be used for functions, i.e. More...
 
Real3d elemMap (const Real coord_local) const
 
Real3d elemMap (const Real2d &coord_local) const
 
Real3d elemMap (const Real3d &coord_local) const
 
virtual const concepts::ElementGraphics< Real > * graphics () const
 Returns element graphics class. More...
 
const concepts::QuadratureRule1dintegration () const
 Returns the integration rule. More...
 
Real jacobianDeterminant (const Real x) const
 Computes the determinant of the Jacobian. More...
 
 NeumannTraceElement (const concepts::EdgeNd &cell, uint p, const mapPart *coarseCell=0)
 Constructor. More...
 
ushort p () const
 
virtual bool quadraturePoint (uint i, intPoint &p, intFormType form=ZERO, bool localCoord=false) const
 Delivers a quadrature point. More...
 
virtual void recomputeShapefunctions ()
 Recompute shape functions, e.g. More...
 
virtual const concepts::ShapeFunction1D< Real > * shpfct () const
 Returns the shape functions. More...
 
virtual const concepts::ShapeFunction1D< Real > * shpfct () const=0
 Returns the shape functions. More...
 
virtual const concepts::Edgesupport () const
 
virtual const concepts::TMatrix< F > & T () const
 
virtual const TMatrixBase< Real > & T () const=0
 Returns the T matrix of the element. More...
 
const concepts::Sequence< UnderlyingElement > & uelm () const
 Returns the Underlying Element(s) More...
 
virtual concepts::Real3d vertex (uint i) const
 
virtual ~NeumannTraceElement ()
 

Static Public Member Functions

static concepts::QuadRuleFactoryrule ()
 Access to the quadrature rule, which is valid for all elements of this type (hp1D::IntegrableElm). More...
 

Protected Member Functions

virtual std::ostream & info (std::ostream &os) const
 

Protected Attributes

const concepts::EdgeNdcell_
 The cell. More...
 
std::unique_ptr< concepts::QuadratureRule1dint_
 The integration rule. More...
 
ushort p_
 T matrix of the element. More...
 

Static Protected Attributes

static std::unique_ptr< LineGraphics > graphics_
 
static concepts::QuadRuleFactory rule_
 

Private Member Functions

void compute_ (const concepts::Array< Real > &q, concepts::Array< Real > &val) const
 Routine calculating the normal derivative of the basis functions (of underlying quads), with normalvector pointing outwards the respective underlying Quad on a given Set of Points. More...
 
void irregularCompute_ (const concepts::Array< Real > &q, concepts::Array< Real > &val) const
 Routine calculating the normal derivative of the basis functions (of underlying quads), with normalvector pointing outwards the respective underlying Quad on a given Set of Points. More...
 

Private Attributes

mapPart coarseCell_
 
bool irregular_
 
std::unique_ptr< concepts::ShapeFunction1D< Real > > shpfct_
 
concepts::TMatrix< Real > T_
 
concepts::Sequence< UnderlyingElementuelm_
 
concepts::Array< Real > values_
 The actual storage of the values of the shape functions along the edge on integration points. More...
 
concepts::Sequence< Real > weights_
 Weight for each side. This defines the TraceType on the Element. More...
 

Detailed Description

template<class F = Real>
class hp2D::NeumannTraceElement< F >

Element on an edge representing the normal derivatives of neighbouring elements, especially their mean or jump.

Also a PLUS and MINUS setting is involded, that is a FIRST trace kind, with respect to the K+ and K- introduced elements via the setted Normalvector.

The number of shape functions are the sum of the number of shape functions of the two neighbouring 2D elements.

Author
Robert Gruhlke 2013, 2016

Definition at line 42 of file neumannTraceElement.hh.

Member Typedef Documentation

◆ FieldT

typedef Real hp1D::BaseElement< Real >::FieldT
inherited

Definition at line 82 of file element.hh.

◆ type

typedef Real concepts::ElementWithCell< Real >::type
inherited

Definition at line 81 of file element.hh.

◆ UnderlyingElement

Definition at line 44 of file neumannTraceElement.hh.

Member Enumeration Documentation

◆ intFormType

Integration form, which determines terms coming from integration over reference element.

Enumerator
ZERO 
ONE 
TWO 
THREE 

Definition at line 29 of file integral.hh.

Constructor & Destructor Documentation

◆ NeumannTraceElement()

template<class F = Real>
hp2D::NeumannTraceElement< F >::NeumannTraceElement ( const concepts::EdgeNd cell,
uint  p,
const mapPart coarseCell = 0 
)

Constructor.

For irregular meshes the underlying elements of a neumanntraceelement is a fine Quad and a neighboured coarse Quad. The coarseCell describes the edgemapping of the coarseedge of which the cell is a subset (part). The coarseCellpart also marks if a neumanntraceElement is build upon regular or irregular mesh

Parameters
cellgeometric cell of this element
phighest Polynomial degree of underlying Elements
coarseCellfor irregular meshes where cell's underlying

◆ ~NeumannTraceElement()

template<class F = Real>
virtual hp2D::NeumannTraceElement< F >::~NeumannTraceElement ( )
virtual

Member Function Documentation

◆ addElement()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::addElement ( const hp2D::Quad< Real > &  quad,
uint  k,
Real  weight = 1.0 
)

Adds the contribution to the Neumann trace from the element on one of the (at most two) sides.

Basis functions between quads do not have to be continous. In this case the highest polynomial degree along edge is set, i.e. for integration rule

Parameters
quadQuad with k-th edge is NeumannTraceElement's topological cell.

◆ addIrrElement()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::addIrrElement ( const hp2D::Quad< Real > &  coarseQuad,
uint  k_coarse,
Real  weight 
)

Basically the same as addElement, but with the difference in the meaning.

In this case the general ancestor edge of this finer Edge is the k-th edge of the coarseQuad.

Parameters
coarseQuadQuad with k-th edge is general ancestor of NeumannTraceElement's topological cell.

◆ buildT()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::buildT ( )

After adding the Elements this routine builds the TMatrix and polynomial informations on the element.

◆ cell()

virtual const concepts::EdgeNd& hp1D::BaseElement< Real >::cell
inlinevirtualinherited

Returns the cell on which the element is built.

Implements concepts::ElementWithCell< Real >.

Definition at line 99 of file element.hh.

◆ chi()

const concepts::Real3d hp1D::IntegrableElm::chi ( const Real  x) const
inlineinherited

Computes the element map.

The reference element is [0,1].

Definition at line 40 of file element.hh.

◆ compute_()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::compute_ ( const concepts::Array< Real > &  q,
concepts::Array< Real > &  val 
) const
private

Routine calculating the normal derivative of the basis functions (of underlying quads), with normalvector pointing outwards the respective underlying Quad on a given Set of Points.

Parameters
qevaluation domain points
valof evaluations of shapefunctions on q values

◆ computeShpfctVals()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::computeShpfctVals ( const concepts::Array< Real > &  xP,
concepts::Array< Real > &  val 
) const

Evaluation method to get basisfunctions evaluated on requested points this should be used for functions, i.e.

hp2D::Value

Parameters
xPrequested Points to evaluate on
valevaluated shpfct, method resizes array.

◆ elemMap() [1/3]

Real3d concepts::ElementWithCell< Real >::elemMap ( const Real  coord_local) const
inlineinherited

Definition at line 86 of file element.hh.

◆ elemMap() [2/3]

Real3d concepts::ElementWithCell< Real >::elemMap ( const Real2d coord_local) const
inlineinherited

Definition at line 90 of file element.hh.

◆ elemMap() [3/3]

Real3d concepts::ElementWithCell< Real >::elemMap ( const Real3d coord_local) const
inlineinherited

Definition at line 94 of file element.hh.

◆ graphics()

virtual const concepts::ElementGraphics<Real>* hp1D::BaseElement< Real >::graphics
virtualinherited

Returns element graphics class.

◆ info()

template<class F = Real>
virtual std::ostream& hp2D::NeumannTraceElement< F >::info ( std::ostream &  os) const
protectedvirtual

Reimplemented from hp1D::BaseElement< Real >.

◆ integration()

const concepts::QuadratureRule1d* hp1D::IntegrableElm::integration ( ) const
inlineinherited

Returns the integration rule.

Definition at line 51 of file element.hh.

◆ irregularCompute_()

template<class F = Real>
void hp2D::NeumannTraceElement< F >::irregularCompute_ ( const concepts::Array< Real > &  q,
concepts::Array< Real > &  val 
) const
private

Routine calculating the normal derivative of the basis functions (of underlying quads), with normalvector pointing outwards the respective underlying Quad on a given Set of Points.

This is for irregular underlying elements

Parameters
qevaluation domain points
valof evaluations of shapefunctions on q values

◆ jacobianDeterminant()

Real hp1D::IntegrableElm::jacobianDeterminant ( const Real  x) const
inlineinherited

Computes the determinant of the Jacobian.

Definition at line 45 of file element.hh.

◆ p()

ushort hp1D::BaseElement< Real >::p
inlineinherited

Definition at line 103 of file element.hh.

◆ quadraturePoint()

virtual bool hp1D::IntegrableElm::quadraturePoint ( uint  i,
intPoint p,
intFormType  form = ZERO,
bool  localCoord = false 
) const
virtualinherited

Delivers a quadrature point.

Quadrature point consists of coordinates (for evaluation of formulas) and intermediate data, consisting of the weight and term coming from mapping.

Returns false, if the number of quadrature points is overstepped.

Parameters
inumber of quadrature point
intPointdata given back
formIntegration form
localCoordIf true, local coordinates are returned. Else physical coordinates.

Implements concepts::IntegrationCell.

◆ recomputeShapefunctions()

template<class F = Real>
virtual void hp2D::NeumannTraceElement< F >::recomputeShapefunctions ( )
virtual

Recompute shape functions, e.g.

for other abscissas redefined through IntegrableElm::rule().set(...)

Implements hp1D::BaseElement< Real >.

◆ rule()

static concepts::QuadRuleFactory& hp1D::IntegrableElm::rule ( )
inlinestaticinherited

Access to the quadrature rule, which is valid for all elements of this type (hp1D::IntegrableElm).

Change of the quadrature rule is put into practice for newly created elements and for already created elements by precomputing the integration points and shape functions on them.

Definition at line 62 of file element.hh.

◆ shpfct() [1/2]

template<class F = Real>
virtual const concepts::ShapeFunction1D<Real>* hp2D::NeumannTraceElement< F >::shpfct ( ) const
virtual

Returns the shape functions.

◆ shpfct() [2/2]

virtual const concepts::ShapeFunction1D<Real>* hp1D::BaseElement< Real >::shpfct
pure virtualinherited

Returns the shape functions.

◆ support()

virtual const concepts::Edge& hp1D::BaseElement< Real >::support
inlinevirtualinherited

Definition at line 91 of file element.hh.

◆ T() [1/2]

template<class F = Real>
virtual const concepts::TMatrix<F>& hp2D::NeumannTraceElement< F >::T ( ) const
inlinevirtual

Definition at line 121 of file neumannTraceElement.hh.

◆ T() [2/2]

virtual const TMatrixBase<Real >& concepts::ElementWithCell< Real >::T
pure virtualinherited

Returns the T matrix of the element.

Implemented in hp3D::Element< Real >, and hp2D::Element< Real >.

◆ uelm()

template<class F = Real>
const concepts::Sequence<UnderlyingElement>& hp2D::NeumannTraceElement< F >::uelm ( ) const
inline

Returns the Underlying Element(s)

Definition at line 151 of file neumannTraceElement.hh.

◆ vertex()

virtual concepts::Real3d hp1D::BaseElement< Real >::vertex ( uint  i) const
inlinevirtualinherited

Definition at line 95 of file element.hh.

Member Data Documentation

◆ cell_

const concepts::EdgeNd& hp1D::IntegrableElm::cell_
protectedinherited

The cell.

Definition at line 70 of file element.hh.

◆ coarseCell_

template<class F = Real>
mapPart hp2D::NeumannTraceElement< F >::coarseCell_
private

Definition at line 190 of file neumannTraceElement.hh.

◆ graphics_

std::unique_ptr<LineGraphics> hp1D::BaseElement< Real >::graphics_
staticprotectedinherited

Definition at line 126 of file element.hh.

◆ int_

std::unique_ptr<concepts::QuadratureRule1d> hp1D::IntegrableElm::int_
protectedinherited

The integration rule.

Definition at line 72 of file element.hh.

◆ irregular_

template<class F = Real>
bool hp2D::NeumannTraceElement< F >::irregular_
private

Definition at line 188 of file neumannTraceElement.hh.

◆ p_

ushort hp1D::BaseElement< Real >::p_
protectedinherited

T matrix of the element.

Definition at line 124 of file element.hh.

◆ rule_

concepts::QuadRuleFactory hp1D::IntegrableElm::rule_
staticprotectedinherited

Definition at line 74 of file element.hh.

◆ shpfct_

template<class F = Real>
std::unique_ptr<concepts::ShapeFunction1D<Real> > hp2D::NeumannTraceElement< F >::shpfct_
mutableprivate

Definition at line 179 of file neumannTraceElement.hh.

◆ T_

template<class F = Real>
concepts::TMatrix<Real> hp2D::NeumannTraceElement< F >::T_
private

Definition at line 182 of file neumannTraceElement.hh.

◆ uelm_

template<class F = Real>
concepts::Sequence<UnderlyingElement> hp2D::NeumannTraceElement< F >::uelm_
private

Definition at line 192 of file neumannTraceElement.hh.

◆ values_

template<class F = Real>
concepts::Array<Real> hp2D::NeumannTraceElement< F >::values_
mutableprivate

The actual storage of the values of the shape functions along the edge on integration points.

The values corresponding to the local shapefunctions of its underlying quad(s) are computed with the normalvector pointing outwards the quad respectivly.

This have to be noticed for defining the weights for trace types. mutable to work with const compute_ as needed for const

Definition at line 201 of file neumannTraceElement.hh.

◆ weights_

template<class F = Real>
concepts::Sequence<Real> hp2D::NeumannTraceElement< F >::weights_
private

Weight for each side. This defines the TraceType on the Element.

Definition at line 204 of file neumannTraceElement.hh.


The documentation for this class was generated from the following file:
Page URL: http://wiki.math.ethz.ch/bin/view/Concepts/WebHome
21 August 2020
© 2020 Eidgenössische Technische Hochschule Zürich