KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
List of all members
Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace > Class Template Reference

Bossak integration scheme (for linear and nonlinear dynamic problems) for displacements adjusted for Material Point Method. More...

#include <mpm_residual_based_bossak_scheme.hpp>

Inheritance diagram for Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >:
Collaboration diagram for Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >:

Type Definitions

typedef Scheme< TSparseSpace, TDenseSpace > BaseType
 
typedef ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace > ImplicitBaseType
 
typedef ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace > BossakBaseType
 
typedef BossakBaseType::TDataType TDataType
 
typedef BossakBaseType::DofsArrayType DofsArrayType
 
typedef Element::DofsVectorType DofsVectorType
 
typedef BossakBaseType::TSystemMatrixType TSystemMatrixType
 
typedef BossakBaseType::TSystemVectorType TSystemVectorType
 
typedef BossakBaseType::LocalSystemVectorType LocalSystemVectorType
 
typedef BossakBaseType::LocalSystemMatrixType LocalSystemMatrixType
 
typedef ModelPart::ElementsContainerType ElementsArrayType
 
typedef ModelPart::ConditionsContainerType ConditionsArrayType
 
typedef BaseType::Pointer BaseTypePointer
 
ModelPartmGridModelPart
 
bool mIsDynamic
 
unsigned int mDomainSize
 
unsigned int mBlockSize
 
MPMBoundaryRotationUtility< LocalSystemMatrixType, LocalSystemVectorTypemRotationTool
 
 KRATOS_CLASS_POINTER_DEFINITION (MPMResidualBasedBossakScheme)
 
 MPMResidualBasedBossakScheme (ModelPart &rGridModelPart, unsigned int DomainSize, unsigned int BlockSize, double Alpha=0.0, double NewmarkBeta=0.25, bool IsDynamic=true)
 Constructor. @detail The MPM bossak method. More...
 
 MPMResidualBasedBossakScheme (MPMResidualBasedBossakScheme &rOther)
 Copy Constructor. More...
 
BaseTypePointer Clone () override
 Clone method. More...
 
virtual ~MPMResidualBasedBossakScheme ()
 
void Update (ModelPart &rModelPart, DofsArrayType &rDofSet, TSystemMatrixType &rA, TSystemVectorType &rDx, TSystemVectorType &rb) override
 Performing the update of the solution. More...
 
void Predict (ModelPart &rModelPart, DofsArrayType &rDofSet, TSystemMatrixType &rA, TSystemVectorType &rDx, TSystemVectorType &rb) override
 Performing the prediction of the solution. More...
 
void FinalizeNonLinIteration (ModelPart &rModelPart, TSystemMatrixType &rA, TSystemVectorType &rDx, TSystemVectorType &rb) override
 Function to be called when it is needed to finalize an iteration. It is designed to be called at the end of each non linear iteration. More...
 
void InitializeSolutionStep (ModelPart &rModelPart, TSystemMatrixType &rA, TSystemVectorType &rDx, TSystemVectorType &rb) override
 It initializes time step solution for MPM simulations. More...
 
void CalculateSystemContributions (Element &rCurrentElement, LocalSystemMatrixType &LHS_Contribution, LocalSystemVectorType &RHS_Contribution, Element::EquationIdVectorType &EquationId, const ProcessInfo &rCurrentProcessInfo) override
 This function is designed to be called in the builder and solver to introduce the selected time integration scheme. More...
 
void CalculateRHSContribution (Element &rCurrentElement, LocalSystemVectorType &RHS_Contribution, Element::EquationIdVectorType &EquationId, const ProcessInfo &rCurrentProcessInfo) override
 This function is designed to calculate just the RHS contribution. More...
 
void CalculateSystemContributions (Condition &rCurrentCondition, LocalSystemMatrixType &LHS_Contribution, LocalSystemVectorType &RHS_Contribution, Element::EquationIdVectorType &EquationId, const ProcessInfo &rCurrentProcessInfo) override
 Functions totally analogous to the precedent but applied to the "condition" objects. More...
 
void CalculateRHSContribution (Condition &rCurrentCondition, LocalSystemVectorType &RHS_Contribution, Element::EquationIdVectorType &EquationId, const ProcessInfo &rCurrentProcessInfo) override
 Functions that calculates the RHS of a "condition" object. More...
 
void ClearReaction () const
 
TSparseSpace::DofUpdaterPointerType mpDofUpdater
 
BossakAlphaMethod mBossak
 Bossak Alpha parameters. More...
 
GeneralMatrices mMatrix
 

Additional Inherited Members

- Public Types inherited from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
using BaseType = Scheme< TSparseSpace, TDenseSpace >
 Base type for the scheme. More...
 
using ImplicitBaseType = ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >
 Implicit base type for the scheme. More...
 
using ClassType = ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
 Class type for the scheme. More...
 
using TDataType = typename ImplicitBaseType::TDataType
 Data type used within the ImplicitBaseType. More...
 
using DofsArrayType = typename ImplicitBaseType::DofsArrayType
 Array type for degrees of freedom within ImplicitBaseType. More...
 
using DofsVectorType = typename Element::DofsVectorType
 Vector type for degrees of freedom within an Element. More...
 
using TSystemMatrixType = typename ImplicitBaseType::TSystemMatrixType
 Type for the system matrix within ImplicitBaseType. More...
 
using TSystemVectorType = typename ImplicitBaseType::TSystemVectorType
 Type for the system vector within ImplicitBaseType. More...
 
using LocalSystemVectorType = typename ImplicitBaseType::LocalSystemVectorType
 Type for local system vectors within ImplicitBaseType. More...
 
using LocalSystemMatrixType = typename ImplicitBaseType::LocalSystemMatrixType
 Type for local system matrices within ImplicitBaseType. More...
 
using NodeIterator = ModelPart::NodeIterator
 Iterator for nodes in a ModelPart. More...
 
using NodesArrayType = ModelPart::NodesContainerType
 Container type for nodes in a ModelPart. More...
 
using ElementsArrayType = ModelPart::ElementsContainerType
 Container type for elements in a ModelPart. More...
 
using ConditionsArrayType = ModelPart::ConditionsContainerType
 Container type for conditions in a ModelPart. More...
 
using BaseTypePointer = typename BaseType::Pointer
 Pointer type for the BaseType. More...
 
using ComponentType = double
 Component type as 'double'. More...
 
- Public Types inherited from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >
typedef Scheme< TSparseSpace, TDenseSpace > BaseType
 Base class definition. More...
 
typedef BaseType::DofsArrayType DofsArrayType
 DoF array type definition. More...
 
typedef Element::DofsVectorType DofsVectorType
 DoF vector type definition. More...
 
typedef BaseType::TDataType TDataType
 Data type definition. More...
 
typedef BaseType::TSystemMatrixType TSystemMatrixType
 Matrix type definition. More...
 
typedef BaseType::TSystemVectorType TSystemVectorType
 Vector type definition. More...
 
typedef BaseType::LocalSystemVectorType LocalSystemVectorType
 Local system matrix type definition. More...
 
typedef BaseType::LocalSystemMatrixType LocalSystemMatrixType
 Local system vector type definition. More...
 
typedef ModelPart::NodesContainerType NodesArrayType
 Nodes containers definition. More...
 
typedef ModelPart::ElementsContainerType ElementsArrayType
 Elements containers definition. More...
 
typedef ModelPart::ConditionsContainerType ConditionsArrayType
 Conditions containers definition. More...
 
typedef std::size_t IndexType
 Index type definition. More...
 
- Public Types inherited from Kratos::Scheme< TSparseSpace, TDenseSpace >
using ClassType = Scheme< TSparseSpace, TDenseSpace >
 The definition of the current class. More...
 
using TDataType = typename TSparseSpace::DataType
 Data type definition. More...
 
using TSystemMatrixType = typename TSparseSpace::MatrixType
 Matrix type definition. More...
 
using TSystemVectorType = typename TSparseSpace::VectorType
 Vector type definition. More...
 
using LocalSystemMatrixType = typename TDenseSpace::MatrixType
 Local system matrix type definition. More...
 
using LocalSystemVectorType = typename TDenseSpace::VectorType
 Local system vector type definition. More...
 
using TDofType = Dof< double >
 DoF type definition. More...
 
using DofsArrayType = ModelPart::DofsArrayType
 DoF array type definition. More...
 
using ElementsArrayType = ModelPart::ElementsContainerType
 Elements containers definition. More...
 
using ConditionsArrayType = ModelPart::ConditionsContainerType
 Conditions containers definition. More...
 
- Public Member Functions inherited from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
 ResidualBasedBossakDisplacementScheme (Parameters ThisParameters)
 Construct from a Parameters object. More...
 
 ResidualBasedBossakDisplacementScheme (const double Alpha=0.0)
 Constructor from a Bossak parameter. More...
 
 ResidualBasedBossakDisplacementScheme (const double Alpha, const double NewmarkBeta)
 Constructor. More...
 
 ResidualBasedBossakDisplacementScheme (ResidualBasedBossakDisplacementScheme &rOther)
 
 ~ResidualBasedBossakDisplacementScheme () override
 
std::string Info () const override
 Return information as a string. More...
 
void PrintInfo (std::ostream &rOStream) const override
 Print information about this object. More...
 
void PrintData (std::ostream &rOStream) const override
 Print the instance's data. More...
 
 KRATOS_CLASS_POINTER_DEFINITION (ResidualBasedBossakDisplacementScheme)
 
BaseType::Pointer Create (Parameters ThisParameters) const override
 Construct a dynamically allocated new scheme from a Parameters object. More...
 
void CalculateBossakCoefficients ()
 Recalculate the Newmark coefficients, taking the alpha parameters into account. More...
 
int Check (const ModelPart &rModelPart) const override
 Check whether the scheme and the provided ModelPart are configured correctly. More...
 
void Clear () override
 Release dynamic memory allocated by this instance. More...
 
Parameters GetDefaultParameters () const override
 This function returns the default Parameters to help avoiding conflicts between different constructors. More...
 
- Public Member Functions inherited from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >
 ResidualBasedImplicitTimeScheme ()
 
 ResidualBasedImplicitTimeScheme (Parameters ThisParameters)
 Constructor. The implicit method method. More...
 
 ResidualBasedImplicitTimeScheme (ResidualBasedImplicitTimeScheme &rOther)
 
 ~ResidualBasedImplicitTimeScheme () override
 
int Check (const ModelPart &rModelPart) const override
 This function is designed to be called once to perform all the checks needed on the input provided. More...
 
Parameters GetDefaultParameters () const override
 This method provides the defaults parameters to avoid conflicts between the different constructors. More...
 
std::string Info () const override
 Turn back information as a string. More...
 
void PrintInfo (std::ostream &rOStream) const override
 Print information about this object. More...
 
void PrintData (std::ostream &rOStream) const override
 Print object's data. More...
 
 KRATOS_CLASS_POINTER_DEFINITION (ResidualBasedImplicitTimeScheme)
 
- Public Member Functions inherited from Kratos::Scheme< TSparseSpace, TDenseSpace >
 Scheme ()
 Default Constructor. More...
 
 Scheme (Parameters ThisParameters)
 Constructor with Parameters. More...
 
 Scheme (Scheme &rOther)
 
virtual ~Scheme ()
 
 KRATOS_CLASS_POINTER_DEFINITION (Scheme)
 Pointer definition of Scheme. More...
 
virtual void Initialize (ModelPart &rModelPart)
 This is the place to initialize the Scheme. More...
 
bool SchemeIsInitialized ()
 This method returns if the scheme is initialized. More...
 
void SetSchemeIsInitialized (bool SchemeIsInitializedFlag=true)
 This method sets if the elements have been initialized or not (true by default) More...
 
bool ElementsAreInitialized ()
 This method returns if the elements are initialized. More...
 
void SetElementsAreInitialized (bool ElementsAreInitializedFlag=true)
 This method sets if the elements have been initialized or not (true by default) More...
 
bool ConditionsAreInitialized ()
 This method returns if the conditions are initialized. More...
 
void SetConditionsAreInitialized (bool ConditionsAreInitializedFlag=true)
 This method sets if the conditions have been initialized or not (true by default) More...
 
virtual void InitializeElements (ModelPart &rModelPart)
 This is the place to initialize the elements. More...
 
virtual void InitializeConditions (ModelPart &rModelPart)
 This is the place to initialize the conditions. More...
 
virtual void FinalizeSolutionStep (ModelPart &rModelPart, TSystemMatrixType &A, TSystemVectorType &Dx, TSystemVectorType &b)
 Function called once at the end of a solution step, after convergence is reached if an iterative process is needed. More...
 
virtual void InitializeNonLinIteration (ModelPart &rModelPart, TSystemMatrixType &A, TSystemVectorType &Dx, TSystemVectorType &b)
 unction to be called when it is needed to initialize an iteration. It is designed to be called at the beginning of each non linear iteration More...
 
virtual void CalculateOutputData (ModelPart &rModelPart, DofsArrayType &rDofSet, TSystemMatrixType &A, TSystemVectorType &Dx, TSystemVectorType &b)
 Functions to be called to prepare the data needed for the output of results. More...
 
virtual void CleanOutputData ()
 Functions that cleans the results data. More...
 
virtual void Clean ()
 This function is intended to be called at the end of the solution step to clean up memory storage not needed after the end of the solution step. More...
 
virtual int Check (ModelPart &rModelPart)
 
virtual void CalculateLHSContribution (Element &rElement, LocalSystemMatrixType &LHS_Contribution, Element::EquationIdVectorType &rEquationIdVector, const ProcessInfo &rCurrentProcessInfo)
 This function is designed to calculate just the LHS contribution. More...
 
virtual void CalculateLHSContribution (Condition &rCondition, LocalSystemMatrixType &LHS_Contribution, Element::EquationIdVectorType &rEquationIdVector, const ProcessInfo &rCurrentProcessInfo)
 Functions totally analogous to the precedent but applied to the "condition" objects. More...
 
virtual void EquationId (const Element &rElement, Element::EquationIdVectorType &rEquationId, const ProcessInfo &rCurrentProcessInfo)
 This method gets the eqaution id corresponding to the current element. More...
 
virtual void EquationId (const Condition &rCondition, Element::EquationIdVectorType &rEquationId, const ProcessInfo &rCurrentProcessInfo)
 Functions totally analogous to the precedent but applied to the "condition" objects. More...
 
virtual void GetDofList (const Element &rElement, Element::DofsVectorType &rDofList, const ProcessInfo &rCurrentProcessInfo)
 Function that returns the list of Degrees of freedom to be assembled in the system for a Given element. More...
 
virtual void GetDofList (const Condition &rCondition, Element::DofsVectorType &rDofList, const ProcessInfo &rCurrentProcessInfo)
 Function that returns the list of Degrees of freedom to be assembled in the system for a Given condition. More...
 
- Static Public Member Functions inherited from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
static std::string Name ()
 Return the name of the class as used in the settings (snake_case). More...
 
- Static Public Member Functions inherited from Kratos::Scheme< TSparseSpace, TDenseSpace >
static std::string Name ()
 Returns the name of the class as used in the settings (snake_case format) More...
 
- Protected Member Functions inherited from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
void UpdateVelocity (array_1d< double, 3 > &rCurrentVelocity, const array_1d< double, 3 > &rDeltaDisplacement, const array_1d< double, 3 > &rPreviousVelocity, const array_1d< double, 3 > &rPreviousAcceleration)
 Update the first time derivative. More...
 
void UpdateAcceleration (array_1d< double, 3 > &rCurrentAcceleration, const array_1d< double, 3 > &rDeltaDisplacement, const array_1d< double, 3 > &rPreviousVelocity, const array_1d< double, 3 > &rPreviousAcceleration)
 Update the second time derivative. More...
 
void AddDynamicsToLHS (LocalSystemMatrixType &LHS_Contribution, LocalSystemMatrixType &D, LocalSystemMatrixType &M, const ProcessInfo &rCurrentProcessInfo) override
 Add dynamic left hand side contribution from Element s. More...
 
void AddDynamicsToRHS (Element &rElement, LocalSystemVectorType &RHS_Contribution, LocalSystemMatrixType &D, LocalSystemMatrixType &M, const ProcessInfo &rCurrentProcessInfo) override
 Add dynamic right hand side contribution from an Element. More...
 
void AddDynamicsToRHS (Condition &rCondition, LocalSystemVectorType &RHS_Contribution, LocalSystemMatrixType &D, LocalSystemMatrixType &M, const ProcessInfo &rCurrentProcessInfo) override
 Add dynamic right hand side contribution of a Condition. More...
 
void AssignSettings (const Parameters ThisParameters) override
 Assign member variables from Parameters. More...
 
Protected Operations
- Protected Member Functions inherited from Kratos::Scheme< TSparseSpace, TDenseSpace >
virtual Parameters ValidateAndAssignParameters (Parameters ThisParameters, const Parameters DefaultParameters) const
 This method validate and assign default parameters. More...
 
- Protected Attributes inherited from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >
TSparseSpace::DofUpdaterPointerType mpDofUpdater = TSparseSpace::CreateDofUpdater()
 
BossakAlphaMethod mBossak
 Bossak Alpha parameters. More...
 
NewmarkMethod mNewmark
 Newmark Beta parameters. More...
 
GeneralVectors mVector
 Aggregate struct for velocities and accelerations. More...
 
- Protected Attributes inherited from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >
GeneralMatrices mMatrix
 
- Protected Attributes inherited from Kratos::Scheme< TSparseSpace, TDenseSpace >
bool mSchemeIsInitialized
 
bool mElementsAreInitialized
 Flag to be used in controlling if the Scheme has been initialized or not. More...
 
bool mConditionsAreInitialized
 Flag taking in account if the elements were initialized correctly or not. More...
 

Detailed Description

template<class TSparseSpace, class TDenseSpace>
class Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >

Bossak integration scheme (for linear and nonlinear dynamic problems) for displacements adjusted for Material Point Method.

This is an implicit scheme based of the Bossak algorithm for displacements suitable for quasi-static and dynamic problems. Furthermore, this scheme has been adjusted for mixed formulation where pressure is also solved as one of the DoFs. The parameter Alpha of Bossak introduces damping, the value of Bossak is from 0 to -0.5 (negative) Implementation according to: "An alpha modification of Newmark's method; W.L. Wood, M. Bossak, O.C. Zienkiewicz; Numerical Methods in Engineering; 1980" MPM implementation according to: "Implicit time integration for the material point method"; J. E. Guilkey, J. A. Weiss The parameter Alpha of Bossak introduces damping, the value of Bossak is from 0 to -0.5 (negative)

Member Typedef Documentation

◆ BaseType

template<class TSparseSpace , class TDenseSpace >
typedef Scheme<TSparseSpace,TDenseSpace> Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::BaseType

◆ BaseTypePointer

template<class TSparseSpace , class TDenseSpace >
typedef BaseType::Pointer Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::BaseTypePointer

◆ BossakBaseType

template<class TSparseSpace , class TDenseSpace >
typedef ResidualBasedBossakDisplacementScheme<TSparseSpace,TDenseSpace> Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::BossakBaseType

◆ ConditionsArrayType

template<class TSparseSpace , class TDenseSpace >
typedef ModelPart::ConditionsContainerType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::ConditionsArrayType

◆ DofsArrayType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::DofsArrayType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::DofsArrayType

◆ DofsVectorType

template<class TSparseSpace , class TDenseSpace >
typedef Element::DofsVectorType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::DofsVectorType

◆ ElementsArrayType

template<class TSparseSpace , class TDenseSpace >
typedef ModelPart::ElementsContainerType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::ElementsArrayType

◆ ImplicitBaseType

template<class TSparseSpace , class TDenseSpace >
typedef ResidualBasedImplicitTimeScheme<TSparseSpace,TDenseSpace> Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::ImplicitBaseType

◆ LocalSystemMatrixType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::LocalSystemMatrixType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::LocalSystemMatrixType

◆ LocalSystemVectorType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::LocalSystemVectorType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::LocalSystemVectorType

◆ TDataType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::TDataType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::TDataType

◆ TSystemMatrixType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::TSystemMatrixType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::TSystemMatrixType

◆ TSystemVectorType

template<class TSparseSpace , class TDenseSpace >
typedef BossakBaseType::TSystemVectorType Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::TSystemVectorType

Constructor & Destructor Documentation

◆ MPMResidualBasedBossakScheme() [1/2]

template<class TSparseSpace , class TDenseSpace >
Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::MPMResidualBasedBossakScheme ( ModelPart rGridModelPart,
unsigned int  DomainSize,
unsigned int  BlockSize,
double  Alpha = 0.0,
double  NewmarkBeta = 0.25,
bool  IsDynamic = true 
)
inline

Constructor. @detail The MPM bossak method.

Constructor. The MPM bossak scheme

Parameters
rGridModelPartThe background grid model part
DomainSizeDomain size or space dimension of the problems
BlockSizeBlock size of DOF
Alphais the Bossak parameter. Default value is 0, which is the Newmark method
NewmarkBetathe Newmark parameter. Default value is 0.25, for mean constant acceleration.

◆ MPMResidualBasedBossakScheme() [2/2]

template<class TSparseSpace , class TDenseSpace >
Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::MPMResidualBasedBossakScheme ( MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace > &  rOther)
inline

Copy Constructor.

◆ ~MPMResidualBasedBossakScheme()

template<class TSparseSpace , class TDenseSpace >
virtual Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::~MPMResidualBasedBossakScheme ( )
inlinevirtual

Destructor.

Member Function Documentation

◆ CalculateRHSContribution() [1/2]

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::CalculateRHSContribution ( Condition rCurrentCondition,
LocalSystemVectorType RHS_Contribution,
Element::EquationIdVectorType EquationId,
const ProcessInfo rCurrentProcessInfo 
)
inlineoverridevirtual

Functions that calculates the RHS of a "condition" object.

Parameters
rCurrentConditionThe condition to compute
rRHSContributionThe RHS vector contribution
rEquationIdThe ID's of the condition degrees of freedom
rCurrentProcessInfoThe current process info instance

Reimplemented from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >.

◆ CalculateRHSContribution() [2/2]

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::CalculateRHSContribution ( Element rCurrentElement,
LocalSystemVectorType RHS_Contribution,
Element::EquationIdVectorType EquationId,
const ProcessInfo rCurrentProcessInfo 
)
inlineoverridevirtual

This function is designed to calculate just the RHS contribution.

Parameters
rCurrentElementThe element to compute
rRHSContributionThe RHS vector contribution
rEquationIdThe ID's of the element degrees of freedom
rCurrentProcessInfoThe current process info instance

Reimplemented from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >.

◆ CalculateSystemContributions() [1/2]

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::CalculateSystemContributions ( Condition rCurrentCondition,
LocalSystemMatrixType LHS_Contribution,
LocalSystemVectorType RHS_Contribution,
Element::EquationIdVectorType EquationId,
const ProcessInfo rCurrentProcessInfo 
)
inlineoverridevirtual

Functions totally analogous to the precedent but applied to the "condition" objects.

Parameters
rCurrentConditionThe condition to compute
rLHSContributionThe LHS matrix contribution
rRHSContributionThe RHS vector contribution
rEquationIdThe ID's of the element degrees of freedom
rCurrentProcessInfoThe current process info instance

Reimplemented from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >.

◆ CalculateSystemContributions() [2/2]

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::CalculateSystemContributions ( Element rCurrentElement,
LocalSystemMatrixType LHS_Contribution,
LocalSystemVectorType RHS_Contribution,
Element::EquationIdVectorType EquationId,
const ProcessInfo rCurrentProcessInfo 
)
inlineoverridevirtual

This function is designed to be called in the builder and solver to introduce the selected time integration scheme.

Parameters
rCurrentElementThe element to compute
LHS_ContributionThe LHS matrix contribution
RHS_ContributionThe RHS vector contribution
EquationIdThe ID's of the element degrees of freedom
rCurrentProcessInfoThe current process info instance

Reimplemented from Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >.

◆ ClearReaction()

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::ClearReaction ( ) const
inlineprotected

◆ Clone()

template<class TSparseSpace , class TDenseSpace >
BaseTypePointer Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::Clone ( )
inlineoverridevirtual

◆ FinalizeNonLinIteration()

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::FinalizeNonLinIteration ( ModelPart rModelPart,
TSystemMatrixType A,
TSystemVectorType Dx,
TSystemVectorType b 
)
inlineoverridevirtual

Function to be called when it is needed to finalize an iteration. It is designed to be called at the end of each non linear iteration.

Parameters
rModelPartThe model part of the problem to solve
ALHS matrix
DxIncremental update of primary variables
bRHS Vector

Reimplemented from Kratos::Scheme< TSparseSpace, TDenseSpace >.

◆ InitializeSolutionStep()

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::InitializeSolutionStep ( ModelPart rModelPart,
TSystemMatrixType rA,
TSystemVectorType rDx,
TSystemVectorType rb 
)
inlineoverridevirtual

It initializes time step solution for MPM simulations.

The initialize solution step here also perform the following procedures:

  1. Loop over the grid nodes performed to clear all historical nodal information
  2. Assign a new nodal variables by means of extrapolation from material points
    Parameters
    rModelPartThe model of the problem to solve
    ALHS matrix
    DxIncremental update of primary variables
    bRHS Vector

Reimplemented from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >.

◆ KRATOS_CLASS_POINTER_DEFINITION()

template<class TSparseSpace , class TDenseSpace >
Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::KRATOS_CLASS_POINTER_DEFINITION ( MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >  )

◆ Predict()

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::Predict ( ModelPart rModelPart,
DofsArrayType rDofSet,
TSystemMatrixType rA,
TSystemVectorType rDx,
TSystemVectorType rb 
)
inlineoverridevirtual

Performing the prediction of the solution.

It predicts the solution for the current step x = xold + vold * Dt

Parameters
rModelPartThe model of the problem to solve
rDofSetset of all primary variables
rALHS matrix
rDxIncremental update of primary variables
rbRHS Vector

Reimplemented from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >.

◆ Update()

template<class TSparseSpace , class TDenseSpace >
void Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::Update ( ModelPart rModelPart,
DofsArrayType rDofSet,
TSystemMatrixType rA,
TSystemVectorType rDx,
TSystemVectorType rb 
)
inlineoverridevirtual

Performing the update of the solution.

Incremental update within newton iteration. It updates the state variables at the end of the time step u_{n+1}^{k+1}= u_{n+1}^{k}+ \Delta u

Parameters
rModelPartThe model of the problem to solve
rDofSetSet of all primary variables
rALHS matrix
rDxincremental update of primary variables
rbRHS Vector

Reimplemented from Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >.

Member Data Documentation

◆ mBlockSize

template<class TSparseSpace , class TDenseSpace >
unsigned int Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::mBlockSize
protected

◆ mBossak

template<class TSparseSpace , class TDenseSpace >
BossakAlphaMethod Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >::mBossak

Bossak Alpha parameters.

◆ mDomainSize

template<class TSparseSpace , class TDenseSpace >
unsigned int Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::mDomainSize
protected

◆ mGridModelPart

template<class TSparseSpace , class TDenseSpace >
ModelPart& Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::mGridModelPart
protected

◆ mIsDynamic

template<class TSparseSpace , class TDenseSpace >
bool Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::mIsDynamic
protected

◆ mMatrix

template<class TSparseSpace , class TDenseSpace >
GeneralMatrices Kratos::ResidualBasedImplicitTimeScheme< TSparseSpace, TDenseSpace >::mMatrix

◆ mpDofUpdater

template<class TSparseSpace , class TDenseSpace >
TSparseSpace::DofUpdaterPointerType Kratos::ResidualBasedBossakDisplacementScheme< TSparseSpace, TDenseSpace >::mpDofUpdater

◆ mRotationTool

template<class TSparseSpace , class TDenseSpace >
MPMBoundaryRotationUtility<LocalSystemMatrixType,LocalSystemVectorType> Kratos::MPMResidualBasedBossakScheme< TSparseSpace, TDenseSpace >::mRotationTool
protected

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