24 #ifndef _DIRECT_LINEAR_SOLVER_H_
25 #define _DIRECT_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>
34 #include <lac/sparse_direct.h>
36 #include <dofs/dof_tools.h>
38 #include <numerics/vector_tools.h>
59 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
72 template<
typename PROBLEM>
93 template<
typename PROBLEM,
typename INTEGRATOR>
94 void Solve(PROBLEM& pde, INTEGRATOR& integr, VECTOR &rhs, VECTOR &solution,
bool force_matrix_build=
false);
99 SPARSITYPATTERN _sparsity_pattern;
102 dealii::SparseDirectUMFPACK* _A_direct;
104 double _linear_global_tol, _linear_tol;
110 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
117 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
126 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
129 if(_A_direct != NULL)
137 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
138 template<
typename PROBLEM>
142 pde.ComputeSparsityPattern(_sparsity_pattern);
143 _matrix.reinit(_sparsity_pattern);
145 if(_A_direct != NULL)
154 template <
typename SPARSITYPATTERN,
typename MATRIX,
typename VECTOR>
155 template<
typename PROBLEM,
typename INTEGRATOR>
160 bool force_matrix_build)
162 if(force_matrix_build)
164 integr.ComputeMatrix (pde,_matrix);
167 integr.ApplyNewtonBoundaryValues(pde,_matrix,rhs,solution);
169 if(_A_direct == NULL)
171 _A_direct =
new dealii::SparseDirectUMFPACK;
172 _A_direct->initialize(_matrix);
174 else if(force_matrix_build)
176 _A_direct->factorize(_matrix);
179 dealii::Vector<double> sol;
181 _A_direct->solve(sol);
184 pde.GetDoFConstraints().distribute(solution);
Definition: parameterreader.h:36
~DirectLinearSolverWithMatrix()
Definition: directlinearsolver.h:127
static void declare_params(ParameterReader ¶m_reader)
Definition: directlinearsolver.h:111
void ReInit(PROBLEM &pde)
Definition: directlinearsolver.h:139
DirectLinearSolverWithMatrix(ParameterReader ¶m_reader)
Definition: directlinearsolver.h:118
void Solve(PROBLEM &pde, INTEGRATOR &integr, VECTOR &rhs, VECTOR &solution, bool force_matrix_build=false)
Definition: directlinearsolver.h:156
Definition: directlinearsolver.h:60