Description
ANCF laminated shell element with four nodes.
This class implements composite material elastic force formulations.
The node numbering is in ccw fashion as in the following scheme:
v ^ | D o-----+-----o C | | | --+-----+-----+----> u | | | A o-----+-----o B
#include <ChElementShellANCF.h>
Classes | |
class | Layer |
Definition of a layer. More... | |
Public Types | |
using | ShapeVector = ChMatrixNM< double, 1, 8 > |
Public Member Functions | |
virtual int | GetNnodes () override |
Get the number of nodes used by this element. | |
virtual int | GetNdofs () override |
Get the number of coordinates in the field used by the referenced nodes. | |
virtual int | GetNodeNdofs (int n) override |
Get the number of coordinates from the n-th node used by this element. | |
void | SetNodes (std::shared_ptr< ChNodeFEAxyzD > nodeA, std::shared_ptr< ChNodeFEAxyzD > nodeB, std::shared_ptr< ChNodeFEAxyzD > nodeC, std::shared_ptr< ChNodeFEAxyzD > nodeD) |
Specify the nodes of this element. | |
void | SetDimensions (double lenX, double lenY) |
Specify the element dimensions. | |
virtual std::shared_ptr< ChNodeFEAbase > | GetNodeN (int n) override |
Access the n-th node of this element. | |
std::shared_ptr< ChNodeFEAxyzD > | GetNodeA () const |
Get a handle to the first node of this element. | |
std::shared_ptr< ChNodeFEAxyzD > | GetNodeB () const |
Get a handle to the second node of this element. | |
std::shared_ptr< ChNodeFEAxyzD > | GetNodeC () const |
Get a handle to the third node of this element. | |
std::shared_ptr< ChNodeFEAxyzD > | GetNodeD () const |
Get a handle to the fourth node of this element. | |
void | AddLayer (double thickness, double theta, std::shared_ptr< ChMaterialShellANCF > material) |
Add a layer. More... | |
size_t | GetNumLayers () const |
Get the number of layers. | |
const Layer & | GetLayer (size_t i) const |
Get a handle to the specified layer. | |
void | SetGravityOn (bool val) |
Turn gravity on/off. | |
void | SetAlphaDamp (double a) |
Set the structural damping. | |
double | GetLengthX () const |
Get the element length in the X direction. | |
double | GetLengthY () const |
Get the element length in the Y direction. | |
double | GetThickness () |
Get the total thickness of the shell element. | |
void | ShapeFunctions (ShapeVector &N, double x, double y, double z) |
Fills the N shape function matrix. More... | |
void | ShapeFunctionsDerivativeX (ShapeVector &Nx, double x, double y, double z) |
Fills the Nx shape function derivative matrix with respect to X. More... | |
void | ShapeFunctionsDerivativeY (ShapeVector &Ny, double x, double y, double z) |
Fills the Ny shape function derivative matrix with respect to Y. More... | |
void | ShapeFunctionsDerivativeZ (ShapeVector &Nz, double x, double y, double z) |
Fills the Nz shape function derivative matrix with respect to Z. More... | |
ChStrainStress3D | EvaluateSectionStrainStress (const ChVector<> &loc, int layer_id) |
Return a struct with 6-component strain and stress vectors evaluated at a given quadrature point and layer number. More... | |
void | EvaluateDeflection (double &defVec) |
virtual void | GetStateBlock (ChVectorDynamic<> &mD) override |
Fill the D vector with the current field values at the nodes of the element, with proper ordering. More... | |
virtual void | ComputeKRMmatricesGlobal (ChMatrixRef H, double Kfactor, double Rfactor=0, double Mfactor=0) override |
Sets H as the stiffness matrix K, scaled by Kfactor. More... | |
virtual void | ComputeMmatrixGlobal (ChMatrixRef M) override |
Set M as the global mass matrix. | |
virtual void | ComputeNodalMass () override |
Add contribution of element inertia to total nodal masses. More... | |
virtual void | ComputeInternalForces (ChVectorDynamic<> &Fi) override |
Computes the internal forces. More... | |
virtual void | Update () override |
Update the state of this element. | |
virtual void | EvaluateSectionDisplacement (const double u, const double v, ChVector<> &u_displ, ChVector<> &u_rotaz) override |
Gets the xyz displacement of a point on the shell, and the rotation RxRyRz of section reference, at parametric coordinates 'u' and 'v'. More... | |
virtual void | EvaluateSectionFrame (const double u, const double v, ChVector<> &point, ChQuaternion<> &rot) override |
Gets the absolute xyz position of a point on the shell, and the absolute rotation of section reference, at parametric coordinates 'u' and 'v'. More... | |
virtual void | EvaluateSectionPoint (const double u, const double v, ChVector<> &point) override |
Gets the absolute xyz position of a point on the shell, at parametric coordinates 'u' and 'v'. More... | |
void | ComputeInternalJacobians (double Kfactor, double Rfactor) |
Compute Jacobians of the internal forces. More... | |
void | ComputeMassMatrix () |
Compute the mass matrix of the element. More... | |
void | ComputeGravityForce (const ChVector<> &g_acc) |
Compute the gravitational forces. | |
void | ShapeFunctionANSbilinearShell (ChMatrixNM< double, 1, 4 > &S_ANS, double x, double y) |
void | CalcStrainANSbilinearShell () |
void | Basis_M (ChMatrixNM< double, 6, 5 > &M, double x, double y, double z) |
double | Calc_detJ0 (double x, double y, double z) |
double | Calc_detJ0 (double x, double y, double z, ShapeVector &Nx, ShapeVector &Ny, ShapeVector &Nz, ChMatrixNM< double, 1, 3 > &Nx_d0, ChMatrixNM< double, 1, 3 > &Ny_d0, ChMatrixNM< double, 1, 3 > &Nz_d0) |
void | CalcCoordMatrix (ChMatrixNM< double, 8, 3 > &d) |
void | CalcCoordDerivMatrix (ChVectorN< double, 24 > &dt) |
virtual int | LoadableGet_ndof_x () override |
Gets the number of DOFs affected by this element (position part). | |
virtual int | LoadableGet_ndof_w () override |
Gets the number of DOFs affected by this element (velocity part). | |
virtual void | LoadableGetStateBlock_x (int block_offset, ChState &mD) override |
Gets all the DOFs packed in a single vector (position part). | |
virtual void | LoadableGetStateBlock_w (int block_offset, ChStateDelta &mD) override |
Gets all the DOFs packed in a single vector (velocity part). | |
virtual void | LoadableStateIncrement (const unsigned int off_x, ChState &x_new, const ChState &x, const unsigned int off_v, const ChStateDelta &Dv) override |
Increment all DOFs using a delta. | |
virtual int | Get_field_ncoords () override |
Number of coordinates in the interpolated field, ex=3 for a tetrahedron finite element or a cable, = 1 for a thermal problem, etc. | |
virtual int | GetSubBlocks () override |
Tell the number of DOFs blocks (ex. =1 for a body, =4 for a tetrahedron, etc.) | |
virtual unsigned int | GetSubBlockOffset (int nblock) override |
Get the offset of the i-th sub-block of DOFs in global vector. | |
virtual unsigned int | GetSubBlockSize (int nblock) override |
Get the size of the i-th sub-block of DOFs in global vector. | |
virtual void | EvaluateSectionVelNorm (double U, double V, ChVector<> &Result) override |
Virtual method to plot velocity field distribution. More... | |
virtual void | LoadableGetVariables (std::vector< ChVariables * > &mvars) override |
Get the pointers to the contained ChVariables, appending to the mvars vector. | |
virtual void | ComputeNF (const double U, const double V, ChVectorDynamic<> &Qi, double &detJ, const ChVectorDynamic<> &F, ChVectorDynamic<> *state_x, ChVectorDynamic<> *state_w) override |
Evaluate N'*F , where N is some type of shape function evaluated at U,V coordinates of the surface, each ranging in -1..+1 F is a load, N'*F is the resulting generalized load Returns also det[J] with J=[dx/du,..], that might be useful in gauss quadrature. More... | |
virtual void | ComputeNF (const double U, const double V, const double W, ChVectorDynamic<> &Qi, double &detJ, const ChVectorDynamic<> &F, ChVectorDynamic<> *state_x, ChVectorDynamic<> *state_w) override |
Evaluate N'*F , where N is some type of shape function evaluated at U,V,W coordinates of the volume, each ranging in -1..+1 F is a load, N'*F is the resulting generalized load Returns also det[J] with J=[dx/du,..], that might be useful in gauss quadrature. More... | |
virtual double | GetDensity () override |
This is needed so that it can be accessed by ChLoaderVolumeGravity. More... | |
virtual ChVector | ComputeNormal (const double U, const double V) override |
Gets the normal to the surface at the parametric coordinate U,V. More... | |
Public Member Functions inherited from chrono::fea::ChElementShell | |
virtual bool | IsTriangleShell () |
Return false if quadrilateral shell - hence u,v parametric coordinates assumed in -1..+1, return true if triangular shell - hence u,v are triangle natural coordinates assumed in 0..+1. | |
Public Member Functions inherited from chrono::fea::ChElementGeneric | |
ChKblockGeneric & | Kstiffness () |
Access the proxy to stiffness, for sparse solver. | |
virtual void | EleIntLoadResidual_F (ChVectorDynamic<> &R, const double c) override |
(This is a default (a bit unoptimal) book keeping so that in children classes you can avoid implementing this EleIntLoadResidual_F function, unless you need faster code) | |
virtual void | EleIntLoadResidual_Mv (ChVectorDynamic<> &R, const ChVectorDynamic<> &w, const double c) override |
(This is a default (VERY UNOPTIMAL) book keeping so that in children classes you can avoid implementing this EleIntLoadResidual_Mv function, unless you need faster code.) | |
virtual void | EleIntLoadResidual_F_gravity (ChVectorDynamic<> &R, const ChVector<> &G_acc, const double c) override |
(This is a default (VERY UNOPTIMAL) book keeping so that in children classes you can avoid implementing this EleIntLoadResidual_F_gravity function, unless you need faster code. More... | |
virtual void | ComputeGravityForces (ChVectorDynamic<> &Fg, const ChVector<> &G_acc) override |
(This is the default implementation (POTENTIALLY INEFFICIENT) so that in children classes you can avoid implementing this ComputeGravityForces() function, unless you need faster code.) This fallback implementation uses a temp ChLoaderGravity that applies the load to elements only if they are inherited by ChLoadableUVW so it can use GetDensity() and Gauss quadrature. | |
virtual void | InjectKRMmatrices (ChSystemDescriptor &mdescriptor) override |
Tell to a system descriptor that there are item(s) of type ChKblock in this object (for further passing it to a solver) | |
virtual void | KRMmatricesLoad (double Kfactor, double Rfactor, double Mfactor) override |
Adds the current stiffness K and damping R and mass M matrices in encapsulated ChKblock item(s), if any. More... | |
virtual void | VariablesFbLoadInternalForces (double factor=1.) override |
Adds the internal forces, expressed as nodal forces, into the encapsulated ChVariables, in the 'fb' part: qf+=forces*factor (This is a default (a bit unoptimal) book keeping so that in children classes you can avoid implementing this VariablesFbLoadInternalForces function, unless you need faster code) | |
virtual void | VariablesFbIncrementMq () override |
Adds M*q (internal masses multiplied current 'qb') to Fb, ex. More... | |
Public Member Functions inherited from chrono::fea::ChElementBase | |
virtual void | EleDoIntegration () |
This is optionally implemented if there is some internal state that requires integration. | |
Public Member Functions inherited from chrono::ChLoadableUV | |
virtual bool | IsTriangleIntegrationNeeded () |
If true, use quadrature over u,v in [0..1] range as triangle area coords (with z=1-u-v) otherwise use default quadrature over u,v in [-1..+1] as rectangular isoparametric coords. | |
Public Member Functions inherited from chrono::ChLoadableUVW | |
virtual bool | IsTetrahedronIntegrationNeeded () |
If true, use quadrature over u,v,w in [0..1] range as tetrahedron volumetric coords (with z=1-u-v-w) otherwise use default quadrature over u,v,w in [-1..+1] as box isoparametric coords. | |
virtual bool | IsTrianglePrismIntegrationNeeded () |
If true, use quadrature over u,v in [0..1] range as triangle natural coords (with z=1-u-v), and use linear quadrature over w in [-1..+1], otherwise use default quadrature over u,v,w in [-1..+1] as box isoparametric coords. | |
Friends | |
class | ShellANCF_Mass |
class | ShellANCF_Gravity |
class | ShellANCF_Force |
class | ShellANCF_Jacobian |
Additional Inherited Members | |
Protected Attributes inherited from chrono::fea::ChElementShell | |
double | mass |
Protected Attributes inherited from chrono::fea::ChElementGeneric | |
ChKblockGeneric | Kmatr |
Member Function Documentation
◆ AddLayer()
void chrono::fea::ChElementShellANCF::AddLayer | ( | double | thickness, |
double | theta, | ||
std::shared_ptr< ChMaterialShellANCF > | material | ||
) |
Add a layer.
- Parameters
-
thickness layer thickness theta fiber angle (radians) material layer material
◆ ComputeInternalForces()
|
overridevirtual |
Computes the internal forces.
(E.g. the actual position of nodes is not in relaxed reference position) and set values in the Fi vector.
Implements chrono::fea::ChElementBase.
◆ ComputeInternalJacobians()
void chrono::fea::ChElementShellANCF::ComputeInternalJacobians | ( | double | Kfactor, |
double | Rfactor | ||
) |
Compute Jacobians of the internal forces.
This function calculates a linear combination of the stiffness (K) and damping (R) matrices, J = Kfactor * K + Rfactor * R for given coefficients Kfactor and Rfactor. This Jacobian will be further combined with the global mass matrix M and included in the global stiffness matrix H in the function ComputeKRMmatricesGlobal().
◆ ComputeKRMmatricesGlobal()
|
overridevirtual |
Sets H as the stiffness matrix K, scaled by Kfactor.
Optionally, also superimposes global damping matrix R, scaled by Rfactor, and mass matrix M, scaled by Mfactor. Matrices are expressed in global reference. Corotational elements can take the local Kl & Rl matrices and rotate them.
Implements chrono::fea::ChElementBase.
◆ ComputeMassMatrix()
void chrono::fea::ChElementShellANCF::ComputeMassMatrix | ( | ) |
Compute the mass matrix of the element.
Note: in this 'basic' implementation, constant section and constant material are assumed
◆ ComputeNF() [1/2]
|
overridevirtual |
Evaluate N'*F , where N is some type of shape function evaluated at U,V coordinates of the surface, each ranging in -1..+1 F is a load, N'*F is the resulting generalized load Returns also det[J] with J=[dx/du,..], that might be useful in gauss quadrature.
- Parameters
-
U parametric coordinate in surface V parametric coordinate in surface Qi Return result of Q = N'*F here detJ Return det[J] here F Input F vector, size is =n. field coords. state_x if != 0, update state (pos. part) to this, then evaluate Q state_w if != 0, update state (speed part) to this, then evaluate Q
Implements chrono::ChLoadableUV.
◆ ComputeNF() [2/2]
|
overridevirtual |
Evaluate N'*F , where N is some type of shape function evaluated at U,V,W coordinates of the volume, each ranging in -1..+1 F is a load, N'*F is the resulting generalized load Returns also det[J] with J=[dx/du,..], that might be useful in gauss quadrature.
- Parameters
-
U parametric coordinate in volume V parametric coordinate in volume W parametric coordinate in volume Qi Return result of N'*F here, maybe with offset block_offset detJ Return det[J] here F Input F vector, size is = n.field coords. state_x if != 0, update state (pos. part) to this, then evaluate Q state_w if != 0, update state (speed part) to this, then evaluate Q
Implements chrono::ChLoadableUVW.
◆ ComputeNodalMass()
|
overridevirtual |
Add contribution of element inertia to total nodal masses.
This class computes and adds corresponding masses to ElementGeneric member m_TotalMass.
Reimplemented from chrono::fea::ChElementBase.
◆ ComputeNormal()
|
overridevirtual |
Gets the normal to the surface at the parametric coordinate U,V.
Each coordinate ranging in -1..+1.
Implements chrono::ChLoadableUV.
◆ EvaluateSectionDisplacement()
|
overridevirtual |
Gets the xyz displacement of a point on the shell, and the rotation RxRyRz of section reference, at parametric coordinates 'u' and 'v'.
Note, u=-1..+1 , v= -1..+1 parametric coordinates, except if triangular shell, where u=0..+1, v=0..+1, natural triangle coords. Results are not corotated.
Implements chrono::fea::ChElementShell.
◆ EvaluateSectionFrame()
|
overridevirtual |
Gets the absolute xyz position of a point on the shell, and the absolute rotation of section reference, at parametric coordinates 'u' and 'v'.
Note, u=-1..+1 , v= -1..+1 parametric coordinates, except if triangular shell, where u=0..+1, v=0..+1, natural triangle coords. Results are corotated.
Implements chrono::fea::ChElementShell.
◆ EvaluateSectionPoint()
|
overridevirtual |
Gets the absolute xyz position of a point on the shell, at parametric coordinates 'u' and 'v'.
Note, u=-1..+1 , v= -1..+1 parametric coordinates, except if triangular shell, where u=0..+1, v=0..+1, natural triangle coords. Results are corotated.
Implements chrono::fea::ChElementShell.
◆ EvaluateSectionStrainStress()
ChStrainStress3D chrono::fea::ChElementShellANCF::EvaluateSectionStrainStress | ( | const ChVector<> & | loc, |
int | layer_id | ||
) |
Return a struct with 6-component strain and stress vectors evaluated at a given quadrature point and layer number.
Beta
◆ EvaluateSectionVelNorm()
|
overridevirtual |
Virtual method to plot velocity field distribution.
Note, u=-1..+1 , v= -1..+1 parametric coordinates, except if triangular shell, where u=0..+1, v=0..+1, natural triangle coords.
Implements chrono::fea::ChElementShell.
◆ GetDensity()
|
overridevirtual |
This is needed so that it can be accessed by ChLoaderVolumeGravity.
Density is mass per unit surface.
Implements chrono::ChLoadableUVW.
◆ GetStateBlock()
|
overridevirtual |
Fill the D vector with the current field values at the nodes of the element, with proper ordering.
If the D vector has not the size of this->GetNdofs(), it will be resized. {x_a y_a z_a Dx_a Dx_a Dx_a x_b y_b z_b Dx_b Dy_b Dz_b}
Implements chrono::fea::ChElementBase.
◆ ShapeFunctions()
void chrono::fea::ChElementShellANCF::ShapeFunctions | ( | ShapeVector & | N, |
double | x, | ||
double | y, | ||
double | z | ||
) |
Fills the N shape function matrix.
NOTE! actually N should be a 3row, 24 column sparse matrix, as N = [s1*eye(3) s2*eye(3) s3*eye(3) s4*eye(3)...]; , but to avoid wasting zero and repeated elements, here it stores only the s1 through s8 values in a 1 row, 8 columns matrix!
◆ ShapeFunctionsDerivativeX()
void chrono::fea::ChElementShellANCF::ShapeFunctionsDerivativeX | ( | ShapeVector & | Nx, |
double | x, | ||
double | y, | ||
double | z | ||
) |
Fills the Nx shape function derivative matrix with respect to X.
NOTE! to avoid wasting zero and repeated elements, here it stores only the four values in a 1 row, 8 columns matrix!
◆ ShapeFunctionsDerivativeY()
void chrono::fea::ChElementShellANCF::ShapeFunctionsDerivativeY | ( | ShapeVector & | Ny, |
double | x, | ||
double | y, | ||
double | z | ||
) |
Fills the Ny shape function derivative matrix with respect to Y.
NOTE! to avoid wasting zero and repeated elements, here it stores only the four values in a 1 row, 8 columns matrix!
◆ ShapeFunctionsDerivativeZ()
void chrono::fea::ChElementShellANCF::ShapeFunctionsDerivativeZ | ( | ShapeVector & | Nz, |
double | x, | ||
double | y, | ||
double | z | ||
) |
Fills the Nz shape function derivative matrix with respect to Z.
NOTE! to avoid wasting zero and repeated elements, here it stores only the four values in a 1 row, 8 columns matrix!
The documentation for this class was generated from the following files:
- /builds/uwsbel/chrono/src/chrono/fea/ChElementShellANCF.h
- /builds/uwsbel/chrono/src/chrono/fea/ChElementShellANCF.cpp