A derived class of linearAlgebra::LinearSolverFunction to encapsulate the Poisson partial differential equation (PDE) discretized in a finite element (FE) basis. The Possion PDE is given as: \(\nabla^2 v(\textbf{r}) = -4 \pi \rho(\textbf{r})$\f
with the boundary condition on
\)@_fakenlv(\textbf{r})|_{\partial \Omega}=g(\textbf{r})$\f ( \(\\partial Omega$\f denoting the boundary of a domain \)\Omega$\f). Here \(v$\f has the physical notion of a potential (e.g.,
Hartree potential, nuclear potential, etc.) arising due to a charge
distributin \)\rho$\f.
More...
|
| | PoissonSolverDealiiMatrixFreeFE (std::shared_ptr< const basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > feBasisManagerField, std::shared_ptr< const basis::FEBasisDataStorage< ValueTypeOperator, memorySpace > > feBasisDataStorageStiffnessMatrix, const std::map< std::string, std::shared_ptr< const basis::FEBasisDataStorage< ValueTypeOperator, memorySpace > > > &feBasisDataStorageRhs, const std::map< std::string, const quadrature::QuadratureValuesContainer< ValueType, memorySpace > & > &inpRhs, const linearAlgebra::PreconditionerType pcType) |
| | This constructor creates an instance of a base LinearSolverFunction called PoissonSolverDealiiMatrixFreeFE. More...
|
| |
| | PoissonSolverDealiiMatrixFreeFE (std::shared_ptr< const basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > feBasisManagerField, std::shared_ptr< const basis::FEBasisDataStorage< ValueTypeOperator, memorySpace > > feBasisDataStorageStiffnessMatrix, std::shared_ptr< const basis::FEBasisDataStorage< ValueTypeOperator, memorySpace > > feBasisDataStorageRhs, const quadrature::QuadratureValuesContainer< ValueType, memorySpace > &inpRhs, const linearAlgebra::PreconditionerType pcType) |
| | This constructor creates an instance of a base LinearSolverFunction called PoissonSolverDealiiMatrixFreeFE. More...
|
| |
| void | reinit (std::shared_ptr< const basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > feBasisManagerField, const std::map< std::string, const quadrature::QuadratureValuesContainer< ValueType, memorySpace > & > &inpRhs) |
| |
| void | reinit (std::shared_ptr< const basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > feBasisManagerField, const quadrature::QuadratureValuesContainer< ValueType, memorySpace > &inpRhs) |
| |
| | ~PoissonSolverDealiiMatrixFreeFE ()=default |
| |
| void | solve (const double absTolerance, const unsigned int maxNumberIterations) |
| |
| void | getSolution (linearAlgebra::MultiVector< ValueType, memorySpace > &solution) |
| |
| const utils::mpi::MPIComm & | getMPIComm () const |
| |
|
| const distributedCPUVec< ValueTypeOperand > & | getRhs () const |
| |
| const distributedCPUVec< ValueType > & | getInitialGuess () const |
| |
| void | setSolution (const distributedCPUVec< ValueType > &x) |
| |
| void | computeRhs (distributedCPUVec< double > &rhs, const std::map< std::string, const quadrature::QuadratureValuesContainer< linearAlgebra::blasLapack::scalar_type< ValueTypeOperator, ValueTypeOperand >, memorySpace > & > &inpRhs) |
| |
| void | vmult (distributedCPUVec< double > &Ax, distributedCPUVec< double > &x) |
| |
| void | precondition_Jacobi (distributedCPUVec< double > &dst, const distributedCPUVec< double > &src) const |
| |
| void | computeDiagonalA () |
| |
| void | AX (const dealii::MatrixFree< dim, double > &matrixFreeData, distributedCPUVec< double > &dst, const distributedCPUVec< double > &src, const std::pair< unsigned int, unsigned int > &cell_range) const |
| |
| void | CGsolve (const double absTolerance, const unsigned int maxNumberIterations, bool distributeFlag) |
| |
|
| size_type | d_numComponents |
| |
| std::shared_ptr< const basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > | d_feBasisManagerField |
| |
| linearAlgebra::PreconditionerType | d_pcType |
| |
| utils::Profiler | d_p |
| |
| std::shared_ptr< basis::FEBasisManager< ValueTypeOperand, ValueTypeOperator, memorySpace, dim > > | d_feBasisManagerHomo |
| |
| distributedCPUVec< ValueTypeOperator > | d_x |
| |
| distributedCPUVec< ValueTypeOperator > | d_rhs |
| |
| distributedCPUVec< ValueTypeOperator > | d_initial |
| |
| distributedCPUVec< ValueTypeOperator > | d_diagonalA |
| |
| std::shared_ptr< dealii::MatrixFree< dim, ValueTypeOperator > > | d_dealiiMatrixFree |
| |
| std::shared_ptr< const dealii::DoFHandler< dim > > | d_dealiiDofHandler |
| |
| const dealii::AffineConstraints< ValueTypeOperand > * | d_dealiiAffineConstraintMatrix |
| |
| const dealii::AffineConstraints< ValueTypeOperand > * | d_constraintsInfo |
| |
| unsigned int | d_num1DQuadPointsStiffnessMatrix |
| |
| std::map< std::string, unsigned int > | d_num1DQuadPointsRhs |
| |
| size_type | d_feOrder |
| |
| unsigned int | d_dofHandlerIndex |
| |
| std::map< std::string, std::shared_ptr< const basis::FEBasisDataStorage< ValueTypeOperator, memorySpace > > > | d_feBasisDataStorageRhs |
| |
| unsigned int | d_matrixFreeQuadCompStiffnessMatrix |
| |
| std::map< dealii::CellId, unsigned int > | d_cellIdToCellIndexMap |
| |
| std::vector< dealii::Quadrature< dim > > | d_dealiiQuadratureRuleVec |
| |
| dealii::MappingQ1< dim, dim > | d_mappingDealii |
| |
| utils::ConditionalOStream | pcout |
| |
template<typename ValueTypeOperator, typename ValueTypeOperand,
utils::MemorySpace memorySpace,
size_type dim>
class dftefe::electrostatics::PoissonSolverDealiiMatrixFreeFE< ValueTypeOperator, ValueTypeOperand, memorySpace, dim >
A derived class of linearAlgebra::LinearSolverFunction to encapsulate the Poisson partial differential equation (PDE) discretized in a finite element (FE) basis. The Possion PDE is given as: \(\nabla^2 v(\textbf{r}) = -4 \pi \rho(\textbf{r})$\f
with the boundary condition on
\)@_fakenlv(\textbf{r})|_{\partial \Omega}=g(\textbf{r})$\f ( \(\\partial Omega$\f denoting the boundary of a domain \)\Omega$\f). Here \(v$\f has the physical notion of a potential (e.g.,
Hartree potential, nuclear potential, etc.) arising due to a charge
distributin \)\rho$\f.
- Template Parameters
-
| ValueTypeOperator | The datatype (float, double, complex<double>, etc.) for the underlying operator |
| ValueTypeOperand | The datatype (float, double, complex<double>, etc.) of the vector, matrices, etc. on which the operator will act |
| memorySpace | The meory sapce (HOST, DEVICE, HOST_PINNES, etc.) in which the data of the operator and its operands reside |
| dim | Dimension of the Poisson problem |