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

Class implementing "simple" cross section model in finite element problem. More...

#include <simplecrosssection.h>

+ Inheritance diagram for oofem::SimpleCrossSection:
+ Collaboration diagram for oofem::SimpleCrossSection:

Public Member Functions

 SimpleCrossSection (int n, Domain *d)
 Constructor. More...
 
virtual void giveRealStress_3d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveRealStress_3dDegeneratedShell (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveRealStress_PlaneStrain (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveRealStress_PlaneStress (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveRealStress_1d (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveRealStress_Warping (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 
virtual void giveStiffnessMatrix_3d (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing the stiffness matrix. More...
 
virtual void giveStiffnessMatrix_PlaneStress (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 
virtual void giveStiffnessMatrix_PlaneStrain (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 
virtual void giveStiffnessMatrix_1d (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 
virtual void giveGeneralizedStress_Beam2d (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 Computes the generalized stress vector for given strain and integration point. More...
 
virtual void giveGeneralizedStress_Beam3d (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
virtual void giveGeneralizedStress_Plate (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
virtual void giveGeneralizedStress_Shell (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
virtual void giveGeneralizedStress_MembraneRot (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
virtual void giveGeneralizedStress_PlateSubSoil (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
virtual void giveCharMaterialStiffnessMatrix (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the stiffness matrix of receiver in given integration point, respecting its history. More...
 
virtual bool isCharacteristicMtrxSymmetric (MatResponseMode mode)
 Check for symmetry of stiffness matrix. More...
 
virtual void give3dDegeneratedShellStiffMtrx (FloatMatrix &answer, MatResponseMode rMode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 3d shell stiffness matrix on degenerated shell elements. More...
 
virtual void give2dBeamStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the stiffness matrix for 2d beams. More...
 
virtual void give3dBeamStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the stiffness matrix for 2d beams. More...
 
virtual void give2dPlateStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 2d plate stiffness matrix. More...
 
virtual void give3dShellStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing 3d shell stiffness matrix. More...
 
virtual void giveMembraneRotStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing membrane stiffness matrix with added drilling stiffness. More...
 
virtual void give2dPlateSubSoilStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing subsoil stiffness matrix for plates. More...
 
virtual IRResultType initializeFrom (InputRecord *ir)
 Initializes receiver acording to object description stored in input record. More...
 
virtual void giveInputRecord (DynamicInputRecord &input)
 Setups the input record string of receiver. More...
 
virtual void createMaterialStatus (GaussPoint &iGP)
 
virtual const char * giveClassName () const
 
virtual const char * giveInputRecordName () const
 
virtual double give (int aProperty, GaussPoint *gp)
 Returns the value of cross section property. More...
 
virtual double give (CrossSectionProperty a, GaussPoint *gp)
 Returns the value of cross section property at given point. More...
 
virtual double give (CrossSectionProperty a, const FloatArray &coords, Element *elem, bool local)
 Returns the value of cross section property at given point (belonging to given element). More...
 
virtual int giveIPValue (FloatArray &answer, GaussPoint *ip, InternalStateType type, TimeStep *tStep)
 Returns the integration point corresponding value in Reduced form. More...
 
virtual MaterialgiveMaterial (IntegrationPoint *ip)
 Returns the material associated with the GP. More...
 
int giveMaterialNumber () const
 
void setMaterialNumber (int matNum)
 
virtual int checkConsistency ()
 Allows programmer to test some internal data, before computation begins. More...
 
virtual InterfacegiveMaterialInterface (InterfaceType t, IntegrationPoint *ip)
 
virtual void giveFirstPKStresses (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedFIncrement, TimeStep *tStep)
 Computes the First Piola-Kirchoff stress vector for a given deformation gradient and integration point. More...
 
virtual void giveCauchyStresses (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedFIncrement, TimeStep *tStep)
 Computes the Cauchy stress vector for a given increment of deformation gradient and given integration point. More...
 
virtual void giveEshelbyStresses (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedvF, TimeStep *tStep)
 Computes the Eshelby stress vector. More...
 
virtual void giveStiffnessMatrix_dPdF (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the material stiffness matrix dPdF of receiver in a given integration point, respecting its history. More...
 
virtual void giveStiffnessMatrix_dCde (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Computes the material stiffness matrix dCde of receiver in a given integration point, respecting its history. More...
 
virtual void giveTemperatureVector (FloatArray &answer, GaussPoint *gp, TimeStep *tStep)
 
virtual int packUnknowns (DataStream &buff, TimeStep *tStep, GaussPoint *gp)
 Pack all necessary data of integration point (according to element parallel_mode) into given communication buffer. More...
 
virtual int unpackAndUpdateUnknowns (DataStream &buff, TimeStep *tStep, GaussPoint *gp)
 Unpack and updates all necessary data of given integration point (according to element parallel_mode) into given communication buffer. More...
 
virtual int estimatePackSize (DataStream &buff, GaussPoint *gp)
 Estimates the necessary pack size to hold all packed data of receiver. More...
 
- Public Member Functions inherited from oofem::StructuralCrossSection
 StructuralCrossSection (int n, Domain *d)
 Constructor. More...
 
virtual ~StructuralCrossSection ()
 Destructor. More...
 
virtual void give3dBeamSubSoilStiffMtrx (FloatMatrix &answer, MatResponseMode mode, GaussPoint *gp, TimeStep *tStep)
 Method for computing subsoil stiffness matrix for 2d beams. More...
 
virtual FloatArrayimposeStressConstrainsOnGradient (GaussPoint *gp, FloatArray *gradientStressVector3d)
 Returns modified gradient of stress vector, which is used to bring stresses back to yield surface. More...
 
virtual FloatArrayimposeStrainConstrainsOnGradient (GaussPoint *gp, FloatArray *gradientStressVector3d)
 Returns modified gradient of strain vector, which is used to compute plastic strain increment. More...
 
virtual int testCrossSectionExtension (CrossSectExtension ext)
 Returns nonzero, if receiver implements required extension. More...
 
void giveRealStresses (FloatArray &answer, GaussPoint *gp, const FloatArray &reducedStrain, TimeStep *tStep)
 Computes the real stress vector for given strain and integration point. More...
 
virtual void giveGeneralizedStress_3dBeamSubSoil (FloatArray &answer, GaussPoint *gp, const FloatArray &generalizedStrain, TimeStep *tStep)
 
- Public Member Functions inherited from oofem::CrossSection
 CrossSection (int n, Domain *d)
 Constructor. More...
 
virtual ~CrossSection ()
 Destructor. More...
 
int giveSetNumber () const
 
virtual bool hasProperty (CrossSectionProperty a)
 Returns true if the dictionary contains the requested property. More...
 
virtual void printYourself ()
 Prints receiver state on stdout. Useful for debugging. More...
 
virtual int setupIntegrationPoints (IntegrationRule &irule, int npoints, Element *element)
 Sets up integration rule for the given element. More...
 
virtual int setupIntegrationPoints (IntegrationRule &irule, int npointsXY, int npointsZ, Element *element)
 Sets up integration rule for the given element. More...
 
virtual double predictRelativeComputationalCost (GaussPoint *ip)
 Returns the weight representing relative computational cost of receiver The reference cross section is integral model in plane stress. More...
 
virtual double giveRelativeSelfComputationalCost ()
 Returns the weight representing relative computational cost of receiver The reference element is integral model in plane stress. More...
 
virtual double predictRelativeRedistributionCost (GaussPoint *gp)
 
virtual contextIOResultType saveIPContext (DataStream &stream, ContextMode mode, GaussPoint *gp)
 Stores integration point state to output stream. More...
 
virtual contextIOResultType restoreIPContext (DataStream &stream, ContextMode mode, GaussPoint *gp)
 Reads integration point state to output stream. 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 updateLocalNumbering (EntityRenumberingFunctor &f)
 Local renumbering support. More...
 
virtual contextIOResultType saveContext (DataStream &stream, ContextMode mode, void *obj=NULL)
 Stores receiver state to output stream. More...
 
virtual contextIOResultType restoreContext (DataStream &stream, ContextMode mode, void *obj=NULL)
 Restores the receiver state previously written in stream. More...
 
virtual void printOutputAt (FILE *file, TimeStep *tStep)
 Prints output of receiver to stream, for given time step. 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 Attributes

int materialNumber
 
int czMaterialNumber
 
- Protected Attributes inherited from oofem::CrossSection
Dictionary propertyDictionary
 Dictionary for storing cross section parameters (like dimensions). More...
 
int setNumber
 
- 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

Class implementing "simple" cross section model in finite element problem.

A cross section is an attribute of a domain. It is usually also attribute of many elements.

The simple cross section implementation does not perform any integration over cross-section volume, it represents a cross section model, where the whole cross section model is represented by single integration point. and therefore all requests for characteristic contributions (stiffness) and for real stress computations are simply passed to parent StructuralCrossSection class, which invokes corresponding material mode services. Please note, that it is assumed that material model will support these material modes and provide corresponding services for characteristic components and stress evaluation. For description, how to incorporate more elaborate models of cross section, please read base CrossSection documentation.

The overloaded methods giveFullCharacteristicVector and giveFullCharacteristicVector add some additional support for integrated cross section models - _3dShell, _3dBeam, _2dPlate and _2dBeam.

This class also reads into its property dictionary necessary geometric cross section characteristics, which are requested by particular material models. These parameters can be requested using get service and include those defined by CrossSectionProperty.

Definition at line 86 of file simplecrosssection.h.

Constructor & Destructor Documentation

oofem::SimpleCrossSection::SimpleCrossSection ( int  n,
Domain d 
)
inline

Constructor.

Parameters
nCross section number.
dAssociated domain.

Definition at line 94 of file simplecrosssection.h.

Member Function Documentation

int oofem::SimpleCrossSection::checkConsistency ( )
virtual

Allows programmer to test some internal data, before computation begins.

For example, one may use this function, to ensure that element has material with required capabilities is assigned to element. This must be done after all mesh components are instanciated.

Returns
Nonzero if receiver is consistent.

Implements oofem::StructuralCrossSection.

Definition at line 665 of file simplecrosssection.C.

References oofem::FEMComponent::giveClassName(), oofem::FEMComponent::giveDomain(), oofem::Domain::giveMaterial(), materialNumber, and OOFEM_WARNING.

void oofem::SimpleCrossSection::createMaterialStatus ( GaussPoint iGP)
virtual
int oofem::SimpleCrossSection::estimatePackSize ( DataStream buff,
GaussPoint ip 
)
virtual

Estimates the necessary pack size to hold all packed data of receiver.

The corresponding material model service is invoked. The nature of packed data is typically material model dependent.

Parameters
buffCommunication buffer.
ipIntegration point.
Returns
Estimate of pack size.

Implements oofem::CrossSection.

Definition at line 843 of file simplecrosssection.C.

References oofem::Material::estimatePackSize(), and giveMaterial().

double oofem::SimpleCrossSection::give ( int  aProperty,
GaussPoint gp 
)
virtual
virtual double oofem::SimpleCrossSection::give ( CrossSectionProperty  a,
GaussPoint gp 
)
inlinevirtual

Returns the value of cross section property at given point.

The default implementation assumes constant properties stored in propertyDictionary.

Parameters
aId of requested property.
gpIntegration point
Returns
Property value.

Reimplemented from oofem::CrossSection.

Reimplemented in oofem::VariableCrossSection.

Definition at line 154 of file simplecrosssection.h.

References oofem::CrossSection::give().

virtual double oofem::SimpleCrossSection::give ( CrossSectionProperty  a,
const FloatArray coords,
Element elem,
bool  local 
)
inlinevirtual

Returns the value of cross section property at given point (belonging to given element).

the point coordinates can be specified using its local element coordinates or global coordinates (one of these two can be set to NULL) The default implementation assumes constant properties stored in propertyDictionary.

Parameters
aId of requested property.
coordslocal or global coordinates (determined by local parameter) of point of interest
elemreference to underlying element containing given point
gpIntegration point
Returns
Property value.

Reimplemented from oofem::CrossSection.

Reimplemented in oofem::VariableCrossSection.

Definition at line 155 of file simplecrosssection.h.

References oofem::CrossSection::give().

void oofem::SimpleCrossSection::give2dBeamStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes the stiffness matrix for 2d beams.

Parameters
answerThe requested matrix.
modeMaterial response mode.
gpIntegration point.
tStepTime step.

Implements oofem::StructuralCrossSection.

Definition at line 314 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::CS_Area, oofem::CS_InertiaMomentY, oofem::CS_ShearAreaZ, oofem::Material::give(), give(), oofem::StructuralMaterial::give1dStressStiffMtrx(), giveMaterial(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by giveCharMaterialStiffnessMatrix(), and giveGeneralizedStress_Beam2d().

void oofem::SimpleCrossSection::give2dPlateStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing 2d plate stiffness matrix.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 372 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::CS_Thickness, give(), giveMaterial(), oofem::StructuralMaterial::givePlaneStressStiffMtrx(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by giveCharMaterialStiffnessMatrix(), and giveGeneralizedStress_Plate().

void oofem::SimpleCrossSection::give2dPlateSubSoilStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing subsoil stiffness matrix for plates.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 457 of file simplecrosssection.C.

References oofem::StructuralMaterial::give2dPlateSubSoilStiffMtrx(), and giveMaterial().

void oofem::SimpleCrossSection::give3dBeamStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes the stiffness matrix for 2d beams.

Parameters
answerThe requested matrix.
modeMaterial response mode.
gpIntegration point.
tStepTime step.
Todo:
Do this by using the general 3d tangent matrix instead somehow!!!

Implements oofem::StructuralCrossSection.

Definition at line 336 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::CS_Area, oofem::CS_InertiaMomentY, oofem::CS_InertiaMomentZ, oofem::CS_ShearAreaY, oofem::CS_ShearAreaZ, oofem::CS_TorsionMomentX, E, oofem::Material::give(), give(), oofem::StructuralMaterial::give1dStressStiffMtrx(), giveMaterial(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by giveCharMaterialStiffnessMatrix(), and giveGeneralizedStress_Beam3d().

void oofem::SimpleCrossSection::give3dDegeneratedShellStiffMtrx ( FloatMatrix answer,
MatResponseMode  rMode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing 3d shell stiffness matrix on degenerated shell elements.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::StructuralCrossSection.

Definition at line 429 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::StructuralMaterial::give3dMaterialStiffnessMatrix(), and giveMaterial().

Referenced by giveCharMaterialStiffnessMatrix().

void oofem::SimpleCrossSection::give3dShellStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing 3d shell stiffness matrix.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 399 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::CS_Thickness, give(), giveMaterial(), oofem::StructuralMaterial::givePlaneStressStiffMtrx(), oofem::FloatMatrix::resize(), and oofem::FloatMatrix::zero().

Referenced by giveCharMaterialStiffnessMatrix(), and giveGeneralizedStress_Shell().

void oofem::SimpleCrossSection::giveCauchyStresses ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedFIncrement,
TimeStep tStep 
)
virtual

Computes the Cauchy stress vector for a given increment of deformation gradient and given integration point.

The service should use previously reached equilibrium history variables. Also it should update temporary history variables in status according to newly reached state. The temporary history variables are moved into equilibrium ones after global structure equilibrium has been reached by iteration process. Elements should always pass their requests to their cross section model, which performs necessary integration over its volume and invokes necessary material services for corresponding material model defined for given integration point.

Parameters
answerContains the Cauchy stress.
gpIntegration point.
reducedFIncrementIncrement of the deformation gradient vector in reduced form.
Todo:
should this then be in a multiplicative way? /JB
Parameters
tStepCurrent time step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 717 of file simplecrosssection.C.

References oofem::StructuralMaterial::giveCauchyStressVector_1d(), oofem::StructuralMaterial::giveCauchyStressVector_3d(), oofem::StructuralMaterial::giveCauchyStressVector_PlaneStrain(), oofem::StructuralMaterial::giveCauchyStressVector_PlaneStress(), giveMaterial(), and oofem::GaussPoint::giveMaterialMode().

void oofem::SimpleCrossSection::giveCharMaterialStiffnessMatrix ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes the stiffness matrix of receiver in given integration point, respecting its history.

The algorithm should use temporary or equilibrium history variables stored in integration point status to compute and return required result. Elements should always pass their requests to their cross section model, which performs necessary integration over its volume and invokes necessary material services for corresponding material model defined for given integration point.

Parameters
answerContains result.
modeMaterial response mode.
gpIntegration point.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 282 of file simplecrosssection.C.

References oofem::StructuralMaterial::give1dStressStiffMtrx(), give2dBeamStiffMtrx(), give2dPlateStiffMtrx(), give3dBeamStiffMtrx(), give3dDegeneratedShellStiffMtrx(), oofem::StructuralMaterial::give3dMaterialStiffnessMatrix(), give3dShellStiffMtrx(), giveMaterial(), oofem::GaussPoint::giveMaterialMode(), oofem::StructuralMaterial::givePlaneStrainStiffMtrx(), oofem::StructuralMaterial::givePlaneStressStiffMtrx(), and oofem::StructuralMaterial::giveStiffnessMatrix().

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

Implements oofem::FEMComponent.

Reimplemented in oofem::VariableCrossSection, and oofem::WarpingCrossSection.

Definition at line 150 of file simplecrosssection.h.

void oofem::SimpleCrossSection::giveEshelbyStresses ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedvF,
TimeStep tStep 
)
virtual

Computes the Eshelby stress vector.

Does not update history variables, this is a postprocesing computation.

Parameters
answerContains the Eshelby stress.
gpIntegration point.
tStepCurrent time step (most models are able to respond only when tStep is current time step).

Reimplemented from oofem::StructuralCrossSection.

Definition at line 738 of file simplecrosssection.C.

References oofem::StructuralMaterial::giveEshelbyStressVector_PlaneStrain(), giveMaterial(), and oofem::GaussPoint::giveMaterialMode().

void oofem::SimpleCrossSection::giveFirstPKStresses ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedFIncrement,
TimeStep tStep 
)
virtual

Computes the First Piola-Kirchoff stress vector for a given deformation gradient and integration point.

The service should use previously reached equilibrium history variables. Also it should update temporary history variables in status according to newly reached state. The temporary history variables are moved into equilibrium ones after global structure equilibrium has been reached by iteration process. Elements should always pass their requests to their cross section model, which performs necessary integration over its volume and invokes necessary material services for corresponding material model defined for given integration point.

Parameters
answerContains the First Piola-Kirchoff stresses.
gpIntegration point.
reducedFIncrementIncrement of the deformation gradient vector in reduced form.
Todo:
should this then be in a multiplicative way? /JB
Parameters
tStepCurrent time step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 693 of file simplecrosssection.C.

References oofem::__MaterialModeToString(), oofem::StructuralMaterial::giveFirstPKStressVector_1d(), oofem::StructuralMaterial::giveFirstPKStressVector_3d(), oofem::StructuralMaterial::giveFirstPKStressVector_PlaneStrain(), oofem::StructuralMaterial::giveFirstPKStressVector_PlaneStress(), giveMaterial(), oofem::GaussPoint::giveMaterialMode(), and OOFEM_ERROR.

void oofem::SimpleCrossSection::giveGeneralizedStress_Beam2d ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual

Computes the generalized stress vector for given strain and integration point.

Parameters
answerContains result.
gpIntegration point.
generalizedStrainStrain vector in reduced generalized form.
tStepCurrent time step (most models are able to respond only when tStep is current time step).

Note: (by bp): This assumes that the behaviour is elastic there exist a nuumber of nonlinear integral material models for beams defined directly in terms of integral forces and moments and corresponding deformations and curvatures. This would require to implement support at material model level. Mikael: That would not be a continuum material model, but it would highly depend on the cross-section shape, thus, it should be a special cross-section model instead. This cross-section assumes you can split the response into inertia moments and pure material response. This is only possible for a constant, elastic response (i.e. elastic). Therefore, this cross-section may only be allowed to give the elastic response.

Implements oofem::StructuralCrossSection.

Definition at line 132 of file simplecrosssection.C.

References oofem::FloatArray::at(), oofem::FloatArray::beProductOf(), oofem::CS_Thickness, give(), give2dBeamStiffMtrx(), giveMaterial(), oofem::StructuralMaterial::giveReferenceTemperature(), oofem::FloatArray::giveSize(), oofem::Material::giveStatus(), giveTemperatureVector(), oofem::StructuralMaterial::giveThermalDilatationVector(), oofem::StructuralMaterialStatus::letTempStrainVectorBe(), and oofem::StructuralMaterialStatus::letTempStressVectorBe().

void oofem::SimpleCrossSection::giveGeneralizedStress_Beam3d ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual

Note: (by bp): This assumes that the behaviour is elastic there exist a nuumber of nonlinear integral material models for beams defined directly in terms of integral forces and moments and corresponding deformations and curvatures. This would require to implement support at material model level. Mikael: See earlier response to comment

Implements oofem::StructuralCrossSection.

Definition at line 163 of file simplecrosssection.C.

References oofem::FloatArray::at(), oofem::FloatArray::beProductOf(), oofem::CS_Thickness, oofem::CS_Width, give(), give3dBeamStiffMtrx(), giveMaterial(), oofem::StructuralMaterial::giveReferenceTemperature(), oofem::FloatArray::giveSize(), oofem::Material::giveStatus(), giveTemperatureVector(), oofem::StructuralMaterial::giveThermalDilatationVector(), oofem::StructuralMaterialStatus::letTempStrainVectorBe(), and oofem::StructuralMaterialStatus::letTempStressVectorBe().

void oofem::SimpleCrossSection::giveGeneralizedStress_MembraneRot ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual
Todo:
We should support nonlinear behavior for the membrane part. In fact, should be even bundle the rotation part with the membrane? We gain nothing from this design anyway as the rotation field is always separate. Separate manual integration by the element would be an option.

Implements oofem::StructuralCrossSection.

Definition at line 258 of file simplecrosssection.C.

References oofem::FloatArray::beProductOf(), giveMaterial(), giveMembraneRotStiffMtrx(), oofem::Material::giveStatus(), oofem::StructuralMaterialStatus::letTempStrainVectorBe(), and oofem::StructuralMaterialStatus::letTempStressVectorBe().

void oofem::SimpleCrossSection::giveGeneralizedStress_Plate ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual

Note: (by bp): This assumes that the behaviour is elastic there exist a number of nonlinear integral material models for beams/plates/shells defined directly in terms of integral forces and moments and corresponding deformations and curvatures. This would require to implement support at material model level. Mikael: See earlier response to comment

Implements oofem::StructuralCrossSection.

Definition at line 196 of file simplecrosssection.C.

References oofem::FloatArray::at(), oofem::FloatArray::beProductOf(), oofem::CS_Thickness, give(), give2dPlateStiffMtrx(), giveMaterial(), oofem::FloatArray::giveSize(), oofem::Material::giveStatus(), giveTemperatureVector(), oofem::StructuralMaterial::giveThermalDilatationVector(), oofem::StructuralMaterialStatus::letTempStrainVectorBe(), and oofem::StructuralMaterialStatus::letTempStressVectorBe().

void oofem::SimpleCrossSection::giveGeneralizedStress_PlateSubSoil ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveGeneralizedStress_Shell ( FloatArray answer,
GaussPoint gp,
const FloatArray generalizedStrain,
TimeStep tStep 
)
virtual

Note: (by bp): This assumes that the behaviour is elastic there exist a nuumber of nonlinear integral material models for beams/plates/shells defined directly in terms of integral forces and moments and corresponding deformations and curvatures. This would require to implement support at material model level. Mikael: See earlier response to comment

Implements oofem::StructuralCrossSection.

Definition at line 226 of file simplecrosssection.C.

References oofem::FloatArray::at(), oofem::FloatArray::beProductOf(), oofem::CS_Thickness, give(), give3dShellStiffMtrx(), giveMaterial(), oofem::StructuralMaterial::giveReferenceTemperature(), oofem::FloatArray::giveSize(), oofem::Material::giveStatus(), giveTemperatureVector(), oofem::StructuralMaterial::giveThermalDilatationVector(), oofem::StructuralMaterialStatus::letTempStrainVectorBe(), and oofem::StructuralMaterialStatus::letTempStressVectorBe().

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

Implements oofem::FEMComponent.

Reimplemented in oofem::VariableCrossSection, and oofem::WarpingCrossSection.

Definition at line 151 of file simplecrosssection.h.

References _IFT_SimpleCrossSection_Name.

int oofem::SimpleCrossSection::giveIPValue ( FloatArray answer,
GaussPoint ip,
InternalStateType  type,
TimeStep tStep 
)
virtual

Returns the integration point corresponding value in Reduced form.

Parameters
answercontain corresponding ip value, zero sized if not available
ipIntegration point.
typeDetermines the type of internal variable.
tStepTime step.
Returns
Nonzero if o.k, zero otherwise.

Reimplemented from oofem::CrossSection.

Definition at line 652 of file simplecrosssection.C.

References oofem::FloatArray::at(), oofem::Material::giveIPValue(), giveMaterial(), oofem::FEMComponent::giveNumber(), and oofem::FloatArray::resize().

Material * oofem::SimpleCrossSection::giveMaterial ( IntegrationPoint ip)
virtual
Interface * oofem::SimpleCrossSection::giveMaterialInterface ( InterfaceType  t,
IntegrationPoint ip 
)
virtual

Reimplemented from oofem::StructuralCrossSection.

Definition at line 684 of file simplecrosssection.C.

References oofem::FEMComponent::giveInterface(), and giveMaterial().

int oofem::SimpleCrossSection::giveMaterialNumber ( ) const
inline

Definition at line 159 of file simplecrosssection.h.

Referenced by giveMaterial(), and isCharacteristicMtrxSymmetric().

void oofem::SimpleCrossSection::giveMembraneRotStiffMtrx ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing membrane stiffness matrix with added drilling stiffness.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 447 of file simplecrosssection.C.

References oofem::FloatMatrix::at(), oofem::CS_DrillingStiffness, give(), giveMaterial(), oofem::StructuralMaterial::givePlaneStressStiffMtrx(), and oofem::FloatMatrix::resizeWithData().

Referenced by giveGeneralizedStress_MembraneRot().

void oofem::SimpleCrossSection::giveRealStress_1d ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveRealStress_3d ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveRealStress_3dDegeneratedShell ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveRealStress_PlaneStrain ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveRealStress_PlaneStress ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveRealStress_Warping ( FloatArray answer,
GaussPoint gp,
const FloatArray reducedStrain,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveStiffnessMatrix_1d ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveStiffnessMatrix_3d ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Method for computing the stiffness matrix.

Parameters
answerStiffness matrix.
modeMaterial response mode.
gpIntegration point, which load history is used.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 100 of file simplecrosssection.C.

References oofem::StructuralMaterial::give3dMaterialStiffnessMatrix(), and giveMaterial().

void oofem::SimpleCrossSection::giveStiffnessMatrix_dCde ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes the material stiffness matrix dCde of receiver in a given integration point, respecting its history.

The algorithm should use temporary or equilibrium history variables stored in integration point status to compute and return required result. Elements should always pass their requests to their cross section model, which performs necessary integration over its volume and invokes necessary material services for corresponding material model defined for given integration point.

Parameters
answerContains result.
modeMaterial response mode.
gpIntegration point.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 778 of file simplecrosssection.C.

References oofem::__MaterialModeToString(), oofem::StructuralMaterial::give1dStressStiffMtrx_dCde(), oofem::StructuralMaterial::give3dMaterialStiffnessMatrix_dCde(), giveMaterial(), oofem::GaussPoint::giveMaterialMode(), oofem::StructuralMaterial::givePlaneStrainStiffMtrx_dCde(), oofem::StructuralMaterial::givePlaneStressStiffMtrx_dCde(), and OOFEM_ERROR.

void oofem::SimpleCrossSection::giveStiffnessMatrix_dPdF ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual

Computes the material stiffness matrix dPdF of receiver in a given integration point, respecting its history.

The algorithm should use temporary or equilibrium history variables stored in integration point status to compute and return required result. Elements should always pass their requests to their cross section model, which performs necessary integration over its volume and invokes necessary material services for corresponding material model defined for given integration point.

Parameters
answerContains result.
modeMaterial response mode.
gpIntegration point.
tStepTime step (most models are able to respond only when tStep is current time step).

Implements oofem::StructuralCrossSection.

Definition at line 756 of file simplecrosssection.C.

References oofem::__MaterialModeToString(), oofem::StructuralMaterial::give1dStressStiffMtrx_dPdF(), oofem::StructuralMaterial::give3dMaterialStiffnessMatrix_dPdF(), giveMaterial(), oofem::GaussPoint::giveMaterialMode(), oofem::StructuralMaterial::givePlaneStrainStiffMtrx_dPdF(), oofem::StructuralMaterial::givePlaneStressStiffMtrx_dPdF(), and OOFEM_ERROR.

void oofem::SimpleCrossSection::giveStiffnessMatrix_PlaneStrain ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual
void oofem::SimpleCrossSection::giveStiffnessMatrix_PlaneStress ( FloatMatrix answer,
MatResponseMode  mode,
GaussPoint gp,
TimeStep tStep 
)
virtual
IRResultType oofem::SimpleCrossSection::initializeFrom ( InputRecord ir)
virtual

Initializes receiver acording to object description stored in input record.

Calls CrossSection initializeFrom service and reads the values of

  • 'thick' thickness
  • 'width' width
  • 'area' area
  • 'iy' Moment of inertia around y
  • 'iz' Moment of inertia around z
  • 'ik' Torsion moment around x
  • 'beamshearcoeff' Beam shear coefficient
    Parameters
    irRecord to read off.

Reimplemented from oofem::CrossSection.

Reimplemented in oofem::VariableCrossSection, and oofem::WarpingCrossSection.

Definition at line 466 of file simplecrosssection.C.

References _IFT_SimpleCrossSection_area, _IFT_SimpleCrossSection_directorx, _IFT_SimpleCrossSection_directory, _IFT_SimpleCrossSection_directorz, _IFT_SimpleCrossSection_drillStiffness, _IFT_SimpleCrossSection_drillType, _IFT_SimpleCrossSection_ik, _IFT_SimpleCrossSection_iy, _IFT_SimpleCrossSection_iz, _IFT_SimpleCrossSection_MaterialNumber, _IFT_SimpleCrossSection_relDrillStiffness, _IFT_SimpleCrossSection_shearareay, _IFT_SimpleCrossSection_shearareaz, _IFT_SimpleCrossSection_shearcoeff, _IFT_SimpleCrossSection_thick, _IFT_SimpleCrossSection_width, oofem::Dictionary::add(), oofem::CS_Area, oofem::CS_BeamShearCoeff, oofem::CS_DirectorVectorX, oofem::CS_DirectorVectorY, oofem::CS_DirectorVectorZ, oofem::CS_DrillingStiffness, oofem::CS_DrillingType, oofem::CS_InertiaMomentY, oofem::CS_InertiaMomentZ, oofem::CS_RelDrillingStiffness, oofem::CS_ShearAreaY, oofem::CS_ShearAreaZ, oofem::CS_Thickness, oofem::CS_TorsionMomentX, oofem::CS_Width, oofem::InputRecord::hasField(), oofem::CrossSection::initializeFrom(), IR_GIVE_FIELD, IR_GIVE_OPTIONAL_FIELD, materialNumber, and oofem::CrossSection::propertyDictionary.

Referenced by oofem::WarpingCrossSection::initializeFrom().

bool oofem::SimpleCrossSection::isCharacteristicMtrxSymmetric ( MatResponseMode  rMode)
virtual

Check for symmetry of stiffness matrix.

Default implementation returns true. It can be moved to base Cross section class in the future.

Parameters
rModeResponse mode of material.
Returns
True if stiffness matrix of receiver is symmetric.

Implements oofem::StructuralCrossSection.

Definition at line 623 of file simplecrosssection.C.

References oofem::FEMComponent::domain, oofem::Domain::giveMaterial(), giveMaterialNumber(), and oofem::Material::isCharacteristicMtrxSymmetric().

int oofem::SimpleCrossSection::packUnknowns ( DataStream buff,
TimeStep tStep,
GaussPoint ip 
)
virtual

Pack all necessary data of integration point (according to element parallel_mode) into given communication buffer.

The corresponding material model service for particular integration point is invoked. The nature of packed data is material model dependent. Typically, for material of "local" response (response depends only on integration point local state) no data are exchanged. For "nonlocal" constitutive models the send/receive of local values which undergo averaging is performed between local and corresponding remote elements.

Parameters
buffCommunication buffer.
tStepSolution step.
ipIntegration point.
Returns
Nonzero if successful.

Implements oofem::CrossSection.

Definition at line 831 of file simplecrosssection.C.

References giveMaterial(), and oofem::Material::packUnknowns().

void oofem::SimpleCrossSection::setMaterialNumber ( int  matNum)
inline

Definition at line 160 of file simplecrosssection.h.

int oofem::SimpleCrossSection::unpackAndUpdateUnknowns ( DataStream buff,
TimeStep tStep,
GaussPoint ip 
)
virtual

Unpack and updates all necessary data of given integration point (according to element parallel_mode) into given communication buffer.

See also
packUnknowns service.
Parameters
buffCommunication buffer.
tStepSolution step.
ipIntegration point.
Returns
Nonzero if successful.

Implements oofem::CrossSection.

Definition at line 837 of file simplecrosssection.C.

References giveMaterial(), and oofem::Material::unpackAndUpdateUnknowns().

Member Data Documentation

int oofem::SimpleCrossSection::czMaterialNumber
protected

Definition at line 179 of file simplecrosssection.h.


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:41 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011