solid_element.hpp

Go to the documentation of this file.

//
//   Project Name:        KratosDamApplication $
//   Created by:          $Author:            JMCarbonell $
//   Last modified by:    $Co-Author:                     $
//   Date:                $Date:                July 2013 $
//   Revision:            $Revision:                  0.0 $
//
//
 
#if !defined(KRATOS_SOLID_ELEMENT_H_INCLUDED)
#define  KRATOS_SOLID_ELEMENT_H_INCLUDED
 
// System includes
 
// External includes
 
// Project includes
#include "includes/kratos_flags.h"
#include "includes/checks.h"
#include "includes/element.h"
#include "custom_utilities/element_utilities.hpp"
 
 
namespace Kratos
{
 
 
class KRATOS_API(DAM_APPLICATION) SolidElement : public Element
{
    public:
 
    typedef ConstitutiveLaw ConstitutiveLawType;
    typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
    typedef ConstitutiveLawType::StressMeasure StressMeasureType;
    typedef GeometryData::IntegrationMethod IntegrationMethod;
    typedef GeometryData::SizeType SizeType;
 
    KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION( SolidElement );
 
    protected:
 
    KRATOS_DEFINE_LOCAL_FLAG( COMPUTE_RHS_VECTOR );
    KRATOS_DEFINE_LOCAL_FLAG( COMPUTE_LHS_MATRIX );
    KRATOS_DEFINE_LOCAL_FLAG( FINALIZED_STEP );
 
    struct ElementData
    {
        private:
 
        //variables including all integration points
        const GeometryType::ShapeFunctionsGradientsType* pDN_De;
        const Matrix* pNcontainer;
        const ProcessInfo* pProcessInfo;
 
        public:
 
        StressMeasureType StressMeasure;
 
        double  Tau;
        double  IntegrationWeight;
 
        //for axisymmetric use only
        double  CurrentRadius;
        double  ReferenceRadius;
 
        //general variables for large displacement use
        double  detF;
        double  detF0;
        double  detH; //Wildcard ( detF(0 to n+1) )
        double  detJ;
        Vector  StrainVector;
        Vector  StressVector;
        Vector  N;
        Matrix  B;
        Matrix  H;    //Wildcard ( Displacement Gradient, F(0 to n+1), B-bar, Velocity Gradient...)
        Matrix  F;    //Incremental Deformation Gradient (n to n+1)
        Matrix  F0;   //Historical Deformation Gradient  (0 to n)
        Matrix  DN_DX;
        Matrix  ConstitutiveMatrix;
 
        //variables including all integration points
        GeometryType::JacobiansType J;
        GeometryType::JacobiansType j;
        Matrix  DeltaPosition;
 
 
        void SetShapeFunctionsGradients(const GeometryType::ShapeFunctionsGradientsType &rDN_De)
        {
            pDN_De=&rDN_De;
        };
 
        void SetShapeFunctions(const Matrix& rNcontainer)
        {
            pNcontainer=&rNcontainer;
        };
 
        void SetProcessInfo(const ProcessInfo& rProcessInfo)
        {
            pProcessInfo=&rProcessInfo;
        };
 
 
 
        const GeometryType::ShapeFunctionsGradientsType& GetShapeFunctionsGradients()
        {
            return *pDN_De;
        };
 
        const Matrix& GetShapeFunctions()
        {
            return *pNcontainer;
        };
 
        const ProcessInfo& GetProcessInfo()
        {
            return *pProcessInfo;
        };
 
        void Initialize(const unsigned int& voigt_size,
                              const unsigned int& dimension,
                              const unsigned int& number_of_nodes)
        {
              StressMeasure = ConstitutiveLaw::StressMeasure_PK2;
 
            //stabilization
            Tau = 0;
 
            //time step
            IntegrationWeight = 1;
 
            //radius
              CurrentRadius = 0;
              ReferenceRadius = 0;
 
            //jacobians
              detF  = 1;
              detF0 = 1;
              detH  = 1;
              detJ  = 1;
 
            //vectors
              StrainVector.resize(voigt_size,false);
            StressVector.resize(voigt_size,false);
              N.resize(number_of_nodes,false);
              noalias(StrainVector) = ZeroVector(voigt_size);
              noalias(StressVector) = ZeroVector(voigt_size);
              noalias(N) = ZeroVector(number_of_nodes);
 
            //matrices
              B.resize(voigt_size, dimension*number_of_nodes,false);
              H.resize(dimension,dimension,false);
              F.resize(dimension,dimension,false);
              F0.resize(dimension,dimension,false);
              DN_DX.resize(number_of_nodes, dimension,false);
              ConstitutiveMatrix.resize(voigt_size, voigt_size,false);
              DeltaPosition.resize(number_of_nodes, dimension,false);
 
              noalias(B)  = ZeroMatrix(voigt_size, dimension*number_of_nodes);
              noalias(H)  = IdentityMatrix(dimension);
              noalias(F)  = IdentityMatrix(dimension);
              noalias(F0) = IdentityMatrix(dimension);
              noalias(DN_DX) = ZeroMatrix(number_of_nodes, dimension);
              noalias(ConstitutiveMatrix) = ZeroMatrix(voigt_size, voigt_size);
              noalias(DeltaPosition) = ZeroMatrix(number_of_nodes, dimension);
 
            //others
              J.resize(1,false);
              j.resize(1,false);
              J[0].resize(dimension,dimension,false);
              j[0].resize(dimension,dimension,false);
              noalias(J[0]) = ZeroMatrix(dimension,dimension);
              noalias(j[0]) = ZeroMatrix(dimension,dimension);
 
            //pointers
              pDN_De = NULL;
              pNcontainer = NULL;
            pProcessInfo = NULL;
          }
 
    };
 
 
    struct LocalSystemComponents
    {
        private:
 
        //for calculation local system with compacted LHS and RHS
        MatrixType *mpLeftHandSideMatrix;
        VectorType *mpRightHandSideVector;
 
        public:
 
        //calculation flags
        Flags  CalculationFlags;
 
        void SetLeftHandSideMatrix( MatrixType& rLeftHandSideMatrix ) { mpLeftHandSideMatrix = &rLeftHandSideMatrix; };
 
        void SetRightHandSideVector( VectorType& rRightHandSideVector ) { mpRightHandSideVector = &rRightHandSideVector; };
 
 
        MatrixType& GetLeftHandSideMatrix() { return *mpLeftHandSideMatrix; };
 
        VectorType& GetRightHandSideVector() { return *mpRightHandSideVector; };
 
    };
 
 
    public:
 
 
    typedef ElementData ElementDataType;
 
 
    SolidElement();
 
    SolidElement(IndexType NewId, GeometryType::Pointer pGeometry);
 
    SolidElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties);
 
    SolidElement(SolidElement const& rOther);
 
    ~SolidElement() override;
 
 
    SolidElement& operator=(SolidElement const& rOther);
 
 
    Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, PropertiesType::Pointer pProperties) const override;
 
    Element::Pointer Clone(IndexType NewId, NodesArrayType const& ThisNodes) const override;
 
 
    //************* GETTING METHODS
 
    IntegrationMethod GetIntegrationMethod() const override;
 
    void GetDofList(DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
 
    void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
 
    void GetValuesVector(Vector& rValues, int Step = 0) const override;
 
    void GetFirstDerivativesVector(Vector& rValues, int Step = 0) const override;
 
    void GetSecondDerivativesVector(Vector& rValues, int Step = 0) const override;
 
 
 
    //on integration points:
    //SET
    void SetValuesOnIntegrationPoints(const Variable<double>& rVariable, const std::vector<double>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
 
    void SetValuesOnIntegrationPoints(const Variable<Vector>& rVariable, const std::vector<Vector>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
 
    void SetValuesOnIntegrationPoints(const Variable<Matrix>& rVariable, const std::vector<Matrix>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
 
    void SetValuesOnIntegrationPoints(const Variable<ConstitutiveLaw::Pointer>& rVariable,
                                      const std::vector<ConstitutiveLaw::Pointer>& rValues,
                                      const ProcessInfo& rCurrentProcessInfo ) override;
 
 
    //************* STARTING - ENDING  METHODS
 
    void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
 
    void InitializeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
 
    void InitializeNonLinearIteration(const ProcessInfo& rCurrentProcessInfo) override;
 
    void FinalizeNonLinearIteration(const ProcessInfo& rCurrentProcessInfo) override;
 
    void FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
 
 
    //************* COMPUTING  METHODS
 
 
    void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,
                              VectorType& rRightHandSideVector,
                              const ProcessInfo& rCurrentProcessInfo) override;
 
 
    void CalculateRightHandSide(VectorType& rRightHandSideVector,
                                const ProcessInfo& rCurrentProcessInfo) override;
 
 
    void CalculateLeftHandSide (MatrixType& rLeftHandSideMatrix,
                                const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateFirstDerivativesContributions(MatrixType& rLeftHandSideMatrix,
                                                VectorType& rRightHandSideVector,
                                                const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateSecondDerivativesContributions(MatrixType& rLeftHandSideMatrix,
                                                 VectorType& rRightHandSideVector,
                                                 const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateSecondDerivativesLHS(MatrixType& rLeftHandSideMatrix,
                                       const ProcessInfo& rCurrentProcessInfo) override;
 
 
    void CalculateSecondDerivativesRHS(VectorType& rRightHandSideVector,
                                       const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateMassMatrix(MatrixType& rMassMatrix,
                             const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateDampingMatrix(MatrixType& rDampingMatrix,
                                const ProcessInfo& rCurrentProcessInfo) override;
 
 
    void AddExplicitContribution(const VectorType& rRHSVector,
                                 const Variable<VectorType>& rRHSVariable,
                                 const Variable<array_1d<double,3> >& rDestinationVariable,
                                 const ProcessInfo& rCurrentProcessInfo) override;
 
    //on integration points:
    void CalculateOnIntegrationPoints(const Variable<double>& rVariable,
                                      std::vector<double>& rOutput,
                                      const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateOnIntegrationPoints(const Variable<Vector>& rVariable,
                                      std::vector<Vector>& rOutput,
                                      const ProcessInfo& rCurrentProcessInfo) override;
 
    void CalculateOnIntegrationPoints(const Variable<Matrix >& rVariable,
                                      std::vector< Matrix >& rOutput,
                                      const ProcessInfo& rCurrentProcessInfo) override;
 
 
    //************************************************************************************
    //************************************************************************************
    int Check(const ProcessInfo& rCurrentProcessInfo) const override;
 
 
 
    std::string Info() const override
    {
        std::stringstream buffer;
        buffer << "Large Displacement Element #" << Id();
        return buffer.str();
    }
 
    void PrintInfo(std::ostream& rOStream) const override
    {
        rOStream << "Large Displacement Element #" << Id();
    }
 
    void PrintData(std::ostream& rOStream) const override
    {
        GetGeometry().PrintData(rOStream);
    }
 
    protected:
 
 
    IntegrationMethod mThisIntegrationMethod;
 
    std::vector<ConstitutiveLaw::Pointer> mConstitutiveLawVector;
 
 
 
    //constexpr const std::size_t& Dimension() const {return GetGeometry().WorkingSpaceDimension();}
 
 
    void IncreaseIntegrationMethod(IntegrationMethod& rThisIntegrationMethod,
                   unsigned int increment) const;
 
    virtual void CalculateElementalSystem(LocalSystemComponents& rLocalSystem,
                                          const ProcessInfo& rCurrentProcessInfo);
 
    virtual void CalculateDynamicSystem(LocalSystemComponents& rLocalSystem,
                                        const ProcessInfo& rCurrentProcessInfo);
 
    void PrintElementCalculation(LocalSystemComponents& rLocalSystem, ElementDataType& rVariables);
 
 
    void CalculatePerturbedLeftHandSide (MatrixType& rLeftHandSideMatrix,
                                         const ProcessInfo& rCurrentProcessInfo);
 
    virtual void CalculateAndAddLHS(LocalSystemComponents& rLocalSystem,
                                    ElementDataType& rVariables,
                                    double& rIntegrationWeight);
 
    virtual void CalculateAndAddRHS(LocalSystemComponents& rLocalSystem,
                                    ElementDataType& rVariables,
                                    Vector& rVolumeForce,
                                    double& rIntegrationWeight);
 
    virtual void CalculateAndAddDynamicLHS(MatrixType& rLeftHandSideMatrix,
                       ElementDataType& rVariables,
                       const ProcessInfo& rCurrentProcessInfo,
                       double& rIntegrationWeight);
 
    virtual void CalculateAndAddDynamicRHS(VectorType& rRightHandSideVector,
                       ElementDataType& rVariables,
                       const ProcessInfo& rCurrentProcessInfo,
                       double& rIntegrationWeight);
 
    virtual void CalculateAndAddKuum(MatrixType& rLeftHandSideMatrix,
                                     ElementDataType& rVariables,
                                     double& rIntegrationWeight);
 
    virtual void CalculateAndAddKuug(MatrixType& rLeftHandSideMatrix,
                                     ElementDataType& rVariables,
                                     double& rIntegrationWeight);
 
    virtual void CalculateAndAddExternalForces(VectorType& rRightHandSideVector,
                           ElementDataType& rVariables,
                           Vector& rVolumeForce,
                           double& rIntegrationWeight);
 
    virtual void CalculateAndAddInternalForces(VectorType& rRightHandSideVector,
                           ElementDataType & rVariables,
                           double& rIntegrationWeight);
 
 
    virtual void SetElementData(ElementDataType& rVariables,
                                ConstitutiveLaw::Parameters& rValues,
                                const int & rPointNumber);
 
    virtual void CalculateMaterialResponse(ElementDataType& rVariables,
                                           ConstitutiveLaw::Parameters& rValues,
                                           const int & rPointNumber);
    virtual SizeType GetDofsSize() const;
 
    bool IsSliver()
    {
      if( this->IsDefined(SELECTED) )
        return this->Is(SELECTED);
      else
        return false;
    }
 
    virtual void InitializeSystemMatrices(MatrixType& rLeftHandSideMatrix,
                                          VectorType& rRightHandSideVector,
                                          Flags& rCalculationFlags);
 
    void InitializeConstitutiveLaw();
 
    void ResetConstitutiveLaw() override;
 
    void InitializeExplicitContributions();
 
    virtual void CalculateKinematics(ElementDataType& rVariables,
                                     const double& rPointNumber);
 
    virtual void CalculateKinetics(ElementDataType& rVariables,
                   const double& rPointNumber);
 
    virtual void InitializeElementData(ElementDataType & rVariables,
                                       const ProcessInfo& rCurrentProcessInfo);
 
    virtual void TransformElementData(ElementDataType & rVariables,
                       const double& rPointNumber);
 
    virtual void FinalizeStepVariables(ElementDataType & rVariables,
                       const double& rPointNumber);
 
    virtual double& CalculateIntegrationWeight(double& rIntegrationWeight);
 
    virtual double& CalculateTotalMass(double& rTotalMass, const ProcessInfo& rCurrentProcessInfo);
 
    virtual Matrix& CalculateDeltaPosition(Matrix & rDeltaPosition);
 
    virtual Matrix& CalculateTotalDeltaPosition(Matrix & rDeltaPosition);
 
    virtual double& CalculateVolumeChange(double& rVolumeChange, ElementDataType& rVariables);
 
    virtual Vector& CalculateVolumeForce(Vector& rVolumeForce, ElementDataType& rVariables);
 
 
    private:
 
 
 
 
 
 
 
 
    friend class Serializer;
 
    // A private default constructor necessary for serialization
 
    void save(Serializer& rSerializer) const override;
 
    void load(Serializer& rSerializer) override;
 
 
 
}; // Class SolidElement
 
 
} // namespace Kratos.
#endif // KRATOS_SOLID_ELEMENT_H_INCLUDED  defined