|
IFEM
90A354
|
Abstract base class representing a system level integrated quantity. More...
#include <Integrand.h>

Public Types | |
| enum | Traits { STANDARD = 0 , NO_DERIVATIVES = 1 , SECOND_DERIVATIVES = 1<< 1 , THIRD_DERIVATIVES = 1<< 2 , AVERAGE = 1<< 3 , ELEMENT_CORNERS = 1<< 4 , ELEMENT_CENTER = 1<< 5 , G_MATRIX = 1<< 6 , NODAL_ROTATIONS = 1<< 7 , XO_ELEMENTS = 1<< 8 , INTERFACE_TERMS = 1<< 9 , NORMAL_DERIVS = 1<<10 , UPDATED_NODES = 1<<11 , PIOLA_MAPPING = 1<<12 , POINT_DEFORMATION = 1<<13 } |
| Enum defining the additional terms that an Integrand may require. More... | |
Public Member Functions | |
| virtual | ~Integrand () |
| Empty destructor. | |
| virtual void | setNeumannOrder (char) |
| Defines the Neumann order that is the subject of integration. More... | |
| virtual void | initPatch (size_t) |
| Define the index of the patch being processed. | |
| virtual LocalIntegral * | getLocalIntegral (size_t nen, size_t iEl, bool neumann=false) const =0 |
| Returns a local integral contribution object for the given element. More... | |
| virtual LocalIntegral * | getLocalIntegral (const std::vector< size_t > &nen, size_t iEl, bool neumann=false) const |
| Returns a local integral contribution object for the given element. More... | |
| virtual bool | initElement (const std::vector< int > &MNPC, const FiniteElement &fe, const Vec3 &X0, size_t nPt, LocalIntegral &elmInt)=0 |
| Initializes current element for numerical integration. More... | |
| virtual bool | initElement (const std::vector< int > &MNPC, LocalIntegral &elmInt)=0 |
| Initializes current element for numerical integration. More... | |
| virtual bool | initElement (const std::vector< int > &MNPC, const std::vector< size_t > &elem_sizes, const std::vector< size_t > &basis_sizes, LocalIntegral &elmInt)=0 |
| Initializes current element for numerical integration (mixed). More... | |
| virtual bool | initElement (const std::vector< int > &MNPC, const MxFiniteElement &fe, const std::vector< size_t > &elem_sizes, const std::vector< size_t > &basis_sizes, LocalIntegral &elmInt)=0 |
| Initializes current element for numerical integration (mixed). More... | |
| virtual bool | initElementBou (const std::vector< int > &MNPC, LocalIntegral &elmInt)=0 |
| Initializes current element for boundary integration. More... | |
| virtual bool | initElementBou (const std::vector< int > &MNPC, const std::vector< size_t > &elem_sizes, const std::vector< size_t > &basis_sizes, LocalIntegral &elmInt)=0 |
| Initializes current element for boundary integration (mixed). More... | |
| virtual SIM::SolutionMode | getMode (bool=false) const =0 |
| Returns current solution mode. | |
| virtual int | getIntegrandType () const |
| Defines which FE quantities are needed by the integrand. | |
| virtual int | getReducedIntegration (int) const |
| Returns the number of reduced-order integration points. | |
| virtual int | getBouIntegrationPoints (int nGP) const |
| Returns the number of boundary integration points. | |
| virtual bool | reducedInt (LocalIntegral &elmInt, const FiniteElement &fe, const Vec3 &X) const |
| Evaluates reduced integration terms at an interior point. More... | |
| virtual bool | evalInt (LocalIntegral &elmInt, const FiniteElement &fe, const TimeDomain &time, const Vec3 &X) const |
| Evaluates the integrand at an interior point. More... | |
| virtual bool | evalIntMx (LocalIntegral &elmInt, const MxFiniteElement &fe, const TimeDomain &time, const Vec3 &X) const |
| Evaluates the integrand at an interior point. More... | |
| virtual bool | evalInt (LocalIntegral &elmInt, const FiniteElement &fe, const TimeDomain &time, const Vec3 &X, const Vec3 &normal) const |
| Evaluates the integrand at an element interface point. More... | |
| virtual bool | evalIntMx (LocalIntegral &elmInt, const MxFiniteElement &fe, const TimeDomain &time, const Vec3 &X, const Vec3 &normal) const |
| Evaluates the integrand at an element interface point. More... | |
| virtual bool | evalPoint (LocalIntegral &elmInt, const FiniteElement &fe, const Vec3 &pval) |
| Evaluates the dirac-delta integrand at a specified point. More... | |
| virtual bool | finalizeElement (LocalIntegral &elmInt, const FiniteElement &fe, const TimeDomain &time, size_t iGP=0) |
| Finalizes the element quantities after the numerical integration. More... | |
| virtual bool | finalizeElementBou (LocalIntegral &elmInt, const FiniteElement &fe, const TimeDomain &time) |
| Finalizes the element quantities after boundary integration. More... | |
| virtual bool | evalBou (LocalIntegral &elmInt, const FiniteElement &fe, const TimeDomain &time, const Vec3 &X, const Vec3 &normal) const |
| Evaluates the integrand at a boundary point. More... | |
| virtual bool | evalBouMx (LocalIntegral &elmInt, const MxFiniteElement &fe, const TimeDomain &time, const Vec3 &X, const Vec3 &normal) const |
| Evaluates the integrand at a boundary point. More... | |
| virtual void | setParam (const std::string &, double) |
| Assigns a parameter value to property functions of the integrand. | |
| virtual void | setParam (const std::string &, const Vec3 &) |
| Assigns parameter values to property functions of the integrand. | |
Protected Member Functions | |
| Integrand () | |
| The default constructor is protected to allow sub-classes only. | |
| virtual bool | evalInt (LocalIntegral &, const FiniteElement &fe, const Vec3 &) const |
| Evaluates the integrand at interior points for stationary problems. | |
| virtual bool | evalIntMx (LocalIntegral &, const MxFiniteElement &fe, const Vec3 &) const |
| Evaluates the integrand at interior points for stationary problems. | |
| virtual bool | evalInt (LocalIntegral &, const FiniteElement &fe, const Vec3 &, const Vec3 &) const |
| Evaluates the integrand at interface points, stationary problems. | |
| virtual bool | evalIntMx (LocalIntegral &, const MxFiniteElement &fe, const Vec3 &, const Vec3 &) const |
| Evaluates the integrand at interface points, stationary problems. | |
| virtual bool | evalBou (LocalIntegral &, const FiniteElement &, const Vec3 &, const Vec3 &) const |
| Evaluates the integrand at boundary points for stationary problems. | |
| virtual bool | evalBouMx (LocalIntegral &, const MxFiniteElement &, const Vec3 &, const Vec3 &) const |
| Evaluates the integrand at boundary points for stationary problems. | |
| virtual bool | finalizeElement (LocalIntegral &elmInt, const TimeDomain &, size_t) |
| Finalizes the element quantities after the numerical integration. More... | |
| virtual bool | finalizeElement (LocalIntegral &) |
| Finalizes the element quantities after the numerical integration. More... | |
Abstract base class representing a system level integrated quantity.
This class defines the interface between the finite element (FE) assembly drivers of the ASM-hierarchy and the problem-dependent classes containing all physical properties for the problem to be solved.
The interface consists of methods for evaluating the integrand at interior integration points (evalInt), and at boundary integration points (evalBou). The latter are used for Neumann boundary conditions, typically. The integrand evaluation methods have access to the FE basis function values and derivatives through the FiniteElement argument. There are also a set of methods dedicated for mixed field interpolation problems, which take a MxFiniteElement object as argument instead.
| enum Integrand::Traits |
Enum defining the additional terms that an Integrand may require.
| Enumerator | |
|---|---|
| STANDARD | Default integrand type (first derivatives only) |
| NO_DERIVATIVES | Integrand don't want any derivatives. |
| SECOND_DERIVATIVES | Integrand wants second derivatives. |
| THIRD_DERIVATIVES | Integrand wants third derivatives. |
| AVERAGE | Integrand wants basis function averages. |
| ELEMENT_CORNERS | Integrand wants element corner coordinates. |
| ELEMENT_CENTER | Integrand wants element center coordinates. |
| G_MATRIX | Integrand wants the G matrix. |
| NODAL_ROTATIONS | Integrand wants nodal rotation tensors. |
| XO_ELEMENTS | Integrand uses extraordinary elements. |
| INTERFACE_TERMS | Integrand has element interface terms. |
| NORMAL_DERIVS | Integrand uses p-order normal derivatives. |
| UPDATED_NODES | Integrand wants updated nodal coordinates. |
| PIOLA_MAPPING | Integrand wants Piola mapping. |
| POINT_DEFORMATION | Integrand wants point-wise deformation. |
|
inlinevirtual |
Evaluates the integrand at a boundary point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
| [in] | normal | Boundary normal vector at current integration point |
The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
Referenced by ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs2D::integrate(), ASMs2DLag::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu3D::integrate(), ASMs3D::integrateEdge(), ASMs3DLag::integrateEdge(), and ASMs3DSpec::integrateEdge().
|
inlinevirtual |
Evaluates the integrand at a boundary point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Mixed finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
| [in] | normal | Boundary normal vector at current integration point |
This interface is used for mixed formulations. The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
Referenced by ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), and ASMu2Dmx::integrate().
|
inlinevirtual |
Evaluates the integrand at an interior point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
Referenced by evalInt(), ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs1DSpec::integrate(), ASMs2D::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu3D::integrate(), and ASMs2DLag::integrateElm().
|
inlinevirtual |
Evaluates the integrand at an element interface point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
| [in] | normal | Interface normal vector at current integration point |
The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
References evalInt().
|
inlinevirtual |
Evaluates the integrand at an interior point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Mixed finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
This interface is used for mixed formulations only. The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
Referenced by evalIntMx(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMu2Dmx::integrate(), and ASMu3Dmx::integrate().
|
inlinevirtual |
Evaluates the integrand at an element interface point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | X | Cartesian coordinates of current integration point |
| [in] | normal | Interface normal vector at current integration point |
The default implementation forwards to the stationary version. Reimplement this method for time-dependent or non-linear problems.
References evalIntMx().
|
inlinevirtual |
Evaluates the dirac-delta integrand at a specified point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | pval | Function value at the specified point |
This interface can be used for implementation calculation of consistent load vectors due to a point load on the element.
Referenced by ASMstruct::diracPoint(), ASMu2D::diracPoint(), and ASMu3D::diracPoint().
|
inlineprotectedvirtual |
Finalizes the element quantities after the numerical integration.
Basic interface. Reimplement this method if no additional parameters are needed.
|
inlinevirtual |
Finalizes the element quantities after the numerical integration.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Nodal and integration point data for current element |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
| [in] | iGP | Global integration point counter of first point in element |
This method is invoked once for each element, after the numerical integration loop over interior points is finished and before the resulting element quantities are assembled into their system level equivalents. It can also be used to implement multiple integration point loops within the same element, provided the necessary integration point values are stored internally in the object during the first integration loop.
The default implementation forwards to the simple interface taking no FiniteElement argument. Reimplement this method if time domain parameters, element quantities and/or the integration point counter are needed.
Referenced by ASMstruct::diracPoint(), finalizeElement(), ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs2D::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), ASMu3D::integrate(), ASMu3Dmx::integrate(), and ASMs2DLag::integrateElm().
|
inlineprotectedvirtual |
Finalizes the element quantities after the numerical integration.
Simplified interface for when finite element data is not needed. The default implementation forwards to the basic interface taking no additional arguments. Reimplement this method if time domain parameters, and/or the integration point counter are needed.
References finalizeElement().
|
inlinevirtual |
Finalizes the element quantities after boundary integration.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Nodal and integration point data for current element |
| [in] | time | Parameters for nonlinear and time-dependent simulations |
This method is invoked once for each element, after the numerical integration loop over boundary points is finished and before the resulting element quantities are assembled into their system level equivalents. It can also be used to implement multiple integration point loops within the same element, provided the necessary integration point values are stored internally in the object during the first integration loop.
Referenced by ASMs2D::integrate(), ASMs2DLag::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), ASMu3D::integrate(), ASMs3D::integrateEdge(), ASMs3DLag::integrateEdge(), and ASMs3DSpec::integrateEdge().
|
inlinevirtual |
Returns a local integral contribution object for the given element.
| [in] | nen | Number of nodes on each basis |
| [in] | iEl | Global element number (1-based) |
| [in] | neumann | Whether or not we are assembling Neumann BCs |
This form is used for mixed formulations only. The default implementation just forwards to the single-basis version. Reimplement this method if your mixed formulation requires specialized local integral objects.
References getLocalIntegral().
|
pure virtual |
Returns a local integral contribution object for the given element.
| [in] | nen | Number of nodes on element |
| [in] | iEl | Global element number (1-based) |
| [in] | neumann | Whether or not we are assembling Neumann BCs |
Implemented in ForceBase, NormBase, IntegrandBase, and GlbL2.
Referenced by ASMstruct::diracPoint(), ASMu2D::diracPoint(), ASMu3D::diracPoint(), getLocalIntegral(), ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs1DSpec::integrate(), ASMs2D::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), ASMu3D::integrate(), ASMu3Dmx::integrate(), ASMs2DLag::integrate(), ASMs3D::integrateEdge(), ASMs3DLag::integrateEdge(), ASMs3DSpec::integrateEdge(), and ASMs2DLag::integrateElm().
|
pure virtual |
Initializes current element for numerical integration.
| [in] | MNPC | Matrix of nodal point correspondance for current element |
| [in] | fe | Nodal and integration point data for current element |
| [in] | X0 | Cartesian coordinates of the element center |
| [in] | nPt | Number of integration points in this element |
| elmInt | Local integral for element |
This method is invoked once before starting the numerical integration loop over the Gaussian quadrature points over an element. It is supposed to perform all the necessary internal initializations needed before the numerical integration is started for current element. Reimplement this method for problems requiring the element center and/or the number of integration points during/before the integrand evaluations.
Implemented in NormBase, IntegrandBase, ForceBase, and GlbL2.
Referenced by ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs1DSpec::integrate(), ASMs2D::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), ASMu3D::integrate(), ASMu3Dmx::integrate(), and ASMs2DLag::integrateElm().
|
pure virtual |
Initializes current element for numerical integration (mixed).
| [in] | MNPC | Nodal point correspondance for the bases |
| [in] | fe | Nodal and integration point data for current element |
| [in] | elem_sizes | Size of each basis on the element |
| [in] | basis_sizes | Size of each basis on the patch level |
| elmInt | Local integral for element |
Implemented in IntegrandBase, NormBase, ForceBase, and GlbL2.
|
pure virtual |
Initializes current element for numerical integration (mixed).
| [in] | MNPC | Nodal point correspondance for the bases |
| [in] | elem_sizes | Size of each basis on the element |
| [in] | basis_sizes | Size of each basis on the patch level |
| elmInt | Local integral for element |
Implemented in NormBase, IntegrandBase, ForceBase, and GlbL2.
|
pure virtual |
Initializes current element for numerical integration.
| [in] | MNPC | Matrix of nodal point correspondance for current element |
| elmInt | Local integral for element |
This method is invoked once before starting the numerical integration loop over the Gaussian quadrature points over an element. It is supposed to perform all the necessary internal initializations needed before the numerical integration is started for current element. Reimplement this method for problems not requiring the the element center nor the number of integration points before the integration loop is started.
Implemented in NormBase, IntegrandBase, ForceBase, and GlbL2.
|
pure virtual |
Initializes current element for boundary integration (mixed).
| [in] | MNPC | Nodal point correspondance for the bases |
| [in] | elem_sizes | Size of each basis on the element |
| [in] | basis_sizes | Size of each basis on the patch |
| elmInt | Local integral for element |
Implemented in ForceBase, NormBase, IntegrandBase, and GlbL2.
|
pure virtual |
Initializes current element for boundary integration.
| [in] | MNPC | Matrix of nodal point correspondance for current element |
| elmInt | Local integral for element |
Implemented in ForceBase, NormBase, IntegrandBase, and GlbL2.
Referenced by ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs2D::integrate(), ASMs2DLag::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DSpec::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMs3DSpec::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), ASMu3D::integrate(), ASMs3D::integrateEdge(), ASMs3DLag::integrateEdge(), and ASMs3DSpec::integrateEdge().
|
inlinevirtual |
Evaluates reduced integration terms at an interior point.
| elmInt | The local integral object to receive the contributions | |
| [in] | fe | Finite element data of current integration point |
| [in] | X | Cartesian coordinates of current integration point |
Reimplement this method to evaluate terms at other points than the regular integration points. This method is invoked in a separate loop prior to the main integration point loop.
Reimplemented in NormBase.
Referenced by ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs2D::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMu2D::integrate(), ASMu3D::integrate(), ASMs2DLag::integrateElm(), and NormBase::reducedInt().
|
inlinevirtual |
Defines the Neumann order that is the subject of integration.
This method is invoked once before the integration loop over a patch boundary for a certain property. Reimplement this method if your integrand has a differential operator of second (or higher) order, for which you need to have two (or more) distinct Neumann boundary conditions.
Referenced by ASMs1D::integrate(), ASMs1DLag::integrate(), ASMs2D::integrate(), ASMs2DLag::integrate(), ASMs2Dmx::integrate(), ASMs2DmxLag::integrate(), ASMs2DTri::integrate(), ASMs3D::integrate(), ASMs3DLag::integrate(), ASMs3Dmx::integrate(), ASMs3DmxLag::integrate(), ASMu2D::integrate(), ASMu2Dmx::integrate(), and ASMu3D::integrate().