24 #ifndef _CG_LINEAR_SOLVER_H_
25 #define _CG_LINEAR_SOLVER_H_
27 #include <lac/vector.h>
28 #include <lac/block_sparsity_pattern.h>
29 #include <lac/block_sparse_matrix.h>
30 #include <lac/compressed_simple_sparsity_pattern.h>
31 #include <lac/solver_cg.h>
32 #include <lac/precondition.h>
33 #include <lac/full_matrix.h>
35 #include <dofs/dof_tools.h>
37 #include <numerics/vector_tools.h>
56 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
69 template<
typename PROBLEM>
90 template<
typename PROBLEM,
typename INTEGRATOR>
91 void Solve(PROBLEM& pde, INTEGRATOR& integr, VECTOR &rhs, VECTOR &solution,
bool force_matrix_build=
false);
96 SPARSITYPATTERN _sparsity_pattern;
99 double _linear_global_tol, _linear_tol;
105 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
108 param_reader.
SetSubsection(
"cglinearsolver_withmatrix parameters");
109 param_reader.
declare_entry(
"linear_global_tol",
"1.e-16",Patterns::Double(0),
"global tolerance for the cg iteration");
110 param_reader.
declare_entry(
"linear_tol",
"1.e-12",Patterns::Double(0),
"relative tolerance for the cg iteration");
111 param_reader.
declare_entry(
"linear_maxiter",
"1000",Patterns::Integer(0),
"maximal number of cg steps");
115 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
119 param_reader.
SetSubsection(
"cglinearsolver_withmatrix parameters");
120 _linear_global_tol = param_reader.
get_double (
"linear_global_tol");
121 _linear_tol = param_reader.
get_double (
"linear_tol");
122 _linear_maxiter = param_reader.
get_integer (
"linear_maxiter");
128 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
135 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
136 template<
typename PROBLEM>
140 pde.ComputeSparsityPattern(_sparsity_pattern);
141 _matrix.reinit(_sparsity_pattern);
145 template <
typename PRECONDITIONER,
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
146 template<
typename PROBLEM,
typename INTEGRATOR>
151 bool force_matrix_build)
153 if(force_matrix_build)
155 integr.ComputeMatrix (pde,_matrix);
158 integr.ApplyNewtonBoundaryValues(pde,_matrix,rhs,solution);
160 dealii::SolverControl solver_control (_linear_maxiter, _linear_global_tol,
false,
false);
161 dealii::SolverCG<VECTOR> cg (solver_control);
162 PRECONDITIONER precondition;
163 precondition.initialize(_matrix);
164 cg.solve (_matrix, solution, rhs,
167 pde.GetDoFConstraints().distribute(solution);
static void declare_params(ParameterReader ¶m_reader)
Definition: cglinearsolver.h:106
double get_double(const std::string &entry_name)
Definition: parameterreader.h:115
Definition: parameterreader.h:36
void ReInit(PROBLEM &pde)
Definition: cglinearsolver.h:137
void declare_entry(const std::string &entry, const std::string &default_value, const Patterns::PatternBase &pattern=Patterns::Anything(), const std::string &documentation=std::string())
Definition: parameterreader.h:98
CGLinearSolverWithMatrix(ParameterReader ¶m_reader)
Definition: cglinearsolver.h:117
Definition: cglinearsolver.h:57
int get_integer(const std::string &entry_name)
Definition: parameterreader.h:126
void SetSubsection(const std::string subsection)
Definition: parameterreader.h:93
~CGLinearSolverWithMatrix()
Definition: cglinearsolver.h:129
void Solve(PROBLEM &pde, INTEGRATOR &integr, VECTOR &rhs, VECTOR &solution, bool force_matrix_build=false)
Definition: cglinearsolver.h:147