LifeV
VenantKirchhoffViscoelasticSolver< Mesh, SolverType > Class Template Reference

SecondOrderProblem this class solver second order problem, waves equation and linear viscoelastic problem. More...

#include <VenantKirchhoffViscoelasticSolver.hpp>

+ Collaboration diagram for VenantKirchhoffViscoelasticSolver< Mesh, SolverType >:

Protected Attributes

std::shared_ptr< data_typeM_data
 data of problem More...
 
std::shared_ptr< FESpace< Mesh, MapEpetra > > M_FESpace
 feSpace More...
 
std::unique_ptr< DisplayerM_displayer
 displayer More...
 
Int M_me
 current rank More...
 
bchandler_type M_BCh
 boundary condition More...
 
std::shared_ptr< const MapEpetraM_localMap
 Epetra map need to define the VectorEpetra;. More...
 
matrix_ptrtype M_matrMass
 Matrix mass. More...
 
matrix_ptrtype M_matrLinearStiffness
 Matrix stiffness linear. More...
 
matrix_ptrtype M_matrSystem
 Matrix System. More...
 
matrix_ptrtype M_matrDamping
 Matrix Damping $\gamma A + \beta M$. More...
 
vector_ptrtype M_solution
 unknowns vector More...
 
vector_ptrtype M_rhs
 right hand side More...
 
vector_ptrtype M_rhsNoBC
 right hand side without boundary condition More...
 
std::shared_ptr< vector_typeM_residual
 residual More...
 
source_type M_source
 source term More...
 
std::shared_ptr< solver_typeM_linearSolver
 data for solving tangent problem with aztec More...
 
bool M_reusePrec
 true if reuse the preconditonar More...
 
UInt M_maxIterForReuse
 max number of iteration before to update preconditonator More...
 
bool M_resetPrec
 reset preconditionator More...
 
UInt M_maxIterSolver
 maximum number of iteration for solver method More...
 
UInt M_offset
 offset More...
 

Private Member Functions

void applyBoundaryConditions (matrix_type &matrix, vector_type &rhs, bchandler_raw_type &BCh, UInt offset=0)
 apply boundary condition More...
 

Public Types

typedef Real(* Function) (const Real &, const Real &, const Real &, const Real &, const ID &)
 
typedef std::function< Real(Real const &, Real const &, Real const &, Real const &, ID const &) > source_type
 
typedef BCHandler bchandler_raw_type
 
typedef std::shared_ptr< bchandler_raw_typebchandler_type
 
typedef SolverType solver_type
 
typedef solver_type::matrix_type matrix_type
 
typedef std::shared_ptr< matrix_typematrix_ptrtype
 
typedef solver_type::vector_type vector_type
 
typedef std::shared_ptr< vector_typevector_ptrtype
 
typedef SolverType::prec_raw_type prec_raw_type
 
typedef SolverType::prec_type prec_type
 
typedef VenantKirchhoffViscoelasticData data_type
 

Constructors & Destructor

 VenantKirchhoffViscoelasticSolver ()
 Constructor. More...
 
virtual ~VenantKirchhoffViscoelasticSolver ()
 virtual Destructor More...
 

Methods

void setup (std::shared_ptr< data_type > data, const std::shared_ptr< FESpace< Mesh, MapEpetra > > &feSpace, std::shared_ptr< Epetra_Comm > &comm)
 class setup More...
 
void setup (std::shared_ptr< data_type > data, const std::shared_ptr< FESpace< Mesh, MapEpetra > > &feSpace, bchandler_type &BCh, std::shared_ptr< Epetra_Comm > &comm)
 class setup More...
 
virtual void setup (std::shared_ptr< data_type > data, const std::shared_ptr< FESpace< Mesh, MapEpetra > > &feSpace, std::shared_ptr< Epetra_Comm > &comm, const std::shared_ptr< const MapEpetra > &epetraMap, UInt offset=0)
 class setup More...
 
void buildSystem (const Real &xi=1, const Real &alpha=0)
 buildSystem More...
 
void buildSystem (matrix_ptrtype matrix, const Real &xi)
 buildSystem More...
 
void buildDamping (matrix_ptrtype damping, const Real &alpha)
 build Damping matrix More...
 
void updateSystem (const vector_type &rhs, const Real &xi=0, const Real &alpha=0)
 updateSystem with coefficients of time advance methods More...
 
void computeStartValue (vector_type &solution, const bchandler_raw_type &bch, const vector_type &rhs)
 compute the start solution More...
 
void iterate (bchandler_raw_type &bch)
 Solve the system. More...
 
void updateSourceTerm (const vector_type &source)
 update source term More...
 
void updateRHS (const vector_type &rhs)
 updateRHS More...
 
void resetPrec ()
 reset the Prec More...
 
void setSourceTerm (source_type const &source)
 Sets source term. More...
 
void setDataFromGetPot (const GetPot &dataFile)
 Set parameters. More...
 
MapEpetra const & map () const
 Return the map. More...
 
Displayer const & displayer () const
 Return the map. More...
 
matrix_type const matrMass () const
 Return the mass matrix. More...
 
matrix_ptrtype const matrDamping () const
 Return the damping matrix. More...
 
matrix_ptrtype const matrSystem () const
 Return the damping matrix. More...
 
matrix_ptrtype const matrLinearStiff () const
 Return the Stiffness matrix. More...
 
FESpace< Mesh, MapEpetra > & feSpace ()
 Return FESpace. More...
 
BCHandler const & BChandler () const
 BCHandler getter and setter. More...
 
vector_ptrtypesolution ()
 Return the solution. More...
 
vector_ptrtyperesidual ()
 Return residual. More...
 
vector_ptrtyperhsContributionSecondDerivativeithoutBC ()
 Return right hand side without boundary conditions. More...
 
std::shared_ptr< Epetra_Comm > const & comm () const
 return the communicator More...
 
const UIntoffset () const
 Return offset;. More...
 
const Realthickness () const
 thickness More...
 
const Realdensity () const
 density More...
 
const Realyoung () const
 young More...
 
const Realpoisson () const
 poisson More...
 

Elementary matrices and vectors:

std::shared_ptr< MatrixElementalM_elmatK
 linear stiffness More...
 
std::shared_ptr< MatrixElementalM_elmatM
 mass More...
 
std::shared_ptr< MatrixElementalM_elmatC
 mass+ linear stiffness More...
 
std::shared_ptr< MatrixElementalM_elmatD
 damping More...
 

Detailed Description

template<typename Mesh, typename SolverType = LifeV::SolverAztecOO>
class LifeV::VenantKirchhoffViscoelasticSolver< Mesh, SolverType >

SecondOrderProblem this class solver second order problem, waves equation and linear viscoelastic problem.

This class solves the following waves equation:

\[M \frac{d^2 u}{d t^2} + D(u,\frac{d u}{d t} ) \frac{d u}{d t} + A(u, \frac{d u}{dt}) = f\]

where $M$ is mass matrix, $A$ stiffness matrix and $D$ is damping matrix given by $ \beta M + \gamma A$.

If $A$ and $D$ depend on $u$ and $dot{u}$ we linearize the problem using suitable extrapolations. we use the time advancing scheme defined in TimeAdvanceBase class for details see TimeAdvanceBase.pdf notes.

Definition at line 87 of file VenantKirchhoffViscoelasticSolver.hpp.

Member Typedef Documentation

◆ Function

typedef Real( * Function) (const Real &, const Real &, const Real &, const Real &, const ID &)

Definition at line 94 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ source_type

typedef std::function<Real ( Real const&, Real const&, Real const&, Real const&, ID const& ) > source_type

Definition at line 95 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ bchandler_raw_type

◆ bchandler_type

typedef std::shared_ptr<bchandler_raw_type> bchandler_type

Definition at line 98 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ solver_type

typedef SolverType solver_type

Definition at line 100 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrix_type

typedef solver_type::matrix_type matrix_type

Definition at line 102 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrix_ptrtype

typedef std::shared_ptr<matrix_type> matrix_ptrtype

Definition at line 103 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ vector_type

typedef solver_type::vector_type vector_type

Definition at line 104 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ vector_ptrtype

typedef std::shared_ptr<vector_type> vector_ptrtype

Definition at line 105 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ prec_raw_type

typedef SolverType::prec_raw_type prec_raw_type

Definition at line 107 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ prec_type

typedef SolverType::prec_type prec_type

Definition at line 108 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ data_type

Constructor & Destructor Documentation

◆ VenantKirchhoffViscoelasticSolver()

Constructor.

Empty Constructor.

Definition at line 536 of file VenantKirchhoffViscoelasticSolver.hpp.

+ Here is the caller graph for this function:

◆ ~VenantKirchhoffViscoelasticSolver()

virtual ~VenantKirchhoffViscoelasticSolver ( )
inlinevirtual

virtual Destructor

Definition at line 121 of file VenantKirchhoffViscoelasticSolver.hpp.

Member Function Documentation

◆ setup() [1/3]

void setup ( std::shared_ptr< data_type data,
const std::shared_ptr< FESpace< Mesh, MapEpetra > > &  feSpace,
std::shared_ptr< Epetra_Comm > &  comm 
)

class setup

Parameters
datafile
feSpacefinite element space
commcommunicator

Definition at line 600 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ setup() [2/3]

void setup ( std::shared_ptr< data_type data,
const std::shared_ptr< FESpace< Mesh, MapEpetra > > &  feSpace,
bchandler_type BCh,
std::shared_ptr< Epetra_Comm > &  comm 
)

class setup

Parameters
datafile
feSpacefinite element space
BChboundary condition
commcommunicator

Definition at line 615 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ setup() [3/3]

void setup ( std::shared_ptr< data_type data,
const std::shared_ptr< FESpace< Mesh, MapEpetra > > &  feSpace,
std::shared_ptr< Epetra_Comm > &  comm,
const std::shared_ptr< const MapEpetra > &  epetraMap,
UInt  offset = 0 
)
virtual

class setup

Parameters
datafile
feSpacefinite element space
commcommunicator
epetraMapepetra vector
offset

Definition at line 572 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ buildSystem() [1/2]

void buildSystem ( const Real xi = 1,
const Real alpha = 0 
)

buildSystem

Parameters
xiis the coefficient of mass term defined as $\xi^0/dt^2$;
alphais the coefficient of damping term defined as $\alpha^0/dt$;
Note
by default $xi$ is equal to 1 and $alpha$ is equal to 0;

Definition at line 637 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ buildSystem() [2/2]

void buildSystem ( matrix_ptrtype  matrix,
const Real xi 
)

buildSystem

Parameters
matrixis the Mass $ + $ Damping $+$ Stiffness
xiis the coefficient of mass term defined as $\xi^0/dt^2$;

Definition at line 660 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ buildDamping()

void buildDamping ( matrix_ptrtype  damping,
const Real alpha 
)

build Damping matrix

building damping matrix defined as $\gamma A + \beta M$

Parameters
dampingmatrix
alphais the coefficient of time advancing scheme defined as $\alpha^0/dt$ alpha can be 1 and we can introduce the time advancing coefficient in the updateSystem

Definition at line 729 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ updateSystem()

void updateSystem ( const vector_type rhs,
const Real xi = 0,
const Real alpha = 0 
)

updateSystem with coefficients of time advance methods

Definition at line 773 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ computeStartValue()

void computeStartValue ( vector_type solution,
const bchandler_raw_type bch,
const vector_type rhs 
)

compute the start solution

we consider the system $M w + D v + A u = f$; in general we knowns the $u^0$, $v^0$ but not $w^0$. This method computes $w^0$ starting by $F = f^0 - D v^0 -A u^0$ and boundary conditions associated to $w$;

Parameters
solutionis the vector where save the solution;
bchthe boundary condition respect to $w$;
rhsis $f^0 - D v^0 -A u^0$

Definition at line 812 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ iterate()

void iterate ( bchandler_raw_type bch)

Solve the system.

solve the system

Parameters
bchis boundary condition;

Definition at line 874 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ updateSourceTerm()

void updateSourceTerm ( const vector_type source)

update source term

Parameters
sourceis the vector where we evaluate the source term;

Definition at line 835 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ updateRHS()

void updateRHS ( const vector_type rhs)
inline

updateRHS

Parameters
rhsis the right and side;

Definition at line 234 of file VenantKirchhoffViscoelasticSolver.hpp.

+ Here is the caller graph for this function:

◆ resetPrec()

void resetPrec ( )
inline

reset the Prec

Definition at line 249 of file VenantKirchhoffViscoelasticSolver.hpp.

+ Here is the caller graph for this function:

◆ setSourceTerm()

void setSourceTerm ( source_type const &  source)
inline

Sets source term.

Parameters
sourceis function

Definition at line 263 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ setDataFromGetPot()

void setDataFromGetPot ( const GetPot dataFile)

Set parameters.

Parameters
dataFileis the file contains the parameters of problem this method must be removed

Definition at line 628 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ map()

MapEpetra const& map ( ) const
inline

Return the map.

Returns
epetraMap

Definition at line 284 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ displayer()

Displayer const& displayer ( ) const
inline

Return the map.

Returns
displayer

Definition at line 293 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrMass()

matrix_type const matrMass ( ) const
inline

Return the mass matrix.

Returns
the mass matrix (it doesn't divided by $\xi/dt^2$)

Definition at line 302 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrDamping()

matrix_ptrtype const matrDamping ( ) const
inline

Return the damping matrix.

Returns
the damping matrix (it doesn't divided by $\alpha/dt$)

Definition at line 311 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrSystem()

matrix_ptrtype const matrSystem ( ) const
inline

Return the damping matrix.

Returns
the system matrix given by Stiffness + Damping $\alpha/dt$ + Mass $\xi/dt^2$

Definition at line 320 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ matrLinearStiff()

matrix_ptrtype const matrLinearStiff ( ) const
inline

Return the Stiffness matrix.

Returns
the stiffness matrix

Definition at line 329 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ feSpace()

FESpace<Mesh, MapEpetra>& feSpace ( )
inline

Return FESpace.

Definition at line 335 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ BChandler()

BCHandler const& BChandler ( ) const
inline

BCHandler getter and setter.

Returns
the boundary conditions

Definition at line 344 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ solution()

vector_ptrtype& solution ( )
inline

Return the solution.

Returns
the solution

Definition at line 353 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ residual()

vector_ptrtype& residual ( )
inline

Return residual.

Definition at line 362 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ rhsContributionSecondDerivativeithoutBC()

vector_ptrtype& rhsContributionSecondDerivativeithoutBC ( )
inline

Return right hand side without boundary conditions.

Definition at line 371 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ comm()

std::shared_ptr<Epetra_Comm> const& comm ( ) const
inline

return the communicator

Returns
the communicator

Definition at line 380 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ offset()

const UInt& offset ( ) const
inline

Return offset;.

Definition at line 389 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ thickness()

const Real& thickness ( ) const
inline

thickness

Returns
the thickness

Definition at line 398 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ density()

const Real& density ( ) const
inline

density

Returns
the density

Definition at line 406 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ young()

const Real& young ( ) const
inline

young

Returns
the young

Definition at line 415 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ poisson()

const Real& poisson ( ) const
inline

poisson

Returns
the poisson

Definition at line 423 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ applyBoundaryConditions()

void applyBoundaryConditions ( matrix_type matrix,
vector_type rhs,
bchandler_raw_type BCh,
UInt  offset = 0 
)
private

apply boundary condition

Parameters
matrixof system;
rhsright hand side without boundaryCondition;
BChboundary condition;
offset

Definition at line 925 of file VenantKirchhoffViscoelasticSolver.hpp.

+ Here is the caller graph for this function:

Field Documentation

◆ M_data

std::shared_ptr<data_type> M_data
protected

data of problem

Definition at line 450 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_FESpace

std::shared_ptr<FESpace<Mesh, MapEpetra> > M_FESpace
protected

feSpace

Definition at line 453 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_displayer

std::unique_ptr<Displayer> M_displayer
protected

displayer

Definition at line 456 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_me

Int M_me
protected

current rank

Definition at line 459 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_BCh

bchandler_type M_BCh
protected

boundary condition

Definition at line 462 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_localMap

std::shared_ptr<const MapEpetra> M_localMap
protected

Epetra map need to define the VectorEpetra;.

Definition at line 465 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_matrMass

matrix_ptrtype M_matrMass
protected

Matrix mass.

Definition at line 468 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_matrLinearStiffness

matrix_ptrtype M_matrLinearStiffness
protected

Matrix stiffness linear.

Definition at line 471 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_matrSystem

matrix_ptrtype M_matrSystem
protected

Matrix System.

Definition at line 474 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_matrDamping

matrix_ptrtype M_matrDamping
protected

Matrix Damping $\gamma A + \beta M$.

Definition at line 477 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_elmatK

std::shared_ptr<MatrixElemental> M_elmatK
protected

linear stiffness

Definition at line 482 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_elmatM

std::shared_ptr<MatrixElemental> M_elmatM
protected

mass

Definition at line 485 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_elmatC

std::shared_ptr<MatrixElemental> M_elmatC
protected

mass+ linear stiffness

Definition at line 487 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_elmatD

std::shared_ptr<MatrixElemental> M_elmatD
protected

damping

Definition at line 490 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_solution

vector_ptrtype M_solution
protected

unknowns vector

Definition at line 494 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_rhs

vector_ptrtype M_rhs
protected

right hand side

Definition at line 497 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_rhsNoBC

vector_ptrtype M_rhsNoBC
protected

right hand side without boundary condition

Definition at line 500 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_residual

std::shared_ptr<vector_type> M_residual
protected

residual

Definition at line 503 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_source

source_type M_source
protected

source term

Definition at line 506 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_linearSolver

std::shared_ptr<solver_type> M_linearSolver
protected

data for solving tangent problem with aztec

Definition at line 509 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_reusePrec

bool M_reusePrec
protected

true if reuse the preconditonar

Definition at line 512 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_maxIterForReuse

UInt M_maxIterForReuse
protected

max number of iteration before to update preconditonator

Definition at line 515 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_resetPrec

bool M_resetPrec
protected

reset preconditionator

Definition at line 518 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_maxIterSolver

UInt M_maxIterSolver
protected

maximum number of iteration for solver method

Definition at line 521 of file VenantKirchhoffViscoelasticSolver.hpp.

◆ M_offset

UInt M_offset
protected

offset

Definition at line 524 of file VenantKirchhoffViscoelasticSolver.hpp.


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