G+Smo  23.12.0
Geometry + Simulation Modules
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
gsThinShellAssemblerBase< T > Class Template Referenceabstract

Detailed Description

template<class T>
class gismo::gsThinShellAssemblerBase< T >

Base class for the gsThinShellAssembler.

Template Parameters
TReal type
+ Inheritance diagram for gsThinShellAssemblerBase< T >:

Public Types

typedef gsBoxTopology::ifContainer ifContainer
 Default deconstructor.
 

Public Member Functions

virtual ThinShellAssemblerStatus assemble ()=0
 Assembles the linear system and corresponding right-hand side.
 
virtual ThinShellAssemblerStatus assemble (const gsFunctionSet< T > &deformed, const bool matrix=true, const bool homogenize=true)=0
 Assembles the tangential stiffness matrix and the residual for an iteration of Newton's method. More...
 
virtual ThinShellAssemblerStatus assemble (const gsMatrix< T > &solVector, const bool matrix=true, const bool homogenize=true)=0
 Assembles the tangential stiffness matrix and the residual for an iteration of Newton's method. More...
 
virtual ThinShellAssemblerStatus assembleFoundation ()=0
 Assembles the elastic foundation matrix.
 
virtual ThinShellAssemblerStatus assembleMass (const bool lumped=false)=0
 Assembles the mass matrix (including density and thickness!); if lumped=true, a lumped mass matrix will be constructed,.
 
virtual ThinShellAssemblerStatus assembleMatrix (const gsFunctionSet< T > &deformed)=0
 Assembles the tangential stiffness matrix (nonlinear) More...
 
virtual ThinShellAssemblerStatus assembleMatrix (const gsMatrix< T > &solVector)=0
 Assembles the tangential stiffness matrix (nonlinear) More...
 
virtual ThinShellAssemblerStatus assembleMatrix (const gsFunctionSet< T > &deformed, const gsFunctionSet< T > &previous, gsMatrix< T > &update)=0
 Assembles the tangential stiffness matrix (nonlinear) using the Mixed Integration Point (MIP) method. More...
 
virtual ThinShellAssemblerStatus assembleMatrix (const gsMatrix< T > &solVector, const gsMatrix< T > &prevVector)=0
 Assembles the tangential stiffness matrix (nonlinear) using the Mixed Integration Point (MIP) method. More...
 
virtual ThinShellAssemblerStatus assemblePressureMatrix (const gsFunction< T > &pressFun)=0
 Assembles the pressure contribution in the system matrix (linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureMatrix (const T pressure)=0
 Assembles the pressure contribution in the system matrix (linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureMatrix (const gsFunction< T > &pressFun, const gsFunctionSet< T > &deformed)=0
 Assembles the pressure contribution in the system matrix (non-linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureMatrix (const T pressure, const gsFunctionSet< T > &deformed)=0
 Assembles the pressure contribution in the system matrix (non-linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureVector (const gsFunction< T > &pressFun)=0
 Assembles the pressure contribution in the system vector (linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureVector (const T pressure)=0
 Assembles the pressure contribution in the system vector (linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureVector (const gsFunction< T > &pressFun, const gsFunctionSet< T > &deformed)=0
 Assembles the pressure contribution in the system vector (non-linear) More...
 
virtual ThinShellAssemblerStatus assemblePressureVector (const T pressure, const gsFunctionSet< T > &deformed)=0
 Assembles the pressure contribution in the system vector (non-linear) More...
 
virtual gsExprAssembler< T > assembler ()=0
 Returns the internal expression assembler.
 
virtual ThinShellAssemblerStatus assembleVector (const gsFunctionSet< T > &deformed, const bool homogenize=true)=0
 Assembles the residual vector. More...
 
virtual ThinShellAssemblerStatus assembleVector (const gsMatrix< T > &solVector, const bool homogenize=true)=0
 Assembles the residual vector. More...
 
virtual gsMatrix< T > boundaryForce (const gsFunctionSet< T > &deformed, const std::vector< patchSide > &patchSides) const =0
 Computes the force on a set of boundaries. More...
 
virtual gsMatrix< T > computePrincipalStresses (const gsMatrix< T > &points, const gsFunctionSet< T > &deformed, const T z=0)=0
 Compute the principal stretches in points given a deformed geometry. Optionally, the stretches can be computed on through-thickness coordinate z.
 
virtual gsMatrix< T > computePrincipalStretches (const gsMatrix< T > &points, const gsFunctionSet< T > &deformed, const T z=0)=0
 Compute the principal stretches in points given a deformed geometry. Optionally, the stretches can be computed on through-thickness coordinate z.
 
virtual gsMultiPatch< T > constructDisplacement (const gsMatrix< T > &solVector) const =0
 Construct displacement field from computed solution vector solVector and returns a multipatch.
 
virtual void constructDisplacement (const gsMatrix< T > &solVector, gsMultiPatch< T > &deformed) const =0
 Construct displacement field from computed solution vector solVector and returns the result in deformed.
 
virtual gsMultiPatch< T > constructMultiPatch (const gsMatrix< T > &solVector) const =0
 Construct solution field from computed solution vector solVector and returns a multipatch.
 
virtual gsMultiPatch< T > constructSolution (const gsMatrix< T > &solVector) const =0
 Construct deformed shell geometry from computed solution vector solVector and returns a multipatch.
 
virtual void constructSolution (const gsMatrix< T > &solVector, gsMultiPatch< T > &deformed) const =0
 Construct deformed shell geometry from computed solution vector solVector and returns the result in deformed.
 
virtual gsVector< T > constructSolutionVector (const gsMultiPatch< T > &deformed) const =0
 Reconstruct the solution vector based on the currently stored boundary conditions (thus the mapper).
 
virtual void constructStress (const gsFunctionSet< T > &deformed, gsPiecewiseFunction< T > &result, stress_type::type type)=0
 Construct Cauchy stress tensor for visualization.
 
virtual gsMatrix< T > fullSolutionVector (const gsMatrix< T > &vector) const =0
 Reconstruct the solution vector based on the currently stored boundary conditions (thus the mapper).
 
virtual const gsMultiPatch< T > & geometry () const =0
 Returns the undeformed geometry.
 
virtual T getArea (const gsFunctionSet< T > &geometry)=0
 Returns the area of geometry.
 
virtual const gsMultiBasis< T > & getBasis ()=0
 Gets the basis.
 
virtual T getDisplacementNorm (const gsFunctionSet< T > &deformed)=0
 Returns the displacement norm, i.e. norm = sqrt(u'*u/area)
 
virtual T getElasticEnergy (const gsFunctionSet< T > &deformed)=0
 Returns the elastic energy norm, i.e. norm = 0.5 * u'*F_int.
 
virtual gsDofMapper getMapper ()=0
 Returns the gsDofMapper.
 
virtual const gsPointLoads< T > & getPointLoads ()=0
 Gets the registered point loads.
 
virtual const gsPointLoads< T > & getPointMass ()=0
 Gets the registered point masses.
 
virtual const gsFunctionSet< T > & getSpaceBasis ()=0
 Get the basis that is used for assembly (but not for quadrature!)
 
 gsThinShellAssemblerBase ()
 Default deconstructor.
 
virtual void homogenizeDirichlet ()=0
 Sets the Dirichlet BCs to zero.
 
virtual T interfaceErrorC0 (const gsFunctionSet< T > &deformed)=0
 Returns the C1 error over the interface.
 
virtual T interfaceErrorG1 (const gsFunctionSet< T > &deformed)=0
 Returns the G1 error over the interface.
 
virtual T interfaceErrorGaussCurvature (const gsFunctionSet< T > &deformed)=0
 Returns the Gaussian curvature error over the interface.
 
virtual T interfaceErrorMeanCurvature (const gsFunctionSet< T > &deformed)=0
 Returns the mean curvature error over the interface.
 
virtual T interfaceErrorNormal (const gsFunctionSet< T > &deformed)=0
 Returns the normal vector error over the interface.
 
virtual gsSparseMatrix< T > & massMatrix ()=0
 Returns a reference to the mass matrix that is assembled.
 
virtual gsMaterialMatrixBase< T > * material (const index_t p) const =0
 Returns the material matrix on patch p used in the class.
 
virtual
gsMaterialMatrixContainer< T > 
materials () const =0
 Returns the material matrices used in the class.
 
virtual const gsSparseMatrix< T > & matrix () const =0
 Returns a reference to the system matrix that is assembled.
 
virtual index_t numDofs () const =0
 Returns the number of degrees of freedom in the assembler.
 
virtual gsOptionListoptions ()=0
 Returns the options of the assembler.
 
virtual gsMatrix< T > projectL2 (const gsFunction< T > &fun)=0
 Projects function fun on the basis and geometry stored in the class and returns the coefficients as a matrix.
 
virtual void projectL2_into (const gsFunction< T > &fun, gsMatrix< T > &result)=0
 Projects function fun on the basis and geometry stored in the class and returns the coefficients in result.
 
virtual void projectL2_into (const gsFunction< T > &fun, gsMultiPatch< T > &result)=0
 Projects function fun on the basis and geometry stored in the class and returns a multipatch in result.
 
virtual const gsMatrix< T > & rhs () const =0
 Returns a reference to the right-hand side vector that is assembled.
 
virtual void setBasis (const gsMultiBasis< T > &basis)=0
 Overwrites the basis. More...
 
virtual void setFoundation (const gsFunction< T > &foundation)=0
 Registers a stiffness function to be used for handling an elastic foundation, only relevant for 3D shells, with out-of-plane deformations. More...
 
virtual void setOptions (gsOptionList &options)=0
 Sets the options of the assembler.
 
virtual void setPointLoads (const gsPointLoads< T > &pLoads)=0
 Registers a gsPointLoads object for point loads acting on the shell.
 
virtual void setPointMass (const gsPointLoads< T > &pLoads)=0
 Registers a gsPointLoads object for a point mass acting on the shell. The point masss must be 1-dimensional.
 
virtual void setPressure (const gsFunction< T > &pressure)=0
 Registers a scalar function acting as pressure (in normal direction) in normal direction. More...
 
virtual void setSpaceBasis (const gsFunctionSet< T > &spaceBasis)=0
 Set the basis that is used for assembly (but not for quadrature!)
 
virtual void setUndeformed (const gsMultiPatch< T > &patches)=0
 Overwrites the undeformed geometry. More...
 
virtual ThinShellAssemblerStatus status () const =0
 Returns the assembler status.
 
virtual void updateBCs (const gsBoundaryConditions< T > &bconditions)=0
 Overwrites the boundary conditions. More...
 
virtual ~gsThinShellAssemblerBase ()
 Default empty constructor.
 

Member Function Documentation

virtual ThinShellAssemblerStatus assemble ( const gsFunctionSet< T > &  deformed,
const bool  matrix = true,
const bool  homogenize = true 
)
pure virtual

Assembles the tangential stiffness matrix and the residual for an iteration of Newton's method.

Note
.rhs() returns the negative residual!
Parameters
[in]deformedThe deformed multipatch
[in]MatrixTrue if the matrix should be assembled

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemble ( const gsMatrix< T > &  solVector,
const bool  matrix = true,
const bool  homogenize = true 
)
pure virtual

Assembles the tangential stiffness matrix and the residual for an iteration of Newton's method.

Note
.rhs() returns the negative residual!
Parameters
[in]deformedThe solution vector
[in]MatrixTrue if the matrix should be assembled

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleMatrix ( const gsFunctionSet< T > &  deformed)
pure virtual

Assembles the tangential stiffness matrix (nonlinear)

Parameters
[in]deformedThe deformed geometry

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleMatrix ( const gsMatrix< T > &  solVector)
pure virtual

Assembles the tangential stiffness matrix (nonlinear)

Parameters
[in]deformedThe solution vector

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleMatrix ( const gsFunctionSet< T > &  deformed,
const gsFunctionSet< T > &  previous,
gsMatrix< T > &  update 
)
pure virtual

Assembles the tangential stiffness matrix (nonlinear) using the Mixed Integration Point (MIP) method.

For more details, see Leonetti, L., Magisano, D., Madeo, A., Garcea, G., Kiendl, J., & Reali, A. (2019). A simplified Kirchhoff–Love large deformation model for elastic shells and its effective isogeometric formulation. Computer Methods in Applied Mechanics and Engineering, 354, 369–396. https://doi.org/10.1016/j.cma.2019.05.025

Parameters
[in]deformedThe deformed geometry
[in]previousThe previous geometry
updateThe update vector

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleMatrix ( const gsMatrix< T > &  solVector,
const gsMatrix< T > &  prevVector 
)
pure virtual

Assembles the tangential stiffness matrix (nonlinear) using the Mixed Integration Point (MIP) method.

For more details, see Leonetti, L., Magisano, D., Madeo, A., Garcea, G., Kiendl, J., & Reali, A. (2019). A simplified Kirchhoff–Love large deformation model for elastic shells and its effective isogeometric formulation. Computer Methods in Applied Mechanics and Engineering, 354, 369–396. https://doi.org/10.1016/j.cma.2019.05.025

Parameters
[in]solVectorThe current solution vector
[in]prevVectorThe previous solution vector

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureMatrix ( const gsFunction< T > &  pressFun)
pure virtual

Assembles the pressure contribution in the system matrix (linear)

Parameters
[in]pressFunThe pressure function

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureMatrix ( const T  pressure)
pure virtual

Assembles the pressure contribution in the system matrix (linear)

Parameters
[in]pressureThe pressure value

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureMatrix ( const gsFunction< T > &  pressFun,
const gsFunctionSet< T > &  deformed 
)
pure virtual

Assembles the pressure contribution in the system matrix (non-linear)

Parameters
[in]pressFunThe pressure function
[in]deformedThe deformed shape

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureMatrix ( const T  pressure,
const gsFunctionSet< T > &  deformed 
)
pure virtual

Assembles the pressure contribution in the system matrix (non-linear)

Parameters
[in]pressureThe pressure value
[in]deformedThe deformed shape

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureVector ( const gsFunction< T > &  pressFun)
pure virtual

Assembles the pressure contribution in the system vector (linear)

Parameters
[in]pressFunThe pressure function

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureVector ( const T  pressure)
pure virtual

Assembles the pressure contribution in the system vector (linear)

Parameters
[in]pressureThe pressure value

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureVector ( const gsFunction< T > &  pressFun,
const gsFunctionSet< T > &  deformed 
)
pure virtual

Assembles the pressure contribution in the system vector (non-linear)

Parameters
[in]pressFunThe pressure value
[in]deformedThe deformed shape

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assemblePressureVector ( const T  pressure,
const gsFunctionSet< T > &  deformed 
)
pure virtual

Assembles the pressure contribution in the system vector (non-linear)

Parameters
[in]pressureThe pressure value
[in]deformedThe deformed shape

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleVector ( const gsFunctionSet< T > &  deformed,
const bool  homogenize = true 
)
pure virtual

Assembles the residual vector.

Parameters
[in]deformedThe deformed geometry

Implemented in gsThinShellAssembler< d, T, bending >.

virtual ThinShellAssemblerStatus assembleVector ( const gsMatrix< T > &  solVector,
const bool  homogenize = true 
)
pure virtual

Assembles the residual vector.

Parameters
[in]deformedThe solution vector

Implemented in gsThinShellAssembler< d, T, bending >.

virtual gsMatrix<T> boundaryForce ( const gsFunctionSet< T > &  deformed,
const std::vector< patchSide > &  patchSides 
) const
pure virtual

Computes the force on a set of boundaries.

        This function is typically used when you want to know the
        load on a boundary on which a displacement is applied
Parameters
[in]deformedThe deformed geometry
patchSidespatch sides
Returns
The sim of the loads on the control points.

Implemented in gsThinShellAssembler< d, T, bending >.

virtual void setBasis ( const gsMultiBasis< T > &  basis)
pure virtual

Overwrites the basis.

Parameters
[in]bconditionsThe basis

Implemented in gsThinShellAssembler< d, T, bending >.

virtual void setFoundation ( const gsFunction< T > &  foundation)
pure virtual

Registers a stiffness function to be used for handling an elastic foundation, only relevant for 3D shells, with out-of-plane deformations.

Parameters
[in]foundationThe foundation stiffness function (3D vector)

Implemented in gsThinShellAssembler< d, T, bending >.

virtual void setPressure ( const gsFunction< T > &  pressure)
pure virtual

Registers a scalar function acting as pressure (in normal direction) in normal direction.

Note
Since the pressure acts in normal direction, and since the normals of the shell change under deformation, this load is nonlinear!
Parameters
[in]pressureThe pressure

Implemented in gsThinShellAssembler< d, T, bending >.

virtual void setUndeformed ( const gsMultiPatch< T > &  patches)
pure virtual

Overwrites the undeformed geometry.

Parameters
[in]bconditionsThe undeformed geometry

Implemented in gsThinShellAssembler< d, T, bending >.

virtual void updateBCs ( const gsBoundaryConditions< T > &  bconditions)
pure virtual

Overwrites the boundary conditions.

Parameters
[in]bconditionsThe boundary conditions

Implemented in gsThinShellAssembler< d, T, bending >.