MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM > Class Template Reference

#include <MatrixBasedBidomainSolver.hpp>

Inheritance diagram for MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >:

Inheritance graph
[legend]
Collaboration diagram for MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >:

Collaboration graph
[legend]

List of all members.

Public Member Functions

 MatrixBasedBidomainSolver (bool bathSimulation, AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > *pMesh, BidomainTissue< SPACE_DIM > *pTissue, BoundaryConditionsContainer< ELEMENT_DIM, SPACE_DIM, 2 > *pBoundaryConditions, unsigned numQuadPoints=2)

Private Member Functions

void InitialiseForSolve (Vec initialSolution)
void SetupLinearSystem (Vec currentSolution, bool computeMatrix)

Private Attributes

Mat mMassMatrix
Vec mVecForConstructingRhs
BidomainCorrectionTermAssembler
< ELEMENT_DIM, SPACE_DIM > * 
mpBidomainCorrectionTermAssembler


Detailed Description

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
class MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >

A better Bidomain solver (better than BasicModomainSolver, from which it inherits), which computes the right-hand-side (RHS) vector of the linear system to be solved using matrix-vector products, rather than assembly. Massively more efficient than BasicBidomainSolver

Definition at line 55 of file MatrixBasedBidomainSolver.hpp.


Constructor & Destructor Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::MatrixBasedBidomainSolver ( bool  bathSimulation,
AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > *  pMesh,
BidomainTissue< SPACE_DIM > *  pTissue,
BoundaryConditionsContainer< ELEMENT_DIM, SPACE_DIM, 2 > *  pBoundaryConditions,
unsigned  numQuadPoints = 2 
) [inline]

Constructor

Parameters:
bathSimulation Whether the simulation involves a perfusing bath
pMesh pointer to the mesh
pTissue pointer to the tissue
pBoundaryConditions pointer to the boundary conditions
numQuadPoints number of quadrature points (defaults to 2)

Definition at line 218 of file MatrixBasedBidomainSolver.cpp.

References HeartConfig::Instance(), AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mNumQuadPoints, MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainCorrectionTermAssembler, AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainTissue, AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >::mpMesh, MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mVecForConstructingRhs, and AbstractCardiacTissue< SPACE_DIM >::SetCacheReplication().


Member Function Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
void MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::InitialiseForSolve ( Vec  initialSolution  )  [inline, private, virtual]

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
void MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::SetupLinearSystem ( Vec  currentSolution,
bool  computeMatrix 
) [inline, private, virtual]

Implementation of SetupLinearSystem() which uses the assembler to compute the LHS matrix, but sets up the RHS vector using the mass-matrix (constructed using a separate assembler) multiplied by a vector

Parameters:
currentSolution Solution at current time
computeMatrix Whether to compute the matrix of the linear system

Implements AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >.

Definition at line 59 of file MatrixBasedBidomainSolver.cpp.

References BoundaryConditionsContainer< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >::ApplyDirichletToLinearProblem(), AbstractFeObjectAssembler< DIM, DIM, 2, false, true, CARDIAC >::Assemble(), PetscMatTools::AssembleFinal(), LinearSystem::AssembleFinalLhsMatrix(), LinearSystem::AssembleIntermediateLhsMatrix(), LinearSystem::AssembleRhsVector(), AbstractFeObjectAssembler< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM, CAN_ASSEMBLE_VECTOR, CAN_ASSEMBLE_MATRIX, INTERPOLATION_LEVEL >::AssembleVector(), DistributedVector::Begin(), GenericEventHandler< 13, HeartEventHandler >::BeginEvent(), DistributedVectorFactory::CreateDistributedVector(), DistributedVector::End(), GenericEventHandler< 13, HeartEventHandler >::EndEvent(), AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::FinaliseForBath(), HeartConfig::GetCapacitance(), PdeSimulationTime::GetPdeTimeStepInverse(), HeartConfig::GetSurfaceAreaToVolumeRatio(), HeartConfig::Instance(), HeartRegionCode::IsRegionBath(), AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mBathSimulation, MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mMassMatrix, AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mNumQuadPoints, AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainAssembler, MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainCorrectionTermAssembler, AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainTissue, AbstractBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBoundaryConditions, AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >::mpLinearSystem, AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >::mpMesh, MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mVecForConstructingRhs, DistributedVector::Restore(), LinearSystem::rGetLhsMatrix(), LinearSystem::rGetRhsVector(), AbstractFeObjectAssembler< DIM, DIM, 2, false, true, CARDIAC >::SetMatrixToAssemble(), and AbstractFeObjectAssembler< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM, CAN_ASSEMBLE_VECTOR, CAN_ASSEMBLE_MATRIX, INTERPOLATION_LEVEL >::SetVectorToAssemble().


Member Data Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
Mat MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mMassMatrix [private]

Mass matrix, used to computing the RHS vector (actually: mass-matrix in voltage-voltage block, zero elsewhere

Definition at line 61 of file MatrixBasedBidomainSolver.hpp.

Referenced by MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::InitialiseForSolve(), and MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::SetupLinearSystem().

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
Vec MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mVecForConstructingRhs [private]

The vector multiplied by the mass matrix. Ie, if the linear system to be solved is Ax=b, this vector is z where b=Mz.

Definition at line 67 of file MatrixBasedBidomainSolver.hpp.

Referenced by MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::InitialiseForSolve(), MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::MatrixBasedBidomainSolver(), and MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::SetupLinearSystem().

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
BidomainCorrectionTermAssembler<ELEMENT_DIM,SPACE_DIM>* MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::mpBidomainCorrectionTermAssembler [private]

If using state variable interpolation, points to an assembler to use in computing the correction term to apply to the RHS.

Definition at line 73 of file MatrixBasedBidomainSolver.hpp.

Referenced by MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::MatrixBasedBidomainSolver(), and MatrixBasedBidomainSolver< ELEMENT_DIM, SPACE_DIM >::SetupLinearSystem().


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

Generated on Tue May 31 14:33:39 2011 for Chaste by  doxygen 1.5.5