Seismo-VLAB  1.3
An Open-Source Finite Element Software for Meso-Scale Simulations
Algorithm Class Referenceabstract

Virtual class for defining how the solution between two steps is going to be carried out. More...

#include <Algorithm.hpp>

Inheritance diagram for Algorithm:

Public Member Functions

 Algorithm (const std::shared_ptr< Mesh > &mesh, unsigned int flag=1, double NormFactor=1.0)
 Creates a Algorithm object. More...
 
virtual ~Algorithm ()=0
 Destroys this Algorithm object. More...
 
virtual bool ComputeNewIncrement (std::shared_ptr< Mesh > &mesh, unsigned int i)=0
 Computes a new incremental solution. More...
 
virtual const Eigen::VectorXd & GetDisplacementIncrement ()=0
 Gets the displacement increment. More...
 
virtual void SetIntegrator (std::shared_ptr< Integrator > &integrator)=0
 Sets the integrator for this algorithm. More...
 
virtual void SetLoadFactor (double factor)=0
 Set the load factor. More...
 
void ReverseStatesIncrements (std::shared_ptr< Mesh > &mesh)
 Reverse the incremental state variables to its pervious converged one in the mesh. More...
 
void UpdateStatesIncrements (std::shared_ptr< Mesh > &mesh, const Eigen::VectorXd &dU)
 Update the incremental state variables in the mesh. More...
 
void ReducedParallelResidual (const Eigen::VectorXd &Feff, double &Residual)
 Construct the residual vector force from each processor. More...
 
double ComputeConvergence (const Eigen::VectorXd &Force, const Eigen::VectorXd &Delta, const Eigen::VectorXd &delta, unsigned int k)
 Computes convergence tests for this algorithm. More...
 

Protected Attributes

Eigen::SparseMatrix< double > Total2FreeMatrix
 Operator that enforced restrain/constraint. More...
 

Private Attributes

unsigned int flag
 Convergence test to be performed. More...
 
double NormFactor
 Relative Unbalanced Factor. More...
 

Detailed Description

Virtual class for defining how the solution between two steps is going to be carried out.

See also
Linear.hpp NewtonRaphson.hpp

Constructor & Destructor Documentation

◆ Algorithm()

Algorithm::Algorithm ( const std::shared_ptr< Mesh > &  mesh,
unsigned int  flag = 1,
double  NormFactor = 1.0 
)

Creates a Algorithm object.

Parameters
meshPointer to the Mesh container to extract Node and Element.
flagThe convergence test to be performed.
NormFactorThe relative unbalanced factor for convergence.
See also
Algorithm::flag, Algorithm::NormFactor.

◆ ~Algorithm()

virtual Algorithm::~Algorithm ( )
pure virtual

Destroys this Algorithm object.

Member Function Documentation

◆ ComputeConvergence()

double Algorithm::ComputeConvergence ( const Eigen::VectorXd &  Force,
const Eigen::VectorXd &  Delta,
const Eigen::VectorXd &  delta,
unsigned int  k 
)

Computes convergence tests for this algorithm.

Parameters
ForceThe residual force vector.
DeltaThe accumulated incremental displacement vector.
deltaThe incremental displacement vector .
isFirstIterationWhether this step in Algorithm is the first iteration.
Note
If Algorithm::flag = 1: Unbalanced Force Norm.
If Algorithm::flag = 2: Increment Displacement Norm.
If Algorithm::flag = 3: Energy Incremental Norm.
If Algorithm::flag = 4: Relative Unbalanced Force Norm.
If Algorithm::flag = 5: Relative Incremental Displacement Norm.
If Algorithm::flag = 6: Relative Energy Incremental Norm.
If Algorithm::flag = 7: Total Relative Increment Displacement Norm.
If Algorithm::flag = 8: Maximum Number of Iterations.

◆ ComputeNewIncrement()

virtual bool Algorithm::ComputeNewIncrement ( std::shared_ptr< Mesh > &  mesh,
unsigned int  i 
)
pure virtual

Computes a new incremental solution.

Parameters
meshThe finite element Mesh object.
iThe time step number to be solved.
Returns
Whether or not the algorithm was successful.

Implemented in NewtonRaphson, and Linear.

◆ GetDisplacementIncrement()

virtual const Eigen::VectorXd& Algorithm::GetDisplacementIncrement ( )
pure virtual

Gets the displacement increment.

Returns
The incremental displacement vector.

Implemented in NewtonRaphson, and Linear.

◆ ReducedParallelResidual()

void Algorithm::ReducedParallelResidual ( const Eigen::VectorXd &  Feff,
double &  Residual 
)

Construct the residual vector force from each processor.

Parameters
FeffThe residual vector in this partition.
ResidualThe value of the residual vector norm.
Note
In parallel execution Feff is gathered across processors.

◆ ReverseStatesIncrements()

void Algorithm::ReverseStatesIncrements ( std::shared_ptr< Mesh > &  mesh)

Reverse the incremental state variables to its pervious converged one in the mesh.

Parameters
meshPointer to the Mesh container to extract Node and Element. note This function zeroes out incremental displacements at each Node in Mesh.
See also
Node::SetIncrementalDisplacements(), Element::ReverseState().

◆ SetIntegrator()

virtual void Algorithm::SetIntegrator ( std::shared_ptr< Integrator > &  integrator)
pure virtual

Sets the integrator for this algorithm.

Parameters
integratorPointer to the Integrator to obtain the effective stiffness and force.

Implemented in NewtonRaphson, and Linear.

◆ SetLoadFactor()

virtual void Algorithm::SetLoadFactor ( double  factor)
pure virtual

Set the load factor.

Parameters
factorThe incremental load factor.

Implemented in NewtonRaphson, and Linear.

◆ UpdateStatesIncrements()

void Algorithm::UpdateStatesIncrements ( std::shared_ptr< Mesh > &  mesh,
const Eigen::VectorXd &  dU 
)

Update the incremental state variables in the mesh.

Parameters
meshPointer to the Mesh container to extract Node and Element.
dUVector with the converged incremental displacements. note This function sets incremental displacements at each Node in Mesh.
See also
Node::SetIncrementalDisplacements() Element::UpdateState().

Member Data Documentation

◆ flag

unsigned int Algorithm::flag
private

Convergence test to be performed.

◆ NormFactor

double Algorithm::NormFactor
private

Relative Unbalanced Factor.

◆ Total2FreeMatrix

Eigen::SparseMatrix<double> Algorithm::Total2FreeMatrix
protected

Operator that enforced restrain/constraint.

See also
Mesh::GetTotalToFreeMatrix().