SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM > Class Template Reference

#include <SimpleLinearParabolicSolver.hpp>

Inherits AbstractAssemblerSolverHybrid< ELEMENT_DIM, SPACE_DIM, 1, NORMAL >, and AbstractDynamicLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, 1 >.

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

List of all members.

Public Member Functions

 SimpleLinearParabolicSolver (AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > *pMesh, AbstractLinearParabolicPde< ELEMENT_DIM, SPACE_DIM > *pPde, BoundaryConditionsContainer< ELEMENT_DIM, SPACE_DIM, 1 > *pBoundaryConditions, unsigned numQuadPoints=2)

Protected Member Functions

virtual c_matrix< double,
1 *(ELEMENT_DIM+1),
1 *(ELEMENT_DIM+1)> 
ComputeMatrixTerm (c_vector< double, ELEMENT_DIM+1 > &rPhi, c_matrix< double, SPACE_DIM, ELEMENT_DIM+1 > &rGradPhi, ChastePoint< SPACE_DIM > &rX, c_vector< double, 1 > &rU, c_matrix< double, 1, SPACE_DIM > &rGradU, Element< ELEMENT_DIM, SPACE_DIM > *pElement)
virtual c_vector< double,
1 *(ELEMENT_DIM+1)> 
ComputeVectorTerm (c_vector< double, ELEMENT_DIM+1 > &rPhi, c_matrix< double, SPACE_DIM, ELEMENT_DIM+1 > &rGradPhi, ChastePoint< SPACE_DIM > &rX, c_vector< double, 1 > &rU, c_matrix< double, 1, SPACE_DIM > &rGradU, Element< ELEMENT_DIM, SPACE_DIM > *pElement)
void SetupLinearSystem (Vec currentSolution, bool computeMatrix)

Protected Attributes

AbstractLinearParabolicPde
< ELEMENT_DIM, SPACE_DIM > * 
mpParabolicPde

Detailed Description

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

SimpleLinearParabolicSolver.

Solver for solving AbstractLinearParabolicPdes

Definition at line 43 of file SimpleLinearParabolicSolver.hpp.


Constructor & Destructor Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::SimpleLinearParabolicSolver ( AbstractTetrahedralMesh< ELEMENT_DIM, SPACE_DIM > *  pMesh,
AbstractLinearParabolicPde< ELEMENT_DIM, SPACE_DIM > *  pPde,
BoundaryConditionsContainer< ELEMENT_DIM, SPACE_DIM, 1 > *  pBoundaryConditions,
unsigned  numQuadPoints = 2 
) [inline]

Constructor.

Parameters:
pMesh pointer to the mesh
pPde pointer to the PDE
pBoundaryConditions pointer to the boundary conditions
numQuadPoints number of quadrature points in each dimension to use per element (defaults to 2)

Definition at line 61 of file SimpleLinearParabolicSolver.cpp.

References AbstractDynamicLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, 1 >::mMatrixIsConstant, and SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::mpParabolicPde.


Member Function Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
c_matrix< double, 1 *(ELEMENT_DIM+1), 1 *(ELEMENT_DIM+1)> SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::ComputeMatrixTerm ( c_vector< double, ELEMENT_DIM+1 > &  rPhi,
c_matrix< double, SPACE_DIM, ELEMENT_DIM+1 > &  rGradPhi,
ChastePoint< SPACE_DIM > &  rX,
c_vector< double, 1 > &  rU,
c_matrix< double, 1, SPACE_DIM > &  rGradU,
Element< ELEMENT_DIM, SPACE_DIM > *  pElement 
) [inline, protected, virtual]

The term to be added to the element stiffness matrix - see AbstractFeVolumeIntegralAssembler

grad_phi[row] . ( pde_diffusion_term * grad_phi[col]) + (1.0/mDt) * pPde->ComputeDuDtCoefficientFunction(rX) * rPhi[row] * rPhi[col]

Parameters:
rPhi The basis functions, rPhi(i) = phi_i, i=1..numBases
rGradPhi Basis gradients, rGradPhi(i,j) = d(phi_j)/d(X_i)
rX The point in space
rU The unknown as a vector, u(i) = u_i
rGradU The gradient of the unknown as a matrix, rGradU(i,j) = d(u_i)/d(X_j)
pElement Pointer to the element

Definition at line 32 of file SimpleLinearParabolicSolver.cpp.

References PdeSimulationTime::GetPdeTimeStepInverse(), and SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::mpParabolicPde.

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
c_vector< double, 1 *(ELEMENT_DIM+1)> SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::ComputeVectorTerm ( c_vector< double, ELEMENT_DIM+1 > &  rPhi,
c_matrix< double, SPACE_DIM, ELEMENT_DIM+1 > &  rGradPhi,
ChastePoint< SPACE_DIM > &  rX,
c_vector< double, 1 > &  rU,
c_matrix< double, 1, SPACE_DIM > &  rGradU,
Element< ELEMENT_DIM, SPACE_DIM > *  pElement 
) [inline, protected, virtual]

The term to be added to the element stiffness vector - see AbstractFeVolumeIntegralAssembler

Parameters:
rPhi The basis functions, rPhi(i) = phi_i, i=1..numBases
rGradPhi Basis gradients, rGradPhi(i,j) = d(phi_j)/d(X_i)
rX The point in space
rU The unknown as a vector, u(i) = u_i
rGradU The gradient of the unknown as a matrix, rGradU(i,j) = d(u_i)/d(X_j)
pElement Pointer to the element

Definition at line 47 of file SimpleLinearParabolicSolver.cpp.

References PdeSimulationTime::GetPdeTimeStepInverse(), and SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::mpParabolicPde.

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

Delegate to AbstractAssemblerSolverHybrid::SetupGivenLinearSystem.

Parameters:
currentSolution The current solution which can be used in setting up the linear system if needed (NULL if there isn't a current solution)
computeMatrix Whether to compute the LHS matrix of the linear system (mainly for dynamic solves).

Implements AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >.

Definition at line 105 of file SimpleLinearParabolicSolver.hpp.

References AbstractLinearPdeSolver< ELEMENT_DIM, SPACE_DIM, PROBLEM_DIM >::mpLinearSystem, and AbstractAssemblerSolverHybrid< ELEMENT_DIM, SPACE_DIM, 1, NORMAL >::SetupGivenLinearSystem().


Member Data Documentation

template<unsigned ELEMENT_DIM, unsigned SPACE_DIM>
AbstractLinearParabolicPde<ELEMENT_DIM,SPACE_DIM>* SimpleLinearParabolicSolver< ELEMENT_DIM, SPACE_DIM >::mpParabolicPde [protected]

The documentation for this class was generated from the following files:
Generated on Thu Dec 22 13:07:42 2011 for Chaste by  doxygen 1.6.3