AbstractNonlinearElasticitySolver< DIM > Class Template Reference

#include <AbstractNonlinearElasticitySolver.hpp>

Inheritance diagram for AbstractNonlinearElasticitySolver< DIM >:

Inheritance graph
[legend]
Collaboration diagram for AbstractNonlinearElasticitySolver< DIM >:

Collaboration graph
[legend]

List of all members.

Public Member Functions

 AbstractNonlinearElasticitySolver (unsigned numDofs, AbstractIncompressibleMaterialLaw< DIM > *pMaterialLaw, c_vector< double, DIM > bodyForce, double density, std::string outputDirectory, std::vector< unsigned > &fixedNodes)
 AbstractNonlinearElasticitySolver (unsigned numDofs, std::vector< AbstractIncompressibleMaterialLaw< DIM > * > &rMaterialLaws, c_vector< double, DIM > bodyForce, double density, std::string outputDirectory, std::vector< unsigned > &fixedNodes)
virtual ~AbstractNonlinearElasticitySolver ()
void Solve (double tol=-1.0, unsigned offset=0, unsigned maxNumNewtonIterations=INT_MAX, bool quitIfNoConvergence=true)
void WriteOutput (unsigned counter)
unsigned GetNumNewtonIterations ()
void SetFunctionalBodyForce (c_vector< double, DIM >(*pFunction)(c_vector< double, DIM > &))
void SetWriteOutput (bool writeOutput=true)
void SetKspAbsoluteTolerance (double kspAbsoluteTolerance)
std::vector< double > & rGetCurrentSolution ()

Protected Member Functions

virtual void FormInitialGuess ()=0
virtual void AssembleSystem (bool assembleResidual, bool assembleJacobian)=0
virtual std::vector< c_vector
< double, DIM > > & 
rGetDeformedPosition ()=0
void ApplyBoundaryConditions (bool applyToMatrix)
double ComputeResidualAndGetNorm (bool allowException)
double CalculateResidualNorm ()
void VectorSum (std::vector< double > &rX, ReplicatableVector &rY, double a, std::vector< double > &rZ)
void PrintLineSearchResult (double s, double residNorm)
double TakeNewtonStep ()
double UpdateSolutionUsingLineSearch (Vec solution)
virtual void PostNewtonStep (unsigned counter, double normResidual)

Protected Attributes

double mKspAbsoluteTol
unsigned mNumDofs
std::vector
< AbstractIncompressibleMaterialLaw
< DIM > * > 
mMaterialLaws
LinearSystemmpLinearSystem
LinearSystemmpPreconditionMatrixLinearSystem
c_vector< double, DIM > mBodyForce
double mDensity
std::string mOutputDirectory
std::vector< unsigned > mFixedNodes
std::vector< c_vector< double,
DIM > > 
mFixedNodeDisplacements
bool mWriteOutput
std::vector< double > mCurrentSolution
FourthOrderTensor< DIM, DIM,
DIM, DIM > 
dTdE
unsigned mNumNewtonIterations
std::vector< c_vector< double,
DIM > > 
mDeformedPosition
std::vector< double > mPressures
std::vector< c_vector< double,
DIM > > 
mSurfaceTractions
c_vector< double, DIM >(* mpBodyForceFunction )(c_vector< double, DIM > &)
c_vector< double, DIM >(* mpTractionBoundaryConditionFunction )(c_vector< double, DIM > &)
bool mUsingBodyForceFunction
bool mUsingTractionBoundaryConditionFunction

Static Protected Attributes

static const double MAX_NEWTON_ABS_TOL = 1e-7
static const double MIN_NEWTON_ABS_TOL = 1e-10
static const double NEWTON_REL_TOL = 1e-4


Detailed Description

template<unsigned DIM>
class AbstractNonlinearElasticitySolver< DIM >

Abstract nonlinear elasticity solver.

Definition at line 63 of file AbstractNonlinearElasticitySolver.hpp.


Constructor & Destructor Documentation

template<unsigned DIM>
AbstractNonlinearElasticitySolver< DIM >::AbstractNonlinearElasticitySolver ( unsigned  numDofs,
AbstractIncompressibleMaterialLaw< DIM > *  pMaterialLaw,
c_vector< double, DIM >  bodyForce,
double  density,
std::string  outputDirectory,
std::vector< unsigned > &  fixedNodes 
) [inline]

Constructor.

Parameters:
numDofs 
pMaterialLaw 
bodyForce 
density 
outputDirectory 
fixedNodes 

Definition at line 786 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mMaterialLaws, AbstractNonlinearElasticitySolver< DIM >::mOutputDirectory, and AbstractNonlinearElasticitySolver< DIM >::mWriteOutput.

template<unsigned DIM>
AbstractNonlinearElasticitySolver< DIM >::AbstractNonlinearElasticitySolver ( unsigned  numDofs,
std::vector< AbstractIncompressibleMaterialLaw< DIM > * > &  rMaterialLaws,
c_vector< double, DIM >  bodyForce,
double  density,
std::string  outputDirectory,
std::vector< unsigned > &  fixedNodes 
) [inline]

Variant constructor taking a vector of material laws.

Parameters:
numDofs 
rMaterialLaws 
bodyForce 
density 
outputDirectory 
fixedNodes 

Definition at line 813 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mMaterialLaws, AbstractNonlinearElasticitySolver< DIM >::mOutputDirectory, and AbstractNonlinearElasticitySolver< DIM >::mWriteOutput.

template<unsigned DIM>
AbstractNonlinearElasticitySolver< DIM >::~AbstractNonlinearElasticitySolver (  )  [inline, virtual]


Member Function Documentation

template<unsigned DIM>
virtual void AbstractNonlinearElasticitySolver< DIM >::FormInitialGuess (  )  [protected, pure virtual]

Set up the current guess to be the solution given no displacement. Must be overridden in concrete derived classes.

Implemented in NonlinearElasticitySolver< DIM >.

template<unsigned DIM>
virtual void AbstractNonlinearElasticitySolver< DIM >::AssembleSystem ( bool  assembleResidual,
bool  assembleJacobian 
) [protected, pure virtual]

Assemble the residual vector (using the current solution stored in mCurrentSolution, output going to mpLinearSystem->rGetRhsVector), or Jacobian matrix (using the current solution stored in mCurrentSolution, output going to mpLinearSystem->rGetLhsMatrix). Must be overridden in concrete derived classes.

Parameters:
assembleResidual A bool stating whether to assemble the residual vector.
assembleJacobian A bool stating whether to assemble the Jacobian matrix.

Implemented in NonlinearElasticitySolver< DIM >.

Referenced by AbstractNonlinearElasticitySolver< DIM >::ComputeResidualAndGetNorm(), and AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep().

template<unsigned DIM>
virtual std::vector<c_vector<double,DIM> >& AbstractNonlinearElasticitySolver< DIM >::rGetDeformedPosition (  )  [protected, pure virtual]

Get the deformed position. Must be overridden in concrete derived classes.

Implemented in NonlinearElasticitySolver< DIM >.

Referenced by AbstractNonlinearElasticitySolver< DIM >::WriteOutput().

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::ApplyBoundaryConditions ( bool  applyToMatrix  )  [inline, protected]

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::ComputeResidualAndGetNorm ( bool  allowException  )  [inline, protected]

Set up the residual vector (using the current solution), and get its scaled norm (Calculate |r|_2 / length(r), where r is residual vector)

Parameters:
allowException Sometimes the current solution solution will be such that the residual vector cannot be computed, as (say) the material law will throw an exception as the strains are too large. If this bool is set to true, the exception will be caught, and DBL_MAX returned as the residual norm

Definition at line 432 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::AssembleSystem(), and AbstractNonlinearElasticitySolver< DIM >::CalculateResidualNorm().

Referenced by AbstractNonlinearElasticitySolver< DIM >::Solve(), and AbstractNonlinearElasticitySolver< DIM >::UpdateSolutionUsingLineSearch().

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::CalculateResidualNorm (  )  [inline, protected]

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::VectorSum ( std::vector< double > &  rX,
ReplicatableVector rY,
double  a,
std::vector< double > &  rZ 
) [inline, protected]

Simple helper function, computes Z = X + aY, where X and Z are std::vectors and Y a ReplicatableVector

Parameters:
rX X
rY Y (replicatable vector)
a a
rZ Z the returned vector

Definition at line 467 of file AbstractNonlinearElasticitySolver.hpp.

References ReplicatableVector::GetSize().

Referenced by AbstractNonlinearElasticitySolver< DIM >::UpdateSolutionUsingLineSearch().

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::PrintLineSearchResult ( double  s,
double  residNorm 
) [inline, protected]

Print to std::cout the residual norm for this s, ie ||f(x+su)|| where f is the residual vector, x the current solution and u the update vector

Parameters:
s s
residNorm residual norm.

Definition at line 626 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::UpdateSolutionUsingLineSearch().

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep (  )  [inline, protected]

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::UpdateSolutionUsingLineSearch ( Vec  solution  )  [inline, protected]

Using the update vector (of Newton's method), choose s such that ||f(x+su)|| is most decreased, where f is the residual vector, x the current solution (mCurrentSolution) and u the update vector. This checks s=1 first (most likely to be the current solution, then 0.9, 0.8.. until ||f|| starts increasing.

Parameters:
solution The solution of the linear solve in newton's method, ie the update vector u.

Definition at line 634 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::CalculateResidualNorm(), AbstractNonlinearElasticitySolver< DIM >::ComputeResidualAndGetNorm(), AbstractNonlinearElasticitySolver< DIM >::mCurrentSolution, AbstractNonlinearElasticitySolver< DIM >::PrintLineSearchResult(), and AbstractNonlinearElasticitySolver< DIM >::VectorSum().

Referenced by AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep().

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::PostNewtonStep ( unsigned  counter,
double  normResidual 
) [inline, protected, virtual]

This function may be overloaded by subclasses. It is called after each Newton iteration.

Parameters:
counter 
normResidual 

Definition at line 780 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::Solve().

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::Solve ( double  tol = -1.0,
unsigned  offset = 0,
unsigned  maxNumNewtonIterations = INT_MAX,
bool  quitIfNoConvergence = true 
) [inline]

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::WriteOutput ( unsigned  counter  )  [inline]

template<unsigned DIM>
unsigned AbstractNonlinearElasticitySolver< DIM >::GetNumNewtonIterations (  )  [inline]

Get number of Newton iterations taken in last solve.

Definition at line 961 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mNumNewtonIterations.

Referenced by ImplicitCardiacMechanicsSolver< DIM >::Solve().

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::SetFunctionalBodyForce ( c_vector< double, DIM >(*)(c_vector< double, DIM > &)  pFunction  )  [inline]

Set a function which gives body force as a function of X (undeformed position) Whatever body force was provided in the constructor will now be ignored.

Parameters:
pFunction 

Definition at line 968 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mpBodyForceFunction, and AbstractNonlinearElasticitySolver< DIM >::mUsingBodyForceFunction.

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::SetWriteOutput ( bool  writeOutput = true  )  [inline]

Set whether to write any output.

Parameters:
writeOutput (defaults to true)

Definition at line 976 of file AbstractNonlinearElasticitySolver.hpp.

References EXCEPTION, AbstractNonlinearElasticitySolver< DIM >::mOutputDirectory, and AbstractNonlinearElasticitySolver< DIM >::mWriteOutput.

template<unsigned DIM>
void AbstractNonlinearElasticitySolver< DIM >::SetKspAbsoluteTolerance ( double  kspAbsoluteTolerance  )  [inline]

Set the absolute tolerance to be used when solving the linear system. If this is not called a relative tolerance is used.

Parameters:
kspAbsoluteTolerance the tolerance

Definition at line 368 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mKspAbsoluteTol.

template<unsigned DIM>
std::vector<double>& AbstractNonlinearElasticitySolver< DIM >::rGetCurrentSolution (  )  [inline]

Get the current solution vector (advanced use only - for getting the deformed position use rGetDeformedPosition()

Definition at line 378 of file AbstractNonlinearElasticitySolver.hpp.

References AbstractNonlinearElasticitySolver< DIM >::mCurrentSolution.


Member Data Documentation

template<unsigned DIM>
const double AbstractNonlinearElasticitySolver< DIM >::MAX_NEWTON_ABS_TOL = 1e-7 [inline, static, protected]

Maximum absolute tolerance for Newton solve. The Newton solver uses the absolute tolerance corresponding to the specified relative tolerance, but has a max and min allowable absolute tolerance. (ie: if max_abs = 1e-7, min_abs = 1e-10, rel=1e-4: then if the norm of the initial_residual (=a) is 1e-2, it will solve with tolerance 1e-7; if a=1e-5, it will solve with tolerance 1e-9; a=1e-9, it will solve with tolerance 1e-10.

Definition at line 72 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::Solve().

template<unsigned DIM>
const double AbstractNonlinearElasticitySolver< DIM >::MIN_NEWTON_ABS_TOL = 1e-10 [inline, static, protected]

Minimum absolute tolerance for Newton solve. See documentation for MAX_NEWTON_ABS_TOL.

Definition at line 75 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::Solve().

template<unsigned DIM>
const double AbstractNonlinearElasticitySolver< DIM >::NEWTON_REL_TOL = 1e-4 [inline, static, protected]

Relative tolerance for Newton solve. See documentation for MAX_NEWTON_ABS_TOL.

Definition at line 78 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::Solve().

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::mKspAbsoluteTol [protected]

Absolute tolerance for linear systems. Can be set by calling SetKspAbsoluteTolerances(), but default to -1, in which case a relative tolerance is used.

Definition at line 83 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::SetKspAbsoluteTolerance(), and AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep().

template<unsigned DIM>
unsigned AbstractNonlinearElasticitySolver< DIM >::mNumDofs [protected]

Number of degrees of freedom (eg equal to DIM*N + M if quadratic-linear bases are used, where there are N total nodes and M vertices).

Definition at line 89 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AllocateMatrixMemory(), NonlinearElasticitySolver< DIM >::AssembleSystem(), AbstractNonlinearElasticitySolver< DIM >::CalculateResidualNorm(), and NonlinearElasticitySolver< DIM >::FormInitialGuess().

template<unsigned DIM>
std::vector<AbstractIncompressibleMaterialLaw<DIM>*> AbstractNonlinearElasticitySolver< DIM >::mMaterialLaws [protected]

template<unsigned DIM>
LinearSystem* AbstractNonlinearElasticitySolver< DIM >::mpLinearSystem [protected]

The linear system where we store all residual vectors which are calculated and the Jacobian. Note we don't actually call Solve but solve using Petsc methods explicitly (in order to easily set number of restarts etc).

Definition at line 103 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AllocateMatrixMemory(), AbstractNonlinearElasticitySolver< DIM >::ApplyBoundaryConditions(), NonlinearElasticitySolver< DIM >::AssembleSystem(), AbstractNonlinearElasticitySolver< DIM >::CalculateResidualNorm(), AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep(), and AbstractNonlinearElasticitySolver< DIM >::~AbstractNonlinearElasticitySolver().

The linear system which stores the matrix used for preconditioning (given the helper functions on LinearSystem it is best to use LinearSystem and use these for assembling the preconditioner, rather than just use a Mat The preconditioner is the petsc LU factorisation of

Jp = [A B] in displacement-pressure block form, [C M]

where the A, B and C are the matrices in the normal jacobian, ie

J = [A B] [C 0]

and M is the MASS MATRIX (ie integral phi_i phi_j dV, where phi_i are the pressure basis functions).

Definition at line 123 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AllocateMatrixMemory(), AbstractNonlinearElasticitySolver< DIM >::ApplyBoundaryConditions(), NonlinearElasticitySolver< DIM >::AssembleSystem(), AbstractNonlinearElasticitySolver< DIM >::TakeNewtonStep(), and AbstractNonlinearElasticitySolver< DIM >::~AbstractNonlinearElasticitySolver().

template<unsigned DIM>
c_vector<double,DIM> AbstractNonlinearElasticitySolver< DIM >::mBodyForce [protected]

Body force vector

Definition at line 126 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnElement().

template<unsigned DIM>
double AbstractNonlinearElasticitySolver< DIM >::mDensity [protected]

Mass density of the undeformed body (equal to the density of deformed body)

Definition at line 129 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnElement().

template<unsigned DIM>
std::string AbstractNonlinearElasticitySolver< DIM >::mOutputDirectory [protected]

template<unsigned DIM>
std::vector<unsigned> AbstractNonlinearElasticitySolver< DIM >::mFixedNodes [protected]

All nodes (including non-vertices) which are fixed

Definition at line 135 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::ApplyBoundaryConditions(), and NonlinearElasticitySolver< DIM >::Initialise().

template<unsigned DIM>
std::vector<c_vector<double,DIM> > AbstractNonlinearElasticitySolver< DIM >::mFixedNodeDisplacements [protected]

The displacements of those nodes with displacement boundary conditions

Definition at line 138 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by AbstractNonlinearElasticitySolver< DIM >::ApplyBoundaryConditions(), and NonlinearElasticitySolver< DIM >::Initialise().

template<unsigned DIM>
bool AbstractNonlinearElasticitySolver< DIM >::mWriteOutput [protected]

template<unsigned DIM>
std::vector<double> AbstractNonlinearElasticitySolver< DIM >::mCurrentSolution [protected]

template<unsigned DIM>
FourthOrderTensor<DIM,DIM,DIM,DIM> AbstractNonlinearElasticitySolver< DIM >::dTdE [protected]

Storage space for a 4th order tensor used in assembling the Jacobian (to avoid repeated memory allocation)

Definition at line 154 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnElement().

template<unsigned DIM>
unsigned AbstractNonlinearElasticitySolver< DIM >::mNumNewtonIterations [protected]

template<unsigned DIM>
std::vector<c_vector<double,DIM> > AbstractNonlinearElasticitySolver< DIM >::mDeformedPosition [protected]

Deformed position: mDeformedPosition[i](j) = x_j for node i

Definition at line 160 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::rGetDeformedPosition().

template<unsigned DIM>
std::vector<double> AbstractNonlinearElasticitySolver< DIM >::mPressures [protected]

The solution pressures. mPressures[i] = pressure at node i (ie vertex i).

Definition at line 166 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::rGetPressures().

template<unsigned DIM>
std::vector<c_vector<double,DIM> > AbstractNonlinearElasticitySolver< DIM >::mSurfaceTractions [protected]

The surface tractions (which should really be non-zero) for the boundary elements in mBoundaryElements.

Definition at line 172 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleSystem(), and NonlinearElasticitySolver< DIM >::SetSurfaceTractionBoundaryConditions().

template<unsigned DIM>
c_vector<double,DIM>(* AbstractNonlinearElasticitySolver< DIM >::mpBodyForceFunction)(c_vector< double, DIM > &) [protected]

An optionally provided (pointer to a) function, giving body force as a function of undeformed position.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnElement(), and AbstractNonlinearElasticitySolver< DIM >::SetFunctionalBodyForce().

template<unsigned DIM>
c_vector<double,DIM>(* AbstractNonlinearElasticitySolver< DIM >::mpTractionBoundaryConditionFunction)(c_vector< double, DIM > &) [protected]

An optionally provided (pointer to a) function, giving the surface traction as a function of undeformed position.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnBoundaryElement(), and NonlinearElasticitySolver< DIM >::SetFunctionalTractionBoundaryCondition().

template<unsigned DIM>
bool AbstractNonlinearElasticitySolver< DIM >::mUsingBodyForceFunction [protected]

Whether the functional version of the body force is being used or not

Definition at line 184 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnElement(), and AbstractNonlinearElasticitySolver< DIM >::SetFunctionalBodyForce().

template<unsigned DIM>
bool AbstractNonlinearElasticitySolver< DIM >::mUsingTractionBoundaryConditionFunction [protected]

Whether the functional version of the surface traction is being used or not

Definition at line 187 of file AbstractNonlinearElasticitySolver.hpp.

Referenced by NonlinearElasticitySolver< DIM >::AssembleOnBoundaryElement(), and NonlinearElasticitySolver< DIM >::SetFunctionalTractionBoundaryCondition().


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

Generated on Mon Nov 1 12:35:50 2010 for Chaste by  doxygen 1.5.5