OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
oofem::LIBeam3dNL2 Class Reference

This class implements a 3-dimensional Linear Isoparametric Mindlin theory beam element, with reduced integration. More...

#include <libeam3dnl2.h>

+ Inheritance diagram for oofem::LIBeam3dNL2:
+ Collaboration diagram for oofem::LIBeam3dNL2:

Public Member Functions

 LIBeam3dNL2 (int n, Domain *d)
 
virtual ~LIBeam3dNL2 ()
 
virtual void computeLumpedMassMatrix (FloatMatrix &answer, TimeStep *tStep)
 Computes lumped mass matrix of receiver. More...
 
virtual void computeConsistentMassMatrix (FloatMatrix &answer, TimeStep *tStep, double &mass, const double *ipDensity=NULL)
 Computes consistent mass matrix of receiver using numerical integration over element volume. More...
 
virtual void computeStrainVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Compute strain vector of receiver evaluated at given integration point at given time step from element displacement vector. More...
 
virtual int computeNumberOfDofs ()
 Computes or simply returns total number of element's local DOFs. More...
 
virtual void giveDofManDofIDMask (int inode, IntArray &answer) const
 Returns dofmanager dof mask for node. More...
 
virtual double computeVolumeAround (GaussPoint *gp)
 Returns volume related to given integration point. More...
 
virtual int computeGlobalCoordinates (FloatArray &answer, const FloatArray &lcoords)
 Computes the global coordinates from given element's local coordinates. More...
 
virtual const char * giveInputRecordName () const
 
virtual const char * giveClassName () const
 
virtual IRResultType initializeFrom (InputRecord *ir)
 Initializes receiver according to object description stored in input record. More...
 
virtual Element_Geometry_Type giveGeometryType () const
 Returns the element geometry type. More...
 
virtual void drawRawGeometry (oofegGraphicContext &gc, TimeStep *tStep)
 
virtual void drawDeformedGeometry (oofegGraphicContext &gc, TimeStep *tStep, UnknownType)
 
virtual void computeStiffnessMatrix (FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
 Computes the stiffness matrix of receiver. More...
 
virtual void giveInternalForcesVector (FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord=0)
 Evaluates nodal representation of real internal forces. More...
 
virtual integrationDomain giveIntegrationDomain () const
 Returns integration domain for receiver, used to initialize integration point over receiver volume. More...
 
virtual MaterialMode giveMaterialMode ()
 Returns material mode for receiver integration points. More...
 
virtual contextIOResultType saveContext (DataStream &stream, ContextMode mode, void *obj)
 Stores receiver state to output stream. More...
 
virtual contextIOResultType restoreContext (DataStream &stream, ContextMode mode, void *obj)
 Restores the receiver state previously written in stream. More...
 
- Public Member Functions inherited from oofem::NLStructuralElement
 NLStructuralElement (int n, Domain *d)
 Constructor. More...
 
virtual ~NLStructuralElement ()
 Destructor. More...
 
int giveGeometryMode ()
 Returns the geometry mode describing the formulation used in the internal work 0 - Engineering (small deformation) stress-strain mode 1 - First Piola-Kirchhoff - Deformation gradient mode, P is defined as FS 2 - Second Piola-Kirchhoff - Green-Lagrange strain mode with deformation gradient as input (deprecated and not supported) More...
 
void computeFirstPKStressVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Computes the first Piola-Kirchhoff stress tensor on Voigt format. More...
 
void computeCauchyStressVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Computes the Cauchy stress tensor on Voigt format. More...
 
virtual void computeInitialStressMatrix (FloatMatrix &answer, TimeStep *tStep)
 Computes the initial stiffness matrix of receiver. More...
 
void computeStiffnessMatrix_withIRulesAsSubcells (FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
 Computes the stiffness matrix of receiver. More...
 
void giveInternalForcesVector_withIRulesAsSubcells (FloatArray &answer, TimeStep *tStep, int useUpdatedGpRecord=0)
 Evaluates nodal representation of real internal forces. More...
 
virtual void computeDeformationGradientVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 Computes the deformation gradient in Voigt form at integration point ip and at time step tStep. More...
 
double computeCurrentVolume (TimeStep *tStep)
 Computes the current volume of element. More...
 
virtual void giveInputRecord (DynamicInputRecord &input)
 Setups the input record string of receiver. More...
 
- Public Member Functions inherited from oofem::StructuralElement
 StructuralElement (int n, Domain *d)
 Constructor. More...
 
virtual ~StructuralElement ()
 Destructor. More...
 
virtual void giveCharacteristicMatrix (FloatMatrix &answer, CharType, TimeStep *tStep)
 Computes characteristic matrix of receiver of requested type in given time step. More...
 
virtual void giveCharacteristicVector (FloatArray &answer, CharType type, ValueModeType mode, TimeStep *tStep)
 Computes characteristic vector of receiver of requested type in given time step. More...
 
virtual void computeMassMatrix (FloatMatrix &answer, TimeStep *tStep)
 Computes mass matrix of receiver. More...
 
virtual void giveMassMtrxIntegrationgMask (IntArray &answer)
 Returns mask indicating, which unknowns (their type and ordering is the same as element unknown vector) participate in mass matrix integration. More...
 
void computeStiffnessMatrix_withIRulesAsSubcells (FloatMatrix &answer, MatResponseMode rMode, TimeStep *tStep)
 
virtual void computeField (ValueModeType mode, TimeStep *tStep, const FloatArray &lcoords, FloatArray &answer)
 Computes the unknown vector interpolated at the specified local coordinates. More...
 
virtual int giveIPValue (FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
 Returns the integration point corresponding value in full form. More...
 
virtual void computeResultingIPTemperatureAt (FloatArray &answer, TimeStep *tStep, GaussPoint *gp, ValueModeType mode)
 Computes at given time (tStep) the the resulting temperature component array. More...
 
virtual void computeResultingIPEigenstrainAt (FloatArray &answer, TimeStep *tStep, GaussPoint *gp, ValueModeType mode)
 Computes at given time the resulting eigenstrain component array. More...
 
virtual int adaptiveUpdate (TimeStep *tStep)
 Updates the internal state variables stored in all IPs according to already mapped state. More...
 
virtual void updateInternalState (TimeStep *tStep)
 Updates element state after equilibrium in time step has been reached. More...
 
virtual int giveInternalStateAtNode (FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *tStep)
 Returns internal state variable (like stress,strain) at node of element in Reduced form, the way how is obtained is dependent on InternalValueType. More...
 
virtual void showSparseMtrxStructure (CharType mtrx, oofegGraphicContext &gc, TimeStep *tStep)
 Shows sparse structure. More...
 
virtual void showExtendedSparseMtrxStructure (CharType mtrx, oofegGraphicContext &gc, TimeStep *tStep)
 Shows extended sparse structure (for example, due to nonlocal interactions for tangent stiffness) More...
 
virtual void computeLoadVector (FloatArray &answer, BodyLoad *load, CharType type, ValueModeType mode, TimeStep *tStep)
 Computes the contribution of the given body load (volumetric). More...
 
virtual void computeBoundarySurfaceLoadVector (FloatArray &answer, BoundaryLoad *load, int boundary, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true)
 Computes the contribution of the given load at the given boundary surface in global coordinate system. More...
 
virtual void computeBoundaryEdgeLoadVector (FloatArray &answer, BoundaryLoad *load, int boundary, CharType type, ValueModeType mode, TimeStep *tStep, bool global=true)
 Computes the contribution of the given load at the given boundary edge. More...
 
virtual void computeEdgeNMatrix (FloatMatrix &answer, int boundaryID, const FloatArray &lcoords)
 computes edge interpolation matrix More...
 
virtual void computeSurfaceNMatrix (FloatMatrix &answer, int boundaryID, const FloatArray &lcoords)
 Computes surface interpolation matrix. More...
 
StructuralCrossSectiongiveStructuralCrossSection ()
 Helper function which returns the structural cross-section for the element. More...
 
virtual void createMaterialStatus ()
 
virtual void updateBeforeNonlocalAverage (TimeStep *tStep)
 Updates internal element state (in all integration points of receiver) before nonlocal averaging takes place. More...
 
virtual void giveNonlocalLocationArray (IntArray &locationArray, const UnknownNumberingScheme &us)
 Returns the "nonlocal" location array of receiver. More...
 
virtual void addNonlocalStiffnessContributions (SparseMtrx &dest, const UnknownNumberingScheme &s, TimeStep *tStep)
 Adds the "nonlocal" contribution to stiffness matrix, to account for nonlocality of material model. More...
 
- Public Member Functions inherited from oofem::Element
 Element (int n, Domain *aDomain)
 Constructor. More...
 
 Element (const Element &src)=delete
 
Elementoperator= (const Element &src)=delete
 
virtual ~Element ()
 Virtual destructor. More...
 
virtual void drawYourself (oofegGraphicContext &gc, TimeStep *tStep)
 
virtual void drawAnnotation (oofegGraphicContext &gc, TimeStep *tStep)
 
virtual void drawScalar (oofegGraphicContext &gc, TimeStep *tStep)
 
virtual void drawSpecial (oofegGraphicContext &gc, TimeStep *tStep)
 
virtual void giveLocalIntVarMaxMin (oofegGraphicContext &gc, TimeStep *tStep, double &emin, double &emax)
 
virtual int giveInternalStateAtSide (FloatArray &answer, InternalStateType type, InternalStateMode mode, int side, TimeStep *tStep)
 Returns internal state variable (like stress,strain) at side of element in Reduced form If side is possessing DOFs, otherwise recover techniques will not work due to absence of side-shape functions. More...
 
int giveLabel () const
 
int giveGlobalNumber () const
 
void setGlobalNumber (int num)
 Sets receiver globally unique number. More...
 
elementParallelMode giveParallelMode () const
 Return elementParallelMode of receiver. More...
 
void setParallelMode (elementParallelMode _mode)
 Sets parallel mode of element. More...
 
virtual elementParallelMode giveKnotSpanParallelMode (int) const
 Returns the parallel mode for particular knot span of the receiver. More...
 
int packUnknowns (DataStream &buff, TimeStep *tStep)
 Pack all necessary data of element (according to its parallel_mode) integration points into given communication buffer. More...
 
int unpackAndUpdateUnknowns (DataStream &buff, TimeStep *tStep)
 Unpack and updates all necessary data of element (according to its parallel_mode) integration points into given communication buffer. More...
 
int estimatePackSize (DataStream &buff)
 Estimates the necessary pack size to hold all packed data of receiver. More...
 
const IntArraygivePartitionList () const
 Returns partition list of receiver. More...
 
void setPartitionList (IntArray &pl)
 Sets partition list of receiver. More...
 
virtual double predictRelativeComputationalCost ()
 Returns the weight representing relative computational cost of receiver The reference element is triangular plane stress element with linear approximation, single integration point and linear isotropic material. More...
 
virtual double giveRelativeSelfComputationalCost ()
 Returns the weight representing relative computational cost of receiver The reference element is triangular plane stress element. More...
 
virtual double predictRelativeRedistributionCost ()
 Returns the relative redistribution cost of the receiver. More...
 
IntArraygiveBodyLoadArray ()
 Returns array containing load numbers of loads acting on element. More...
 
IntArraygiveBoundaryLoadArray ()
 Returns array containing load numbers of boundary loads acting on element. More...
 
virtual void printOutputAt (FILE *file, TimeStep *tStep)
 Prints output of receiver to stream, for given time step. More...
 
void giveLocationArray (IntArray &locationArray, const UnknownNumberingScheme &s, IntArray *dofIds=NULL) const
 Returns the location array (array of code numbers) of receiver for given numbering scheme. More...
 
void giveLocationArray (IntArray &locationArray, const IntArray &dofIDMask, const UnknownNumberingScheme &s, IntArray *dofIds=NULL) const
 
virtual void giveBoundaryLocationArray (IntArray &locationArray, const IntArray &bNodes, const UnknownNumberingScheme &s, IntArray *dofIds=NULL)
 Returns the location array for the boundary of the element. More...
 
virtual void giveBoundaryLocationArray (IntArray &locationArray, const IntArray &bNodes, const IntArray &dofIDMask, const UnknownNumberingScheme &s, IntArray *dofIds=NULL)
 
virtual int giveNumberOfDofs ()
 
virtual int giveNumberOfInternalDofManagers () const
 
virtual DofManagergiveInternalDofManager (int i) const
 Returns i-th internal element dof manager of the receiver. More...
 
virtual double giveCharacteristicValue (CharType type, TimeStep *tStep)
 Computes characteristic value of receiver of requested type in given time step. More...
 
virtual void computeTangentFromSurfaceLoad (FloatMatrix &answer, SurfaceLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep)
 Computes the tangent contribution of the given load at the given boundary. More...
 
virtual void computeTangentFromEdgeLoad (FloatMatrix &answer, EdgeLoad *load, int boundary, MatResponseMode rmode, TimeStep *tStep)
 Computes the tangent contribution of the given load at the given boundary. More...
 
const IntArraygiveBodyLoadList () const
 Returns receiver list of bodyloads. More...
 
const IntArraygiveBoundaryLoadList () const
 Returns receiver list of boundary loads. More...
 
void computeVectorOf (ValueModeType u, TimeStep *tStep, FloatArray &answer)
 Returns local vector of unknowns. More...
 
void computeVectorOf (const IntArray &dofIDMask, ValueModeType u, TimeStep *tStep, FloatArray &answer, bool padding=false)
 
void computeBoundaryVectorOf (const IntArray &bNodes, const IntArray &dofIDMask, ValueModeType u, TimeStep *tStep, FloatArray &answer, bool padding=false)
 Boundary version of computeVectorOf. More...
 
void computeVectorOf (PrimaryField &field, const IntArray &dofIDMask, ValueModeType u, TimeStep *tStep, FloatArray &answer, bool padding=false)
 Returns local vector of unknowns. More...
 
void computeVectorOfPrescribed (ValueModeType u, TimeStep *tStep, FloatArray &answer)
 Returns local vector of prescribed unknowns. More...
 
void computeVectorOfPrescribed (const IntArray &dofIDMask, ValueModeType type, TimeStep *tStep, FloatArray &answer)
 Returns local vector of prescribed unknowns. More...
 
virtual int computeNumberOfGlobalDofs ()
 Computes the total number of element's global dofs. More...
 
int computeNumberOfPrimaryMasterDofs ()
 Computes the total number of element's primary master DOFs. More...
 
virtual bool computeGtoLRotationMatrix (FloatMatrix &answer)
 Returns transformation matrix from global c.s. More...
 
virtual bool giveRotationMatrix (FloatMatrix &answer)
 Transformation matrices updates rotation matrix between element-local and primary DOFs, taking into account nodal c.s. More...
 
virtual bool computeDofTransformationMatrix (FloatMatrix &answer, const IntArray &nodes, bool includeInternal)
 Returns transformation matrix for DOFs from global coordinate system to local coordinate system in nodes. More...
 
virtual void giveInternalDofManDofIDMask (int inode, IntArray &answer) const
 Returns internal dofmanager dof mask for node. More...
 
virtual void giveElementDofIDMask (IntArray &answer) const
 Returns element dof mask for node. More...
 
virtual double computeVolumeAreaOrLength ()
 Computes the volume, area or length of the element depending on its spatial dimension. More...
 
double computeMeanSize ()
 Computes the size of the element defined as its length. More...
 
virtual double computeVolume ()
 Computes the volume. More...
 
virtual double computeArea ()
 Computes the area (zero for all but 2d geometries). More...
 
virtual void giveBoundaryEdgeNodes (IntArray &bNodes, int boundary)
 Returns list of receiver boundary nodes for given edge. More...
 
virtual void giveBoundarySurfaceNodes (IntArray &bNodes, int boundary)
 Returns list of receiver boundary nodes for given surface. More...
 
virtual IntegrationRulegiveBoundaryEdgeIntegrationRule (int order, int boundary)
 Returns boundary edge integration rule. More...
 
virtual IntegrationRulegiveBoundarySurfaceIntegrationRule (int order, int boundary)
 Returns boundary surface integration rule. More...
 
int giveDofManagerNumber (int i) const
 Translates local to global indices for dof managers. More...
 
const IntArraygiveDofManArray () const
 
void addDofManager (DofManager *dMan)
 
DofManagergiveDofManager (int i) const
 
NodegiveNode (int i) const
 Returns reference to the i-th node of element. More...
 
virtual ElementSidegiveSide (int i) const
 Returns reference to the i-th side of element. More...
 
virtual FEInterpolationgiveInterpolation () const
 
virtual FEInterpolationgiveInterpolation (DofIDItem id) const
 Returns the interpolation for the specific dof id. More...
 
virtual MaterialgiveMaterial ()
 
int giveMaterialNumber () const
 
CrossSectiongiveCrossSection ()
 
void setMaterial (int matIndx)
 Sets the material of receiver. More...
 
virtual void setCrossSection (int csIndx)
 Sets the cross section model of receiver. More...
 
virtual int giveNumberOfDofManagers () const
 
virtual int giveNumberOfNodes () const
 Returns number of nodes of receiver. More...
 
void setDofManagers (const IntArray &dmans)
 Sets receiver dofManagers. More...
 
void setBodyLoads (const IntArray &bodyLoads)
 Sets receiver bodyLoadArray. More...
 
void setIntegrationRules (std::vector< std::unique_ptr< IntegrationRule > > irlist)
 Sets integration rules. More...
 
virtual int giveIntegrationRuleLocalCodeNumbers (IntArray &answer, IntegrationRule &ie)
 Assembles the code numbers of given integration element (sub-patch) This is done by obtaining list of nonzero shape functions and by collecting the code numbers of nodes corresponding to these shape functions. More...
 
int giveRegionNumber ()
 
virtual void postInitialize ()
 Performs post initialization steps. More...
 
virtual void initializeYourself (TimeStep *timeStepWhenICApply)
 Initialization according to state given by initial conditions. More...
 
virtual bool isActivated (TimeStep *tStep)
 
virtual bool isCast (TimeStep *tStep)
 
virtual int giveSpatialDimension ()
 Returns the element spatial dimension (1, 2, or 3). More...
 
virtual int giveNumberOfBoundarySides ()
 
virtual int giveDefaultIntegrationRule () const
 Returns id of default integration rule. More...
 
virtual IntegrationRulegiveDefaultIntegrationRulePtr ()
 Access method for default integration rule. More...
 
int giveNumberOfIntegrationRules ()
 
virtual IntegrationRulegiveIntegrationRule (int i)
 
std::vector< std::unique_ptr< IntegrationRule > > & giveIntegrationRulesArray ()
 
virtual int testElementExtension (ElementExtension ext)
 Tests if the element implements required extension. More...
 
int giveGlobalIPValue (FloatArray &answer, GaussPoint *gp, InternalStateType type, TimeStep *tStep)
 
virtual double giveLengthInDir (const FloatArray &normalToCrackPlane)
 Default implementation returns length of element projection into specified direction. More...
 
virtual double giveCharacteristicLength (const FloatArray &normalToCrackPlane)
 Returns the size of element in the given direction, in some cases adjusted (e.g. More...
 
double giveCharacteristicLengthForPlaneElements (const FloatArray &normalToCrackPlane)
 Returns the size of element in the given direction if the direction is in the XY plane, otherwise gives the mean size defined as the square root of the element area. More...
 
double giveCharacteristicLengthForAxisymmElements (const FloatArray &normalToCrackPlane)
 Returns the size of an axisymmetric element in the given direction if the direction is in the XY plane, otherwise gives the mean distance vrom the symmetry axis multiplied by pi. More...
 
virtual double giveCharacteristicSize (GaussPoint *gp, FloatArray &normalToCrackPlane, ElementCharSizeMethod method)
 Returns characteristic element size for a given integration point and given direction. More...
 
virtual double giveParentElSize () const
 Returns the size (length, area or volume depending on element type) of the parent element. More...
 
virtual bool computeLocalCoordinates (FloatArray &answer, const FloatArray &gcoords)
 Computes the element local coordinates from given global coordinates. More...
 
virtual void computeMidPlaneNormal (FloatArray &answer, const GaussPoint *gp)
 Computes mid-plane normal of receiver at integration point. More...
 
virtual int adaptiveMap (Domain *oldd, TimeStep *tStep)
 Initializes the internal state variables stored in all IPs according to state in given domain. More...
 
virtual int mapStateVariables (Domain &iOldDom, const TimeStep &iTStep)
 Maps the internal state variables stored in all IPs from the old domain to the new domain. More...
 
virtual int adaptiveFinish (TimeStep *tStep)
 Finishes the mapping for given time step. More...
 
virtual void updateLocalNumbering (EntityRenumberingFunctor &f)
 Local renumbering support. More...
 
template<class T >
void ipEvaluator (T *src, void(T::*f)(GaussPoint *gp))
 Integration point evaluator, loops over receiver IP's and calls given function (passed as f parameter) on them. The IP is parameter to function f. More...
 
template<class T , class S >
void ipEvaluator (T *src, void(T::*f)(GaussPoint *, S &), S &_val)
 Integration point evaluator, loops over receiver IP's and calls given function (passed as f parameter) on them. The IP is parameter to function f as well as additional array. More...
 
- Public Member Functions inherited from oofem::FEMComponent
 FEMComponent (int n, Domain *d)
 Regular constructor, creates component with given number and belonging to given domain. More...
 
virtual ~FEMComponent ()
 Virtual destructor. More...
 
DomaingiveDomain () const
 
virtual void setDomain (Domain *d)
 Sets associated Domain. More...
 
int giveNumber () const
 
void setNumber (int num)
 Sets number of receiver. More...
 
virtual void printYourself ()
 Prints receiver state on stdout. Useful for debugging. More...
 
virtual InterfacegiveInterface (InterfaceType t)
 Interface requesting service. More...
 
std::string errorInfo (const char *func) const
 Returns string for prepending output (used by error reporting macros). More...
 

Protected Member Functions

virtual void giveEdgeDofMapping (IntArray &answer, int iEdge) const
 Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element DOFs. More...
 
virtual double computeEdgeVolumeAround (GaussPoint *gp, int iEdge)
 Computes volume related to integration point on local edge. More...
 
virtual int computeLoadLEToLRotationMatrix (FloatMatrix &answer, int iEdge, GaussPoint *gp)
 Returns transformation matrix from local edge c.s to element local coordinate system of load vector components. More...
 
virtual int computeLoadGToLRotationMtrx (FloatMatrix &answer)
 Returns transformation matrix from global coordinate system to local element coordinate system for element load vector components. More...
 
virtual void computeBodyLoadVectorAt (FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode)
 Computes the load vector due to body load acting on receiver, at given time step. More...
 
virtual void updateYourself (TimeStep *tStep)
 Updates element state after equilibrium in time step has been reached. More...
 
virtual void initForNewStep ()
 Initializes receivers state to new time step. More...
 
virtual void computeBmatrixAt (GaussPoint *gp, FloatMatrix &answer, int, int)
 Computes the geometrical matrix of receiver in given integration point. More...
 
virtual void computeNmatrixAt (const FloatArray &iLocCoord, FloatMatrix &answer)
 Computes interpolation matrix for element unknowns. More...
 
virtual void computeGaussPoints ()
 Initializes the array of integration rules member variable. More...
 
virtual void computeConstitutiveMatrixAt (FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
 Computes constitutive matrix of receiver. More...
 
virtual void computeStressVector (FloatArray &answer, const FloatArray &strain, GaussPoint *gp, TimeStep *tStep)
 Computes the stress vector of receiver at given integration point, at time step tStep. More...
 
virtual double computeLength ()
 Computes the length (zero for all but 1D geometries) More...
 
virtual int giveLocalCoordinateSystem (FloatMatrix &answer)
 Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy. More...
 
void updateTempQuaternion (TimeStep *tStep)
 Updates the temporary triad at the centre to the state identified by given solution step. More...
 
void computeTempCurv (FloatArray &answer, TimeStep *tStep)
 Compute the temporary curvature at the centre to the state identified by given solution step. More...
 
void computeSMtrx (FloatMatrix &answer, FloatArray &vec)
 Evaluates the S matrix from given vector vec. More...
 
void computeRotMtrx (FloatMatrix &answer, FloatArray &psi)
 Evaluates the rotation matrix for large rotations according to Rodrigues formula for given pseudovector psi. More...
 
void computeXMtrx (FloatMatrix &answer, TimeStep *tStep)
 Computes X matrix at given solution state. More...
 
void computeRotMtrxFromQuaternion (FloatMatrix &answer, FloatArray &q)
 Computes rotation matrix from given quaternion. More...
 
void computeQuaternionFromRotMtrx (FloatArray &answer, FloatMatrix &R)
 Computes the normalized quaternion from the given rotation matrix. More...
 
void computeXdVector (FloatArray &answer, TimeStep *tStep)
 Computes x_21' vector for given solution state. More...
 
- Protected Member Functions inherited from oofem::NLStructuralElement
virtual int checkConsistency ()
 Performs consistency check. More...
 
virtual void computeBHmatrixAt (GaussPoint *gp, FloatMatrix &answer)
 Computes a matrix which, multiplied by the column matrix of nodal displacements, gives the displacement gradient stored by columns. More...
 
- Protected Member Functions inherited from oofem::StructuralElement
virtual int giveNumberOfIPForMassMtrxIntegration ()
 Return desired number of integration points for consistent mass matrix computation, if required. More...
 
void condense (FloatMatrix *stiff, FloatMatrix *mass, FloatArray *load, IntArray *what)
 General service for condensation of stiffness and optionally load vector and mass or initial stress matrices of receiver. More...
 
virtual void setupIRForMassMtrxIntegration (IntegrationRule &iRule)
 Setup Integration Rule Gauss Points for Mass Matrix integration. More...
 
virtual void computePointLoadVectorAt (FloatArray &answer, Load *load, TimeStep *tStep, ValueModeType mode, bool global=true)
 Computes point load vector contribution of receiver for given load (should has BoundaryLoad Base). More...
 
virtual void giveSurfaceDofMapping (IntArray &answer, int iSurf) const
 Assembles surface dof mapping mask, which provides mapping between surface local DOFs and "global" element DOFs. More...
 
virtual IntegrationRuleGetSurfaceIntegrationRule (int order)
 
virtual double computeSurfaceVolumeAround (GaussPoint *gp, int iSurf)
 Computes volume related to integration point on local surface. More...
 
virtual int computeLoadLSToLRotationMatrix (FloatMatrix &answer, int iSurf, GaussPoint *gp)
 Returns transformation matrix from local surface c.s to element local coordinate system of load vector components. More...
 

Private Attributes

double l0
 Initial length. More...
 
FloatArray q
 Quaternion at the center (last equilibrated). More...
 
FloatArray tempQ
 Temporary quaternion at the center. More...
 
StateCounterType tempQCounter
 Time stamp of temporary centre quaternion. More...
 
int referenceNode
 Reference node. More...
 

Additional Inherited Members

- Protected Attributes inherited from oofem::NLStructuralElement
int nlGeometry
 Flag indicating if geometrical nonlinearities apply. More...
 
- Protected Attributes inherited from oofem::StructuralElement
std::unique_ptr< FloatArrayinitialDisplacements
 Initial displacement vector, describes the initial nodal displacements when element has been casted. More...
 
- Protected Attributes inherited from oofem::Element
int numberOfDofMans
 Number of dofmanagers. More...
 
IntArray dofManArray
 Array containing dofmanager numbers. More...
 
int material
 Number of associated material. More...
 
int crossSection
 Number of associated cross section. More...
 
IntArray bodyLoadArray
 Array containing indexes of loads (body loads and boundary loads are kept separately), that apply on receiver. More...
 
IntArray boundaryLoadArray
 
std::vector< std::unique_ptr< IntegrationRule > > integrationRulesArray
 List of integration rules of receiver (each integration rule contains associated integration points also). More...
 
FloatMatrix elemLocalCS
 Transformation material matrix, used in orthotropic and anisotropic materials, global->local transformation. More...
 
int activityTimeFunction
 Element activity time function. If defined, nonzero value indicates active receiver, zero value inactive element. More...
 
int globalNumber
 In parallel mode, globalNumber contains globally unique DoFManager number. More...
 
int numberOfGaussPoints
 Number of integration points as specified by nip. More...
 
elementParallelMode parallel_mode
 Determines the parallel mode of the element. More...
 
IntArray partitions
 List of partition sharing the shared element or remote partition containing remote element counterpart. More...
 
- Protected Attributes inherited from oofem::FEMComponent
int number
 Component number. More...
 
Domaindomain
 Link to domain object, useful for communicating with other FEM components. More...
 

Detailed Description

This class implements a 3-dimensional Linear Isoparametric Mindlin theory beam element, with reduced integration.

Geometric nonlinearities are taken into account. Based on Element due to Simo and Vu-Quoc, description taken from Crisfield monograph. Similar to Libeam3dNL, but rotational update is done using quaternions.

Definition at line 55 of file libeam3dnl2.h.

Constructor & Destructor Documentation

oofem::LIBeam3dNL2::LIBeam3dNL2 ( int  n,
Domain d 
)

Definition at line 57 of file libeam3dnl2.C.

References l0, oofem::Element::numberOfDofMans, referenceNode, and tempQCounter.

virtual oofem::LIBeam3dNL2::~LIBeam3dNL2 ( )
inlinevirtual

Definition at line 73 of file libeam3dnl2.h.

References computeLumpedMassMatrix().

Member Function Documentation

virtual void oofem::LIBeam3dNL2::computeBmatrixAt ( GaussPoint gp,
FloatMatrix answer,
int  lowerIndx,
int  upperIndx 
)
inlineprotectedvirtual

Computes the geometrical matrix of receiver in given integration point.

The product of this matrix (assembled at given integration point) and element displacement vector is element strain vector. If lowerIndx and upperIndx parameters are specified, answer is formed only for strains within this interval. This will affects the size of answer.

Parameters
gpIntegration point for which answer is computed.
answerGeometric matrix of receiver.
lowerIndxIf specified, answer is formed only for strain with index equal and greater than lowerIndx. This parameter has default value 1 (answer is formed from first strain).
upperIndxIf specified, answer is formed only for strain with index less and equal than upperIndx. This parameter has default value ALL_STRAINS (answer is formed for all strains).

Implements oofem::StructuralElement.

Definition at line 116 of file libeam3dnl2.h.

References computeConstitutiveMatrixAt(), computeGaussPoints(), computeLength(), computeNmatrixAt(), computeQuaternionFromRotMtrx(), computeRotMtrx(), computeRotMtrxFromQuaternion(), computeSMtrx(), computeStressVector(), computeTempCurv(), computeXdVector(), computeXMtrx(), giveLocalCoordinateSystem(), OOFEM_ERROR, and updateTempQuaternion().

void oofem::LIBeam3dNL2::computeBodyLoadVectorAt ( FloatArray answer,
Load load,
TimeStep tStep,
ValueModeType  mode 
)
protectedvirtual

Computes the load vector due to body load acting on receiver, at given time step.

Default implementation computes body load vector numerically as $ l=\int_V N^{\mathrm{T}} f \rho\;\mathrm{d}V $ using default integration rule. Result is transformed to global c.s.

Parameters
answerComputed load vector due to body load
loadBody load which contribution is computed.
tStepTime step.
modedetermines the response mode

Reimplemented from oofem::StructuralElement.

Definition at line 700 of file libeam3dnl2.C.

References oofem::StructuralElement::computeBodyLoadVectorAt(), oofem::CS_Area, oofem::Element::giveCrossSection(), and oofem::FloatArray::times().

Referenced by giveMaterialMode().

virtual void oofem::LIBeam3dNL2::computeConsistentMassMatrix ( FloatMatrix answer,
TimeStep tStep,
double &  mass,
const double *  ipDensity = NULL 
)
inlinevirtual

Computes consistent mass matrix of receiver using numerical integration over element volume.

Mass matrix is computed as $ M = \int_V N^{\mathrm{T}} \rho N dV $, where $ N $ is displacement approximation matrix. The number of necessary integration points is determined using this->giveNumberOfIPForMassMtrxIntegration service. Only selected degrees of freedom participate in integration of mass matrix. This is described using dof mass integration mask. This mask is obtained from this->giveMassMtrxIntegrationgMask service. The nonzero mask value at i-th position indicates that i-th element DOF participates in mass matrix computation. The result is in element local coordinate system.

Parameters
answerMass matrix.
tStepTime step.
massTotal mass of receiver.

Reimplemented from oofem::StructuralElement.

Definition at line 76 of file libeam3dnl2.h.

References computeLumpedMassMatrix(), and computeStrainVector().

void oofem::LIBeam3dNL2::computeConstitutiveMatrixAt ( FloatMatrix answer,
MatResponseMode  rMode,
GaussPoint gp,
TimeStep tStep 
)
protectedvirtual

Computes constitutive matrix of receiver.

Default implementation uses element cross section giveCharMaterialStiffnessMatrix service.

Parameters
answerConstitutive matrix.
rModeMaterial response mode of answer.
gpIntegration point for which constitutive matrix is computed.
tStepTime step.

Implements oofem::StructuralElement.

Definition at line 434 of file libeam3dnl2.C.

References oofem::StructuralCrossSection::give3dBeamStiffMtrx(), and oofem::StructuralElement::giveStructuralCrossSection().

Referenced by computeBmatrixAt(), and computeStiffnessMatrix().

double oofem::LIBeam3dNL2::computeEdgeVolumeAround ( GaussPoint gp,
int  iEdge 
)
protectedvirtual

Computes volume related to integration point on local edge.

Parameters
gpedge integration point
iEdgeedge number

Reimplemented from oofem::StructuralElement.

Definition at line 607 of file libeam3dnl2.C.

References computeLength(), oofem::GaussPoint::giveWeight(), and OOFEM_ERROR.

Referenced by giveMaterialMode().

void oofem::LIBeam3dNL2::computeGaussPoints ( )
protectedvirtual

Initializes the array of integration rules member variable.

Element can have multiple integration rules for different tasks. For example structural element family class uses this feature to implement transparent support for reduced and selective integration of some strain components. Must be defined by terminator classes.

See also
IntegrationRule

Reimplemented from oofem::Element.

Definition at line 422 of file libeam3dnl2.C.

References oofem::Element::giveCrossSection(), oofem::Element::integrationRulesArray, and oofem::CrossSection::setupIntegrationPoints().

Referenced by computeBmatrixAt().

int oofem::LIBeam3dNL2::computeGlobalCoordinates ( FloatArray answer,
const FloatArray lcoords 
)
virtual

Computes the global coordinates from given element's local coordinates.

Parameters
answerRequested global coordinates.
lcoordsLocal coordinates.
Returns
Nonzero if successful, zero otherwise.

Reimplemented from oofem::Element.

Definition at line 571 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::Node::giveCoordinate(), oofem::Element::giveNode(), and oofem::FloatArray::resize().

Referenced by computeNumberOfDofs().

double oofem::LIBeam3dNL2::computeLength ( )
protectedvirtual

Computes the length (zero for all but 1D geometries)

Returns
Element length.

Reimplemented from oofem::Element.

Definition at line 485 of file libeam3dnl2.C.

References oofem::Node::giveCoordinate(), oofem::Element::giveNode(), and l0.

Referenced by computeBmatrixAt(), computeEdgeVolumeAround(), computeLumpedMassMatrix(), computeVolumeAround(), and giveLocalCoordinateSystem().

int oofem::LIBeam3dNL2::computeLoadGToLRotationMtrx ( FloatMatrix answer)
protectedvirtual

Returns transformation matrix from global coordinate system to local element coordinate system for element load vector components.

If no transformation is necessary, answer is empty matrix (default);

Parameters
answerTransformation matrix.
Returns
Nonzero if transformation matrix is not empty matrix, zero otherwise.

Reimplemented from oofem::StructuralElement.

Definition at line 656 of file libeam3dnl2.C.

References oofem::FloatMatrix::at(), giveLocalCoordinateSystem(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by giveMaterialMode().

int oofem::LIBeam3dNL2::computeLoadLEToLRotationMatrix ( FloatMatrix answer,
int  iEdge,
GaussPoint gp 
)
protectedvirtual

Returns transformation matrix from local edge c.s to element local coordinate system of load vector components.

Necessary, because integration must be done in local coordinate system of entity (edge or surface). If no transformation is necessary, answer is empty matrix (default);

Parameters
answerComputed rotation matrix.
iEdgeEdge number.
gpIntegration point (point, where transformation is computed, useful for curved edges).
Returns
Nonzero if transformation matrix is not empty matrix, zero otherwise.

Reimplemented from oofem::StructuralElement.

Definition at line 685 of file libeam3dnl2.C.

References oofem::FloatMatrix::clear().

Referenced by giveMaterialMode().

void oofem::LIBeam3dNL2::computeLumpedMassMatrix ( FloatMatrix answer,
TimeStep tStep 
)
virtual

Computes lumped mass matrix of receiver.

Default implementation returns lumped consistent mass matrix. Then returns lumped mass transformed into nodal coordinate system. The lumping procedure zeroes all off-diagonal members and zeroes also all diagonal members corresponding to non-displacement DOFs. Such diagonal matrix is then rescaled, to preserve the element mass. Requires the computeNmatrixAt and giveMassMtrxIntegrationgMask services to be implemented.

Parameters
answerLumped mass matrix.
tStepTime step.

Reimplemented from oofem::StructuralElement.

Definition at line 505 of file libeam3dnl2.C.

References oofem::FloatMatrix::at(), computeLength(), oofem::CS_Area, oofem::CrossSection::give(), oofem::Element::giveCrossSection(), oofem::StructuralElement::giveStructuralCrossSection(), oofem::Element::integrationRulesArray, oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by computeConsistentMassMatrix(), and ~LIBeam3dNL2().

void oofem::LIBeam3dNL2::computeNmatrixAt ( const FloatArray iLocCoord,
FloatMatrix answer 
)
protectedvirtual

Computes interpolation matrix for element unknowns.

The order and meaning of unknowns is element dependent.

Parameters
iLocCoordLocal coordinates.
answerInterpolation matrix evaluated at gp.

Reimplemented from oofem::StructuralElement.

Definition at line 520 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by computeBmatrixAt().

virtual int oofem::LIBeam3dNL2::computeNumberOfDofs ( )
inlinevirtual

Computes or simply returns total number of element's local DOFs.

Must be defined by particular element.

Returns
Number of local DOFs of element.

Reimplemented from oofem::Element.

Definition at line 80 of file libeam3dnl2.h.

References computeGlobalCoordinates(), computeVolumeAround(), and giveDofManDofIDMask().

void oofem::LIBeam3dNL2::computeQuaternionFromRotMtrx ( FloatArray answer,
FloatMatrix R 
)
protected

Computes the normalized quaternion from the given rotation matrix.

Parameters
answerComputed quaternion.
RSource rotation matrix.

Definition at line 179 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), and oofem::FloatArray::resize().

Referenced by computeBmatrixAt(), and initializeFrom().

void oofem::LIBeam3dNL2::computeRotMtrx ( FloatMatrix answer,
FloatArray psi 
)
protected

Evaluates the rotation matrix for large rotations according to Rodrigues formula for given pseudovector psi.

Parameters
answerResult.
psiPseudovector.

Definition at line 88 of file libeam3dnl2.C.

References oofem::FloatMatrix::add(), oofem::FloatMatrix::at(), oofem::FloatMatrix::beProductOf(), oofem::FloatArray::computeNorm(), computeSMtrx(), oofem::FloatArray::giveSize(), OOFEM_ERROR, oofem::FloatMatrix::resize(), S, oofem::FloatMatrix::times(), and oofem::FloatMatrix::zero().

Referenced by computeBmatrixAt(), and computeTempCurv().

void oofem::LIBeam3dNL2::computeRotMtrxFromQuaternion ( FloatMatrix answer,
FloatArray q 
)
protected

Computes rotation matrix from given quaternion.

Parameters
answerReturned rotation matrix.
qInput quaternion

Definition at line 158 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::times().

Referenced by computeBmatrixAt(), computeStiffnessMatrix(), computeStrainVector(), computeTempCurv(), drawDeformedGeometry(), and giveInternalForcesVector().

void oofem::LIBeam3dNL2::computeSMtrx ( FloatMatrix answer,
FloatArray vec 
)
protected

Evaluates the S matrix from given vector vec.

Parameters
answerAssembled result.
vecSource vector.

Definition at line 69 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::FloatArray::giveSize(), OOFEM_ERROR, and oofem::FloatMatrix::resize().

Referenced by computeBmatrixAt(), computeRotMtrx(), computeStiffnessMatrix(), and computeXMtrx().

void oofem::LIBeam3dNL2::computeStiffnessMatrix ( FloatMatrix answer,
MatResponseMode  rMode,
TimeStep tStep 
)
virtual

Computes the stiffness matrix of receiver.

The response is evaluated using $ \int B_{\mathrm{H}}^{\mathrm{T}} D B_{\mathrm{H}} \;\mathrm{d}v $, where $ B_{\mathrm{H}} $ is the B-matrix which produces the displacement gradient vector $ H_{\mathrm{V}} $ when multiplied with the solution vector a. Reduced integration are taken into account.

Parameters
answerComputed stiffness matrix.
rModeResponse mode.
tStepTime step.
Todo:
We probably need overloaded function (like above) here as well.
Todo:
We probably need overloaded function (like above) here as well.
Todo:
Verify that it works with large deformations

Reimplemented from oofem::NLStructuralElement.

Definition at line 335 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::FloatMatrix::beProductOf(), oofem::FloatMatrix::beProductTOf(), oofem::FloatMatrix::clear(), computeConstitutiveMatrixAt(), computeRotMtrxFromQuaternion(), computeSMtrx(), computeStrainVector(), computeStressVector(), computeXdVector(), computeXMtrx(), oofem::IntegrationRule::getIntegrationPoint(), oofem::Element::giveDefaultIntegrationRulePtr(), l0, tempQ, oofem::FloatMatrix::times(), and updateTempQuaternion().

Referenced by giveGeometryType().

void oofem::LIBeam3dNL2::computeStrainVector ( FloatArray answer,
GaussPoint gp,
TimeStep tStep 
)
virtual

Compute strain vector of receiver evaluated at given integration point at given time step from element displacement vector.

The nature of strain vector depends on the element type.

Parameters
answerRequested strain vector.
gpIntegration point where to calculate the strain.
tStepTime step.

Reimplemented from oofem::StructuralElement.

Definition at line 227 of file libeam3dnl2.C.

References oofem::FloatArray::at(), computeRotMtrxFromQuaternion(), computeTempCurv(), computeXdVector(), l0, oofem::FloatArray::resize(), tempQ, and updateTempQuaternion().

Referenced by computeConsistentMassMatrix(), computeStiffnessMatrix(), and giveInternalForcesVector().

void oofem::LIBeam3dNL2::computeStressVector ( FloatArray answer,
const FloatArray strain,
GaussPoint gp,
TimeStep tStep 
)
protectedvirtual

Computes the stress vector of receiver at given integration point, at time step tStep.

The nature of these stresses depends on the element's type.

Parameters
answerStress vector.
strainStrain vector.
gpIntegration point.
tStepTime step.

Implements oofem::StructuralElement.

Definition at line 441 of file libeam3dnl2.C.

References oofem::StructuralCrossSection::giveGeneralizedStress_Beam3d(), and oofem::StructuralElement::giveStructuralCrossSection().

Referenced by computeBmatrixAt(), computeStiffnessMatrix(), and giveInternalForcesVector().

double oofem::LIBeam3dNL2::computeVolumeAround ( GaussPoint gp)
virtual

Returns volume related to given integration point.

Used typically in subroutines, that perform integration over element volume. Should be implemented by particular elements.

See also
GaussPoint
Parameters
gpIntegration point for which volume is computed.
Returns
Volume for integration point.

Reimplemented from oofem::Element.

Definition at line 555 of file libeam3dnl2.C.

References computeLength(), and oofem::GaussPoint::giveWeight().

Referenced by computeNumberOfDofs().

void oofem::LIBeam3dNL2::computeXdVector ( FloatArray answer,
TimeStep tStep 
)
protected

Computes x_21' vector for given solution state.

Parameters
answerReturned x_21'.
tStepDetermines solution state.

Definition at line 317 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::Element::computeVectorOf(), oofem::Node::giveCoordinate(), oofem::Element::giveNode(), and oofem::FloatArray::resize().

Referenced by computeBmatrixAt(), computeStiffnessMatrix(), computeStrainVector(), and computeXMtrx().

void oofem::LIBeam3dNL2::computeXMtrx ( FloatMatrix answer,
TimeStep tStep 
)
protected

Computes X matrix at given solution state.

Parameters
answerReturned X matrix.
tStepDetermines solution state.

Definition at line 257 of file libeam3dnl2.C.

References oofem::FloatMatrix::at(), computeSMtrx(), computeXdVector(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by computeBmatrixAt(), computeStiffnessMatrix(), and giveInternalForcesVector().

virtual const char* oofem::LIBeam3dNL2::giveClassName ( ) const
inlinevirtual
Returns
Class name of the receiver.

Reimplemented from oofem::NLStructuralElement.

Definition at line 88 of file libeam3dnl2.h.

References initializeFrom().

void oofem::LIBeam3dNL2::giveDofManDofIDMask ( int  inode,
IntArray answer 
) const
virtual

Returns dofmanager dof mask for node.

This mask defines the dofs which are used by element in node. Mask influences the code number ordering for particular node. Code numbers are ordered according to node order and dofs belonging to particular node are ordered according to this mask. If element requests dofs using node mask which are not in node then error is generated. This masking allows node to be shared by different elements with different dofs in same node. Elements local code numbers are extracted from node using this mask. Must be defined by particular element.

Parameters
inodeMask is computed for local dofmanager with inode number.
answerMask for node.

Reimplemented from oofem::Element.

Definition at line 565 of file libeam3dnl2.C.

Referenced by computeNumberOfDofs().

void oofem::LIBeam3dNL2::giveEdgeDofMapping ( IntArray answer,
int  iEdge 
) const
protectedvirtual

Assembles edge dof mapping mask, which provides mapping between edge local DOFs and "global" element DOFs.

Mask can be imagined as local edge code numbers used to localize local edge DOFs to element DOFs.

Parameters
answerEdge DOF mapping.
iEdgeEdge number.

Reimplemented from oofem::StructuralElement.

Definition at line 589 of file libeam3dnl2.C.

References oofem::IntArray::at(), OOFEM_ERROR, and oofem::IntArray::resize().

Referenced by giveMaterialMode().

virtual Element_Geometry_Type oofem::LIBeam3dNL2::giveGeometryType ( ) const
inlinevirtual

Returns the element geometry type.

This information is assumed to be of general interest, but it is required only for some specialized tasks.

Returns
Geometry type of element.

Reimplemented from oofem::Element.

Definition at line 90 of file libeam3dnl2.h.

References computeStiffnessMatrix(), drawDeformedGeometry(), drawRawGeometry(), gc, and giveInternalForcesVector().

virtual const char* oofem::LIBeam3dNL2::giveInputRecordName ( ) const
inlinevirtual
Returns
Input record name of the receiver.

Implements oofem::FEMComponent.

Definition at line 87 of file libeam3dnl2.h.

References _IFT_LIBeam3dNL2_Name.

virtual integrationDomain oofem::LIBeam3dNL2::giveIntegrationDomain ( ) const
inlinevirtual

Returns integration domain for receiver, used to initialize integration point over receiver volume.

Default behavior is taken from the default interpolation.

Returns
Integration domain of element.

Reimplemented from oofem::Element.

Definition at line 100 of file libeam3dnl2.h.

References oofem::_Line.

void oofem::LIBeam3dNL2::giveInternalForcesVector ( FloatArray answer,
TimeStep tStep,
int  useUpdatedGpRecord = 0 
)
virtual

Evaluates nodal representation of real internal forces.

Necessary transformations are taken into account.

Todo:
what is meant?
Parameters
answerEquivalent nodal forces vector.
tStepTime step
useUpdatedGpRecordIf equal to zero, the stresses in integration points are computed (slow but safe).
Todo:
Is this really what we should do for inactive elements?
Todo:
This is actaully inefficient since it constructs B and twice and collects the nodal unknowns over and over.
Todo:
is this check really necessary?

Reimplemented from oofem::NLStructuralElement.

Definition at line 282 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), oofem::FloatArray::beProductOf(), computeRotMtrxFromQuaternion(), computeStrainVector(), computeStressVector(), computeXMtrx(), oofem::IntegrationRule::getIntegrationPoint(), oofem::Element::giveDefaultIntegrationRulePtr(), oofem::GaussPoint::giveMaterialStatus(), tempQ, and updateTempQuaternion().

Referenced by giveGeometryType().

int oofem::LIBeam3dNL2::giveLocalCoordinateSystem ( FloatMatrix answer)
protectedvirtual

Returns local coordinate system of receiver Required by material models with ortho- and anisotrophy.

Returns a unit vectors of local coordinate system at element stored row-wise. If local system is equal to global one, set answer to empty matrix and return zero value.

Returns
nonzero if answer computed, zero value if answer is empty, i.e. no transformation is necessary.

Reimplemented from oofem::Element.

Definition at line 619 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::FloatMatrix::at(), computeLength(), oofem::Node::giveCoordinate(), oofem::FEMComponent::giveDomain(), oofem::Domain::giveNode(), oofem::Element::giveNode(), referenceNode, oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by computeBmatrixAt(), computeLoadGToLRotationMtrx(), and initializeFrom().

virtual MaterialMode oofem::LIBeam3dNL2::giveMaterialMode ( )
inlinevirtual

Returns material mode for receiver integration points.

Should be specialized.

Returns
Material mode of element.

Reimplemented from oofem::Element.

Definition at line 101 of file libeam3dnl2.h.

References computeBodyLoadVectorAt(), computeEdgeVolumeAround(), computeLoadGToLRotationMtrx(), computeLoadLEToLRotationMatrix(), giveEdgeDofMapping(), initForNewStep(), restoreContext(), saveContext(), and updateYourself().

void oofem::LIBeam3dNL2::initForNewStep ( )
protectedvirtual

Initializes receivers state to new time step.

It can be used also if current time step must be re-started. Default implementation invokes initForNewStep member function for all defined integrationRules. Thus all state variables in all defined integration points are re initialized.

See also
IntegrationRule::initForNewStep.

Reimplemented from oofem::Element.

Definition at line 725 of file libeam3dnl2.C.

References oofem::Element::initForNewStep(), q, and tempQ.

Referenced by giveMaterialMode().

IRResultType oofem::LIBeam3dNL2::initializeFrom ( InputRecord ir)
virtual

Initializes receiver according to object description stored in input record.

This function is called immediately after creating object using constructor. Input record can be imagined as data record in component database belonging to receiver. Receiver may use value-name extracting functions to extract particular field from record.

See also
IR_GIVE_FIELD
IR_GIVE_OPTIONAL_FIELD
Parameters
irInput record to initialize from.
Returns
IRResultType
Todo:
Move this to postInitialize?

Reimplemented from oofem::NLStructuralElement.

Definition at line 448 of file libeam3dnl2.C.

References _IFT_LIBeam3dNL2_refnode, oofem::FloatMatrix::beTranspositionOf(), computeQuaternionFromRotMtrx(), giveLocalCoordinateSystem(), oofem::NLStructuralElement::initializeFrom(), IR_GIVE_FIELD, oofem::IRRT_OK, OOFEM_ERROR, q, and referenceNode.

Referenced by giveClassName().

contextIOResultType oofem::LIBeam3dNL2::restoreContext ( DataStream stream,
ContextMode  mode,
void *  obj 
)
virtual

Restores the receiver state previously written in stream.

See also
saveContext
Parameters
streamInput stream.
modeDetermines amount of info available in stream (state, definition, ...).
objSpecial parameter for sending extra information.
Returns
contextIOResultType.
Exceptions
throwsan ContextIOERR exception if error encountered.

Reimplemented from oofem::Element.

Definition at line 862 of file libeam3dnl2.C.

References oofem::CIO_OK, q, oofem::Element::restoreContext(), oofem::FloatArray::restoreYourself(), and THROW_CIOERR.

Referenced by giveMaterialMode().

contextIOResultType oofem::LIBeam3dNL2::saveContext ( DataStream stream,
ContextMode  mode,
void *  obj 
)
virtual

Stores receiver state to output stream.

Parameters
streamOutput stream.
modeDetermines amount of info required in stream (state, definition, ...).
objSpecial parameter, used only to send particular integration point to material class version of this method.
Returns
contextIOResultType.
Exceptions
throwsan ContextIOERR exception if error encountered.

Reimplemented from oofem::Element.

Definition at line 842 of file libeam3dnl2.C.

References oofem::CIO_OK, q, oofem::Element::saveContext(), oofem::FloatArray::storeYourself(), and THROW_CIOERR.

Referenced by giveMaterialMode().

void oofem::LIBeam3dNL2::updateTempQuaternion ( TimeStep tStep)
protected

Updates the temporary triad at the centre to the state identified by given solution step.

The attribute tempQ is changed to reflect new state and tempQCounter is set to solution step counter to avoid multiple updates.

Parameters
tStepSolution step identifying reached state.

Definition at line 118 of file libeam3dnl2.C.

References oofem::FloatArray::at(), oofem::Element::computeVectorOf(), oofem::TimeStep::giveSolutionStateCounter(), q, tempQ, and tempQCounter.

Referenced by computeBmatrixAt(), computeStiffnessMatrix(), computeStrainVector(), giveInternalForcesVector(), and updateYourself().

void oofem::LIBeam3dNL2::updateYourself ( TimeStep tStep)
protectedvirtual

Updates element state after equilibrium in time step has been reached.

Default implementation updates all integration rules defined by integrationRulesArray member variable. Doing this, all integration points and their material statuses are updated also. All temporary history variables, which now describe equilibrium state are copied into equilibrium ones. The existing internal state is used for update.

Parameters
tStepTime step for newly reached state.
See also
Material::updateYourself
IntegrationRule::updateYourself
GaussPoint::updateYourself
Element::updateInternalState

Reimplemented from oofem::StructuralElement.

Definition at line 709 of file libeam3dnl2.C.

References q, tempQ, updateTempQuaternion(), and oofem::StructuralElement::updateYourself().

Referenced by giveMaterialMode().

Member Data Documentation

double oofem::LIBeam3dNL2::l0
private
FloatArray oofem::LIBeam3dNL2::q
private

Quaternion at the center (last equilibrated).

Definition at line 61 of file libeam3dnl2.h.

Referenced by computeTempCurv(), drawDeformedGeometry(), initForNewStep(), initializeFrom(), restoreContext(), saveContext(), updateTempQuaternion(), and updateYourself().

int oofem::LIBeam3dNL2::referenceNode
private

Reference node.

Definition at line 69 of file libeam3dnl2.h.

Referenced by giveLocalCoordinateSystem(), initializeFrom(), and LIBeam3dNL2().

FloatArray oofem::LIBeam3dNL2::tempQ
private

Temporary quaternion at the center.

Definition at line 63 of file libeam3dnl2.h.

Referenced by computeStiffnessMatrix(), computeStrainVector(), giveInternalForcesVector(), initForNewStep(), updateTempQuaternion(), and updateYourself().

StateCounterType oofem::LIBeam3dNL2::tempQCounter
private

Time stamp of temporary centre quaternion.

Definition at line 67 of file libeam3dnl2.h.

Referenced by LIBeam3dNL2(), and updateTempQuaternion().


The documentation for this class was generated from the following files:

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:37 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011