18 #ifndef _FCST_LINEAR_SOLVERS_H_
19 #define _FCST_LINEAR_SOLVERS_H_
21 #include <lac/vector.h>
22 #include <lac/block_vector.h>
24 #include <lac/sparse_matrix.h>
25 #include <lac/block_sparse_matrix.h>
27 #include <lac/block_matrix_array.h>
28 #include <lac/vector_memory.h>
30 #include <lac/sparse_direct.h>
31 #include <lac/solver_gmres.h>
33 #include <lac/precondition.h>
34 #include <lac/sparse_ilu.h>
36 using namespace dealii;
57 namespace LinearSolvers
106 template<
typename MATRIX,
typename VECTOR>
109 const VECTOR& right_hand_side)
const
112 tmp = right_hand_side;
114 SparseDirectUMFPACK solver;
115 solver.initialize(matrix);
166 template<
typename MATRIX,
typename VECTOR,
typename PRECONDITIONER>
169 const VECTOR& right_hand_side,
170 const unsigned int& n_iter,
171 const double& tolerance,
172 const PRECONDITIONER& preconditioner)
const
174 SolverControl solver_control(n_iter, tolerance);
176 typename SolverGMRES< VECTOR >::AdditionalData additional_data(100,
false,
true);
178 SolverGMRES< VECTOR > solver(solver_control, additional_data);
198 template<
typename MATRIX,
typename VECTOR,
typename PRECONDITIONER>
199 void solve(SolverControl& solver_control,
200 const MATRIX& matrix,
202 const VECTOR& right_hand_side,
203 const PRECONDITIONER& preconditioner)
const
205 SolverGMRES< VECTOR > solver(solver_control);
237 preconditioner_matrix(matrix.n_block_rows()),
238 preconditioner_pointer(matrix.n_block_rows()),
239 preconditioner(matrix.n_block_rows())
243 for(
unsigned int i = 0; i < matrix.n_block_rows(); ++i)
245 preconditioner_matrix[i].initialize(matrix.block(i,i),
246 SparseILU<double>::AdditionalData());
247 preconditioner_pointer[i].set_memory(&
vector_pool);
248 preconditioner_pointer[i] = &preconditioner_matrix[i];
251 for(
unsigned int i = 0; i < matrix.n_block_rows(); ++i)
252 for(
unsigned int j = 0; j < matrix.n_block_cols(); ++j)
255 preconditioner.enter(preconditioner_pointer[i], i, i);
257 preconditioner.enter(matrix.block(i,j), i, j);
This class implements ILU preconditioner.
Definition: linear_solvers.h:225
SparseDirectUMFPACKSolver()
Constructor.
Definition: linear_solvers.h:85
~GMRESSolver()
Destructor.
Definition: linear_solvers.h:148
void solve(const MATRIX &matrix, VECTOR &solution, const VECTOR &right_hand_side) const
This function solves the linear system Ax = b where.
Definition: linear_solvers.h:107
This class implements GMRES solver.
Definition: linear_solvers.h:133
std::vector< PointerMatrixAux< SparseILU< double >, Vector< double > > > preconditioner_pointer
Preconditioner pointer.
Definition: linear_solvers.h:285
static GrowingVectorMemory< BlockVector< double > > block_vector_pool
Definition: linear_solvers.h:64
void solve(SolverControl &solver_control, const MATRIX &matrix, VECTOR &solution, const VECTOR &right_hand_side, const PRECONDITIONER &preconditioner) const
This routine is used when solver_control is an object already initialized.
Definition: linear_solvers.h:199
FCSTLogStream log
Object used to output data to file and, if file attached recorded to a file as well.
BlockTrianglePrecondition< double > preconditioner
Preconditioner.
Definition: linear_solvers.h:274
This class implements an interface to the sparse direct solver UMFPACK, see the link below http://www...
Definition: linear_solvers.h:75
~ILUPreconditioner()
Destructor.
Definition: linear_solvers.h:264
void solve(const MATRIX &matrix, VECTOR &solution, const VECTOR &right_hand_side, const unsigned int &n_iter, const double &tolerance, const PRECONDITIONER &preconditioner) const
This function solves the linear system Ax = b where.
Definition: linear_solvers.h:167
ILUPreconditioner(const BlockSparseMatrix< double > &matrix)
Constructor.
Definition: linear_solvers.h:235
static GrowingVectorMemory< Vector< double > > vector_pool
Memory management objects.
Definition: linear_solvers.h:63
std::vector< SparseILU< double > > preconditioner_matrix
Preconditioner matrix.
Definition: linear_solvers.h:280
GMRESSolver()
Constructor.
Definition: linear_solvers.h:143
~SparseDirectUMFPACKSolver()
Destructor.
Definition: linear_solvers.h:90