LifeV
MonolithicBlockMatrix Class Reference

MonolithicBlockMatrix - class which handles matrices with a block structure. More...

#include <MonolithicBlockMatrix.hpp>

+ Inheritance diagram for MonolithicBlockMatrix:
+ Collaboration diagram for MonolithicBlockMatrix:

Public Types

typedef MonolithicBlock super_Type
 
typedef super_Type::fespacePtr_Type fespacePtr_Type
 
typedef super_Type::vector_Type vector_Type
 
typedef super_Type::vectorPtr_Type vectorPtr_Type
 
typedef super_Type::solverPtr_Type solverPtr_Type
 
typedef super_Type::matrix_Type matrix_Type
 
typedef super_Type::matrixPtr_Type matrixPtr_Type
 
typedef super_Type::epetraOperatorPtr_Type epetraOperatorPtr_Type
 
typedef super_Type::mapPtr_Type mapPtr_Type
 
typedef FactorySingleton< Factory< MonolithicBlockMatrix, std::string > > Factory_Type
 

Constructor & Destructor

 MonolithicBlockMatrix (const std::vector< Int > &flags)
 
 ~MonolithicBlockMatrix ()
 

Virtual methods

virtual void setDataFromGetPot (const GetPot &data, const std::string &section)
 Sets the parameters needed by the class from data file. More...
 
virtual void setupSolver (solver_Type &solver, const GetPot &data)
 Sets the parameters needed by the preconditioner from data file. More...
 
virtual void GlobalAssemble ()
 runs GlobalAssemble on the blocks More...
 
virtual void coupler (mapPtr_Type &map, const std::map< ID, ID > &locDofMap, const vectorPtr_Type &numerationInterface, const Real &timeStep, const Real &coefficient, const Real &rescaleFactor)
 Computes the coupling matrix. More...
 
void coupler (mapPtr_Type &, const std::map< ID, ID > &locDofMap, const vectorPtr_Type &numerationInterface, const Real &timeStep, const Real &coefficient, const Real &rescaleFactor, UInt couplingBlock)
 Adds a coupling part to the already existing coupling matrix. More...
 
virtual bool set ()
 returns true if the operator has at least one block More...
 

Public methods

int solveSystem (const vector_Type &rhs, vector_Type &step, solverPtr_Type &linearSolver)
 Solves the preconditioned linear system (used only when blockMatrix is used as a preconditioner) More...
 
void push_back_matrix (const matrixPtr_Type &Mat, bool)
 pushes a block at the end of the vector More...
 
void replace_matrix (const matrixPtr_Type &Mat, UInt index)
 replaces a block More...
 
void replace_coupling (const matrixPtr_Type &, UInt)
 replaces the coupling matrix without copies More...
 
void replace_precs (const epetraOperatorPtr_Type &Mat, UInt index)
 never used within this class More...
 
void blockAssembling ()
 sums the coupling matrix in the specified position with the global matrix More...
 
matrixPtr_Typematrix ()
 returns the global matrix, with all the blocks and the coupling parts More...
 
void applyPreconditioner (matrixPtr_Type robinCoupling, matrixPtr_Type prec, vectorPtr_Type &rhs)
 multiplies the whole system times a matrix More...
 
void applyPreconditioner (const matrixPtr_Type robinCoupling, matrixPtr_Type &prec)
 multiplies two matrices More...
 
void applyPreconditioner (const matrixPtr_Type matrix, vectorPtr_Type &rhsFull)
 applies a matrix to the system More...
 
void createInterfaceMap (const MapEpetra &interfaceMap, const std::map< ID, ID > &locDofMap, const UInt subdomainMaxId, const std::shared_ptr< Epetra_Comm > epetraWorldComm)
 creates the map for the coupling More...
 
void applyBoundaryConditions (const Real &time)
 applies the b.c. to every block More...
 
void applyBoundaryConditions (const Real &time, const UInt block)
 applies the b.c. to a specific block More...
 
void applyBoundaryConditions (const Real &time, vectorPtr_Type &rhs)
 applies all the b.c.s to the global matrix, updates the rhs More...
 
void applyBoundaryConditions (const Real &time, vectorPtr_Type &rhs, const UInt block)
 applies the b.c. to the a specific block, updates the rhs More...
 
void addToCoupling (const matrixPtr_Type &Mat, UInt)
 adds a block to the coupling matrix More...
 
void addToCoupling (const Real &entry, UInt row, UInt col, UInt position)
 
void addToGlobalMatrix (const matrixPtr_Type &Mat)
 adds a block to the coupling matrix More...
 
void push_back_coupling (matrixPtr_Type &coupling)
 adds a coupling block to the coupling matrix More...
 

Get Methods

UInt interface ()
 returns the dimension of the interface More...
 
mapPtr_Type interfaceMap () const
 returns the map built for theLagrange multipliers More...
 
matrixPtr_Type coupling () const
 
void numerationInterface (vectorPtr_Type &numeration)
 returns the numeration of the interface More...
 
const UInt whereIsBlock (UInt) const
 

Factory Method

static MonolithicBlockMatrixcreateAdditiveSchwarz ()
 

Protected Members

matrixPtr_Type M_globalMatrix
 
matrixPtr_Type M_coupling
 
mapPtr_Type M_interfaceMap
 
UInt M_interface
 

Private Members

std::unique_ptr< std::vector< Int > > M_couplingFlags
 
vectorPtr_Type M_numerationInterface
 

Additional Inherited Members

- Public Types inherited from MonolithicBlock
typedef VectorEpetra vector_Type
 
typedef std::shared_ptr< vector_TypevectorPtr_Type
 
typedef MatrixEpetra< Realmatrix_Type
 
typedef std::shared_ptr< matrix_TypematrixPtr_Type
 
typedef std::shared_ptr< Epetra_OperatorepetraOperatorPtr_Type
 
typedef std::shared_ptr< Preconditionerepetra_preconditioner_ptrtype
 
typedef matrix_Type::matrix_type epetraMatrix_Type
 
typedef SolverAztecOO solver_Type
 
typedef std::shared_ptr< SolverAztecOOsolverPtr_Type
 
typedef FESpace< RegionMesh< LinearTetra >, MapEpetrafespace_Type
 
typedef std::shared_ptr< fespace_TypefespacePtr_Type
 
typedef std::shared_ptr< MapEpetramapPtr_Type
 
typedef std::shared_ptr< BCHandlerbchandlerPtr_Type
 
- Public Member Functions inherited from MonolithicBlock
virtual void replace_bch (bchandlerPtr_Type &, UInt)
 replaces a BCHandler More...
 
virtual void resetBlocks ()
 resets the blocks (frees the shared pointers) More...
 
virtual void setComm (std::shared_ptr< Epetra_Comm > comm)
 sets the communicator More...
 
virtual void reset ()
 resets the blocks, boundary conditions, FE spaces. More...
 
virtual void setRobin (matrixPtr_Type &, vectorPtr_Type &)
 Applies the robin preconditioners. More...
 
void couplingMatrix (matrixPtr_Type &bigMatrix, Int flag, const std::vector< fespacePtr_Type > &problem, const std::vector< UInt > &offset, const std::map< ID, ID > &locDofMap, const vectorPtr_Type &numerationInterface, const Real &timeStep=1.e-3, const Real &value=1., const Real &coefficient=1., const Real &rescaleFactor=1.)
 builds the coupling matrix. More...
 
void setConditions (std::vector< bchandlerPtr_Type > &vec)
 sets the vector of raw pointer to the BCHandler More...
 
void setSpaces (std::vector< fespacePtr_Type > &vec)
 sets the vector of raw pointer to the FESpaces More...
 
void setOffsets (UInt blocks,...)
 sets the vector of raw pointer to the offsets of the different blocks More...
 
void robinCoupling (MonolithicBlock::matrixPtr_Type &matrix, Real &alphaf, Real &alphas, UInt coupling, const MonolithicBlock::fespacePtr_Type &FESpace1, const UInt &offset1, const MonolithicBlock::fespacePtr_Type &FESpace2, const UInt &offset2, const std::map< ID, ID > &locDofMap, const MonolithicBlock::vectorPtr_Type &numerationInterface)
 computes the Robin coupling matrix More...
 
virtual void addToBlock (const matrixPtr_Type &Mat, UInt position)
 
virtual void push_back_oper (MonolithicBlock &Oper)
 Pushes a new block. More...
 
 MonolithicBlock ()
 Empty Constructor. More...
 
virtual ~MonolithicBlock ()
 Destructor. More...
 
virtual void setRecompute (UInt, bool)
 If not present in the derived class it must not be called (gives an assertion fail) More...
 
virtual void push_back_precs (const epetraOperatorPtr_Type &)
 pushes back a block preconditioner More...
 
const std::vector< matrixPtr_Type > & blockVector ()
 returns the vector of pointers to the blocks (by const reference). More...
 
const std::vector< bchandlerPtr_Type > & BChVector ()
 returns the vector of pointers to the BCHandlers (by const reference). More...
 
const std::vector< fespacePtr_Type > & FESpaceVector ()
 returns the vector of pointers to the FE spaces (by const reference). More...
 
const std::vector< UInt > & offsetVector ()
 returns the vector of the offsets (by const reference). More...
 
- Protected Member Functions inherited from MonolithicBlock
virtual void blockAssembling (const UInt)
 sums the coupling matrix in the specified position with the corresponding block More...
 
template<typename Operator >
void swap (std::shared_ptr< Operator > &operFrom, std::shared_ptr< Operator > &OperTo)
 swaps two std::shared_ptr. The tamplate argument of the shared_ptr is templated More...
 
template<typename Operator >
void insert (std::vector< Operator > &operFrom, std::vector< Operator > &OperTo)
 swaps two std::shared_ptr. The tamplate argument of the shared_ptr is templated More...
 
- Protected Attributes inherited from MonolithicBlock
std::vector< bchandlerPtr_TypeM_bch
 
std::vector< matrixPtr_TypeM_blocks
 
std::vector< fespacePtr_TypeM_FESpace
 
std::vector< UIntM_offset
 
vectorPtr_Type M_numerationInterface
 
std::shared_ptr< Epetra_Comm > M_comm
 

Detailed Description

MonolithicBlockMatrix - class which handles matrices with a block structure.

Author
Paolo Crosetto

This class handles a matrix with block structure. In particular a vector of shared_ptr points to the different blocks, while one matrix contains the coupling part. all the blocks and the coupling are summed into the matrix M_globalMatrix for the solution of the linear system. If a MonolithicBlockMatrix is used as a preconditioner by default an algebraic additive Schwarz preconditioner is built on the matrix M_globalMatrix.

Definition at line 55 of file MonolithicBlockMatrix.hpp.

Member Typedef Documentation

◆ super_Type

Definition at line 61 of file MonolithicBlockMatrix.hpp.

◆ fespacePtr_Type

◆ vector_Type

Definition at line 63 of file MonolithicBlockMatrix.hpp.

◆ vectorPtr_Type

◆ solverPtr_Type

◆ matrix_Type

Definition at line 66 of file MonolithicBlockMatrix.hpp.

◆ matrixPtr_Type

◆ epetraOperatorPtr_Type

◆ mapPtr_Type

Definition at line 69 of file MonolithicBlockMatrix.hpp.

◆ Factory_Type

Definition at line 70 of file MonolithicBlockMatrix.hpp.

Constructor & Destructor Documentation

◆ MonolithicBlockMatrix()

MonolithicBlockMatrix ( const std::vector< Int > &  flags)
inline

Definition at line 77 of file MonolithicBlockMatrix.hpp.

◆ ~MonolithicBlockMatrix()

~MonolithicBlockMatrix ( )
inline

Definition at line 87 of file MonolithicBlockMatrix.hpp.

Member Function Documentation

◆ setDataFromGetPot()

void setDataFromGetPot ( const GetPot data,
const std::string &  section 
)
virtual

Sets the parameters needed by the class from data file.

Public Methods.

Parameters
dataGetPot object reading the text data file
sectionstring specifying the path in the data file where to find the options for the operator

Implements MonolithicBlock.

Reimplemented in MonolithicBlockMatrixRN.

Definition at line 42 of file MonolithicBlockMatrix.cpp.

◆ setupSolver()

void setupSolver ( solver_Type ,
const GetPot  
)
virtual

Sets the parameters needed by the preconditioner from data file.

Parameters
dataGetPot object reading the text data file
sectionstring specifying the path in the data file where to find the options for the operator

Reimplemented from MonolithicBlock.

Definition at line 46 of file MonolithicBlockMatrix.cpp.

◆ GlobalAssemble()

void GlobalAssemble ( )
virtual

runs GlobalAssemble on the blocks

closes and distributes all the blocks, sums all the blocks together into M_globalMatrix.

Implements MonolithicBlock.

Reimplemented in MonolithicBlockMatrixRN.

Definition at line 52 of file MonolithicBlockMatrix.cpp.

◆ coupler() [1/2]

void coupler ( mapPtr_Type map,
const std::map< ID, ID > &  locDofMap,
const vectorPtr_Type numerationInterface,
const Real timeStep,
const Real coefficient = 1.,
const Real rescaleFactor = 1. 
)
virtual

Computes the coupling matrix.

computes the coupling matrix for the system (or for the preconditioner). The coupling is handled through an augmented formulation, introducing new variables (multipliers). Needs as input: the global map of the problem, the two FESpaces of the subproblems to be coupled with their offsets, a std::map holding the two numerations of the interface between the two subproblems (the numeration can be different but the nodes must be matching in each subdomain), an EpetraVector defined on the multipliers map containing the corresponding dof at the interface (NB: the multipliers map should be constructed from the second numeration in the std::map)

Parameters
mapthe map of the global problem
FESpace1FESpace of the first problem
offset1offset for the first block in the global matrix
FESpace2FESpace of the second problem
offset2offset for the second block in the global matrix
locDofMapstd::map with the correspondence between the interface dofs for the two different maps in the subproblems
numerationInterfacevector containing the correspondence of the Lagrange multipliers with the interface dofs

Implements MonolithicBlock.

Reimplemented in MonolithicBlockMatrixRN.

Definition at line 58 of file MonolithicBlockMatrix.cpp.

◆ coupler() [2/2]

void coupler ( mapPtr_Type ,
const std::map< ID, ID > &  locDofMap,
const vectorPtr_Type numerationInterface,
const Real timeStep,
const Real coefficient,
const Real rescaleFactor,
UInt  couplingBlock 
)
virtual

Adds a coupling part to the already existing coupling matrix.

The coupling is handled through an augmented formulation, introducing new variables (multipliers). Needs as input: the global map of the problem, the two FESpaces of the subproblems to be coupled with their offsets, a std::map holding the two numerations of the interface between the two subproblems (the numeration can be different but the nodes must be matching in each subdomain), an EpetraVector defined on the multipliers map containing the corresponding dof at the interface (NB: the multipliers map should be constructed from the second numeration in the std::map)

Parameters
mapthe map of the global problem
locDofMapstd::map with the correspondence between the interface dofs for the two different maps in the subproblems
numerationInterfacevector containing the correspondence of the Lagrange multipliers with the interface dofs
unusedflag kept for compliance with the base class
couplingFlagflag specifying which block must be coupled whith which block.

Implements MonolithicBlock.

Reimplemented in MonolithicBlockMatrixRN.

Definition at line 72 of file MonolithicBlockMatrix.cpp.

◆ set()

virtual bool set ( )
inlinevirtual

returns true if the operator has at least one block

Implements MonolithicBlock.

Definition at line 158 of file MonolithicBlockMatrix.hpp.

◆ solveSystem()

int solveSystem ( const vector_Type rhs,
vector_Type step,
solverPtr_Type linearSolver 
)
virtual

Solves the preconditioned linear system (used only when blockMatrix is used as a preconditioner)

Provided the linear solver and the right hand side this method computes the solution and returns it into the result vector.

Parameters
rhsright hand side of the linear system
resultoutput result
linearSolverthe linear system

Implements MonolithicBlock.

Definition at line 85 of file MonolithicBlockMatrix.cpp.

◆ push_back_matrix()

void push_back_matrix ( const matrixPtr_Type Mat,
bool   
)
virtual

pushes a block at the end of the vector

adds a new block

Parameters
Matblock matrix to push
recomputeflag stating wether the preconditioner for this block have to be recomputed at every time step. In this case it is not used since it is equal to the boolean specifying wether the whole preconditioner must be recomputed or not.

Implements MonolithicBlock.

Definition at line 90 of file MonolithicBlockMatrix.cpp.

◆ replace_matrix()

void replace_matrix ( const matrixPtr_Type Mat,
UInt  index 
)
virtual

replaces a block

replaces a block on a specified position in the vector

Parameters
Matblock matrix to push
indexposition in the vector

Implements MonolithicBlock.

Definition at line 96 of file MonolithicBlockMatrix.cpp.

◆ replace_coupling()

void replace_coupling ( const matrixPtr_Type ,
UInt   
)
inlinevirtual

replaces the coupling matrix without copies

this method is deprecated, it is implemented for compatibility with the base class

Parameters
Matreplacing matrix

Implements MonolithicBlock.

Definition at line 202 of file MonolithicBlockMatrix.hpp.

◆ replace_precs()

void replace_precs ( const epetraOperatorPtr_Type Mat,
UInt  index 
)
virtual

never used within this class

runs assert(false) when called. This method is used only when the preconditioner is a composed operator

Reimplemented from MonolithicBlock.

Definition at line 102 of file MonolithicBlockMatrix.cpp.

◆ blockAssembling()

void blockAssembling ( )
virtual

sums the coupling matrix in the specified position with the global matrix

Everything (but the boundary conditions) must have been set before calling this

Reimplemented from MonolithicBlock.

Reimplemented in MonolithicBlockMatrixRN.

Definition at line 108 of file MonolithicBlockMatrix.cpp.

◆ matrix()

matrixPtr_Type& matrix ( )
inline

returns the global matrix, with all the blocks and the coupling parts

Definition at line 223 of file MonolithicBlockMatrix.hpp.

◆ applyPreconditioner() [1/3]

void applyPreconditioner ( matrixPtr_Type  robinCoupling,
matrixPtr_Type  prec,
vectorPtr_Type rhs 
)

multiplies the whole system times a matrix

applies a matrix robinCoupling to the system matrix M_globalMatrix, to the rpeconditioner prec passed as input, to the rhs passed as input. This is useful e.g. when we want to substitute the Dirichlet-Neumann coupling with a Robin-Robin one.

Parameters
robinCouplingmatrix that multiplies the system
precpreconditioner matrix
rhsright hand side of the system

Definition at line 129 of file MonolithicBlockMatrix.cpp.

◆ applyPreconditioner() [2/3]

void applyPreconditioner ( const matrixPtr_Type  robinCoupling,
matrixPtr_Type prec 
)

multiplies two matrices

multiplies on the left the matrix prec times the matrix robinCoupling

Parameters
robinCouplingthe matrix that multiplies
precthe matrix multiplied (modified)

Definition at line 137 of file MonolithicBlockMatrix.cpp.

◆ applyPreconditioner() [3/3]

void applyPreconditioner ( const matrixPtr_Type  matrix,
vectorPtr_Type rhsFull 
)

applies a matrix to the system

multiplies on the left the system matrix M_globalMatrix and the rhs times the matrix matrix.

Parameters
matrixthe preconditioning matrix
rhsFullthe right hand side of the system

Definition at line 122 of file MonolithicBlockMatrix.cpp.

◆ createInterfaceMap()

void createInterfaceMap ( const MapEpetra interfaceMap,
const std::map< ID, ID > &  locDofMap,
const UInt  subdomainMaxId,
const std::shared_ptr< Epetra_Comm >  epetraWorldComm 
)

creates the map for the coupling

This method create the coupling map needed for two blocks. The coupling is handled adding a number of Lagrange multipliers equal to the interface degrees of freedom. This method generates the map for the vector of the Lagrange multipliers.

Parameters
interfaceMapthe map of the interface (not contiguous, the numeration of this map has "holes")
locDofMapthe map between the two (e. g. in FSI the fluid (first) and the solid (second)) numerations of the interface
subdomainMaxIdThe maximum ID of the sub-map built from interfaceMap (is the number of dof in the subdomain used to buid the interfaceMap divided by the number of dimensions). The interface map used to build the vector M_numerationInterface contains once every node, not every interface degree of freedom (that is nDimensions * interfaceNodes)
epetraWorldCommThe communicator

Definition at line 149 of file MonolithicBlockMatrix.cpp.

◆ applyBoundaryConditions() [1/4]

void applyBoundaryConditions ( const Real time)
virtual

applies the b.c. to every block

Reimplemented from MonolithicBlock.

Definition at line 216 of file MonolithicBlockMatrix.cpp.

◆ applyBoundaryConditions() [2/4]

void applyBoundaryConditions ( const Real time,
const UInt  block 
)
virtual

applies the b.c. to a specific block

Reimplemented from MonolithicBlock.

Definition at line 238 of file MonolithicBlockMatrix.cpp.

◆ applyBoundaryConditions() [3/4]

void applyBoundaryConditions ( const Real time,
vectorPtr_Type rhs 
)

applies all the b.c.s to the global matrix, updates the rhs

Definition at line 224 of file MonolithicBlockMatrix.cpp.

◆ applyBoundaryConditions() [4/4]

void applyBoundaryConditions ( const Real time,
vectorPtr_Type rhs,
const UInt  block 
)

applies the b.c. to the a specific block, updates the rhs

Definition at line 232 of file MonolithicBlockMatrix.cpp.

◆ addToCoupling() [1/2]

void addToCoupling ( const matrixPtr_Type Mat,
UInt   
)
virtual

adds a block to the coupling matrix

Parameters
Matblock added
positionof the block, unused here because there is only one coupling matrix

Implements MonolithicBlock.

Definition at line 244 of file MonolithicBlockMatrix.cpp.

◆ addToCoupling() [2/2]

void addToCoupling ( const Real entry,
UInt  row,
UInt  col,
UInt  position 
)
virtual

adds an entry to the coupling matrix

Parameters
entryentry
rowrow for the insertion
colcolon for the insertion

Implements MonolithicBlock.

Definition at line 260 of file MonolithicBlockMatrix.cpp.

◆ addToGlobalMatrix()

void addToGlobalMatrix ( const matrixPtr_Type Mat)
inline

adds a block to the coupling matrix

This method is specific for the MonolithicBlockMatrix class, it is used e.g. to add the shape derivatives block in FSI to the global matrix

Parameters
Matblock matrix to be added.

Definition at line 311 of file MonolithicBlockMatrix.hpp.

◆ push_back_coupling()

void push_back_coupling ( matrixPtr_Type coupling)
inlinevirtual

adds a coupling block to the coupling matrix

This method is kept for compatibility with the base class. It calls the method addToCoupling.

Implements MonolithicBlock.

Definition at line 331 of file MonolithicBlockMatrix.hpp.

◆ interface()

UInt interface ( )
inline

returns the dimension of the interface

NOTE: it has to be multiplied times nDimensions to get the number of interface dofs

Definition at line 343 of file MonolithicBlockMatrix.hpp.

◆ interfaceMap()

mapPtr_Type interfaceMap ( ) const
inline

returns the map built for theLagrange multipliers

Definition at line 349 of file MonolithicBlockMatrix.hpp.

◆ coupling()

matrixPtr_Type coupling ( ) const
inline

Definition at line 354 of file MonolithicBlockMatrix.hpp.

◆ numerationInterface()

void numerationInterface ( vectorPtr_Type numeration)
inline

returns the numeration of the interface

Parameters
numerationoutput vector

Definition at line 363 of file MonolithicBlockMatrix.hpp.

◆ whereIsBlock()

const UInt whereIsBlock ( UInt  ) const
inlinevirtual

Implements MonolithicBlock.

Definition at line 368 of file MonolithicBlockMatrix.hpp.

◆ createAdditiveSchwarz()

static MonolithicBlockMatrix* createAdditiveSchwarz ( )
inlinestatic

Definition at line 377 of file MonolithicBlockMatrix.hpp.

Field Documentation

◆ M_globalMatrix

matrixPtr_Type M_globalMatrix
protected

Definition at line 391 of file MonolithicBlockMatrix.hpp.

◆ M_coupling

matrixPtr_Type M_coupling
protected

Definition at line 392 of file MonolithicBlockMatrix.hpp.

◆ M_interfaceMap

mapPtr_Type M_interfaceMap
protected

Definition at line 393 of file MonolithicBlockMatrix.hpp.

◆ M_interface

UInt M_interface
protected

Definition at line 394 of file MonolithicBlockMatrix.hpp.

◆ M_couplingFlags

std::unique_ptr<std::vector<Int> > M_couplingFlags
private

Definition at line 401 of file MonolithicBlockMatrix.hpp.

◆ M_numerationInterface

vectorPtr_Type M_numerationInterface
private

Definition at line 402 of file MonolithicBlockMatrix.hpp.


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