LifeV
|
SecondOrderProblem this class solver second order problem, waves equation and linear viscoelastic problem. More...
#include <VenantKirchhoffViscoelasticSolver.hpp>
Protected Attributes | |
std::shared_ptr< data_type > | M_data |
data of problem More... | |
std::shared_ptr< FESpace< Mesh, MapEpetra > > | M_FESpace |
feSpace More... | |
std::unique_ptr< Displayer > | M_displayer |
displayer More... | |
Int | M_me |
current rank More... | |
bchandler_type | M_BCh |
boundary condition More... | |
std::shared_ptr< const MapEpetra > | M_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 . 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_type > | M_residual |
residual More... | |
source_type | M_source |
source term More... | |
std::shared_ptr< solver_type > | M_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_type > | bchandler_type |
typedef SolverType | solver_type |
typedef solver_type::matrix_type | matrix_type |
typedef std::shared_ptr< matrix_type > | matrix_ptrtype |
typedef solver_type::vector_type | vector_type |
typedef std::shared_ptr< vector_type > | vector_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_ptrtype & | solution () |
Return the solution. More... | |
vector_ptrtype & | residual () |
Return residual. More... | |
vector_ptrtype & | rhsContributionSecondDerivativeithoutBC () |
Return right hand side without boundary conditions. More... | |
std::shared_ptr< Epetra_Comm > const & | comm () const |
return the communicator More... | |
const UInt & | offset () const |
Return offset;. More... | |
const Real & | thickness () const |
thickness More... | |
const Real & | density () const |
density More... | |
const Real & | young () const |
young More... | |
const Real & | poisson () const |
poisson More... | |
Elementary matrices and vectors: | |
std::shared_ptr< MatrixElemental > | M_elmatK |
linear stiffness More... | |
std::shared_ptr< MatrixElemental > | M_elmatM |
mass More... | |
std::shared_ptr< MatrixElemental > | M_elmatC |
mass+ linear stiffness More... | |
std::shared_ptr< MatrixElemental > | M_elmatD |
damping More... | |
SecondOrderProblem this class solver second order problem, waves equation and linear viscoelastic problem.
This class solves the following waves equation:
where is mass matrix, stiffness matrix and is damping matrix given by .
If and depend on and 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.
Definition at line 94 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef std::function<Real ( Real const&, Real const&, Real const&, Real const&, ID const& ) > source_type |
Definition at line 95 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef BCHandler bchandler_raw_type |
Definition at line 97 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef std::shared_ptr<bchandler_raw_type> bchandler_type |
Definition at line 98 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef SolverType solver_type |
Definition at line 100 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef solver_type::matrix_type matrix_type |
Definition at line 102 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef std::shared_ptr<matrix_type> matrix_ptrtype |
Definition at line 103 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef solver_type::vector_type vector_type |
Definition at line 104 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef std::shared_ptr<vector_type> vector_ptrtype |
Definition at line 105 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef SolverType::prec_raw_type prec_raw_type |
Definition at line 107 of file VenantKirchhoffViscoelasticSolver.hpp.
typedef SolverType::prec_type prec_type |
Definition at line 108 of file VenantKirchhoffViscoelasticSolver.hpp.
Definition at line 110 of file VenantKirchhoffViscoelasticSolver.hpp.
Constructor.
Empty Constructor.
Definition at line 536 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inlinevirtual |
virtual Destructor
Definition at line 121 of file VenantKirchhoffViscoelasticSolver.hpp.
void setup | ( | std::shared_ptr< data_type > | data, |
const std::shared_ptr< FESpace< Mesh, MapEpetra > > & | feSpace, | ||
std::shared_ptr< Epetra_Comm > & | comm | ||
) |
class setup
data | file |
feSpace | finite element space |
comm | communicator |
Definition at line 600 of file VenantKirchhoffViscoelasticSolver.hpp.
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
data | file |
feSpace | finite element space |
BCh | boundary condition |
comm | communicator |
Definition at line 615 of file VenantKirchhoffViscoelasticSolver.hpp.
|
virtual |
class setup
data | file |
feSpace | finite element space |
comm | communicator |
epetraMap | epetra vector |
offset |
Definition at line 572 of file VenantKirchhoffViscoelasticSolver.hpp.
buildSystem
xi | is the coefficient of mass term defined as ; |
alpha | is the coefficient of damping term defined as ; |
Definition at line 637 of file VenantKirchhoffViscoelasticSolver.hpp.
void buildSystem | ( | matrix_ptrtype | matrix, |
const Real & | xi | ||
) |
buildSystem
matrix | is the Mass Damping Stiffness |
xi | is the coefficient of mass term defined as ; |
Definition at line 660 of file VenantKirchhoffViscoelasticSolver.hpp.
void buildDamping | ( | matrix_ptrtype | damping, |
const Real & | alpha | ||
) |
build Damping matrix
building damping matrix defined as
damping | matrix |
alpha | is the coefficient of time advancing scheme defined as alpha can be 1 and we can introduce the time advancing coefficient in the updateSystem |
Definition at line 729 of file VenantKirchhoffViscoelasticSolver.hpp.
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.
void computeStartValue | ( | vector_type & | solution, |
const bchandler_raw_type & | bch, | ||
const vector_type & | rhs | ||
) |
compute the start solution
we consider the system ; in general we knowns the , but not . This method computes starting by and boundary conditions associated to ;
solution | is the vector where save the solution; |
bch | the boundary condition respect to ; |
rhs | is |
Definition at line 812 of file VenantKirchhoffViscoelasticSolver.hpp.
void iterate | ( | bchandler_raw_type & | bch | ) |
Solve the system.
solve the system
bch | is boundary condition; |
Definition at line 874 of file VenantKirchhoffViscoelasticSolver.hpp.
void updateSourceTerm | ( | const vector_type & | source | ) |
update source term
source | is the vector where we evaluate the source term; |
Definition at line 835 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
updateRHS
rhs | is the right and side; |
Definition at line 234 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
reset the Prec
Definition at line 249 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Sets source term.
source | is function |
Definition at line 263 of file VenantKirchhoffViscoelasticSolver.hpp.
void setDataFromGetPot | ( | const GetPot & | dataFile | ) |
Set parameters.
dataFile | is the file contains the parameters of problem this method must be removed |
Definition at line 628 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the map.
Definition at line 284 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the map.
Definition at line 293 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the mass matrix.
Definition at line 302 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the damping matrix.
Definition at line 311 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the damping matrix.
Definition at line 320 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the Stiffness matrix.
Definition at line 329 of file VenantKirchhoffViscoelasticSolver.hpp.
Return FESpace.
Definition at line 335 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
BCHandler getter and setter.
Definition at line 344 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return the solution.
Definition at line 353 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return residual.
Definition at line 362 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return right hand side without boundary conditions.
Definition at line 371 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
return the communicator
Definition at line 380 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
Return offset;.
Definition at line 389 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
thickness
Definition at line 398 of file VenantKirchhoffViscoelasticSolver.hpp.
|
inline |
|
inline |
|
inline |
|
private |
apply boundary condition
matrix | of system; |
rhs | right hand side without boundaryCondition; |
BCh | boundary condition; |
offset |
Definition at line 925 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
data of problem
Definition at line 450 of file VenantKirchhoffViscoelasticSolver.hpp.
feSpace
Definition at line 453 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
displayer
Definition at line 456 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
current rank
Definition at line 459 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
boundary condition
Definition at line 462 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
Epetra map need to define the VectorEpetra;.
Definition at line 465 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
Matrix mass.
Definition at line 468 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
Matrix stiffness linear.
Definition at line 471 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
Matrix System.
Definition at line 474 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
Matrix Damping .
Definition at line 477 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
linear stiffness
Definition at line 482 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
mass
Definition at line 485 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
mass+ linear stiffness
Definition at line 487 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
damping
Definition at line 490 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
unknowns vector
Definition at line 494 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
right hand side
Definition at line 497 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
right hand side without boundary condition
Definition at line 500 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
residual
Definition at line 503 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
source term
Definition at line 506 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
data for solving tangent problem with aztec
Definition at line 509 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
true if reuse the preconditonar
Definition at line 512 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
max number of iteration before to update preconditonator
Definition at line 515 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
reset preconditionator
Definition at line 518 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
maximum number of iteration for solver method
Definition at line 521 of file VenantKirchhoffViscoelasticSolver.hpp.
|
protected |
offset
Definition at line 524 of file VenantKirchhoffViscoelasticSolver.hpp.