16 #ifndef dealii_solver_cg_h 17 #define dealii_solver_cg_h 20 #include <deal.II/base/config.h> 22 #include <deal.II/base/exceptions.h> 23 #include <deal.II/base/logstream.h> 24 #include <deal.II/base/subscriptor.h> 26 #include <deal.II/lac/solver.h> 27 #include <deal.II/lac/solver_control.h> 28 #include <deal.II/lac/tridiagonal_matrix.h> 32 DEAL_II_NAMESPACE_OPEN
95 template <
typename VectorType = Vector<
double>>
133 template <
typename MatrixType,
typename PreconditionerType>
135 solve(
const MatrixType & A,
137 const VectorType & b,
138 const PreconditionerType &preconditioner);
146 boost::signals2::connection
148 const std::function<
void(
typename VectorType::value_type,
149 typename VectorType::value_type)> &slot);
157 boost::signals2::connection
159 const bool every_iteration =
false);
167 boost::signals2::connection
169 const std::function<
void(
const std::vector<double> &)> &slot,
170 const bool every_iteration =
false);
180 const VectorType & x,
181 const VectorType & r,
182 const VectorType & d)
const;
191 const std::vector<typename VectorType::value_type> &diagonal,
192 const std::vector<typename VectorType::value_type> &offdiagonal,
193 const boost::signals2::signal<
void(
const std::vector<double> &)>
195 const boost::signals2::signal<
void(
double)> &cond_signal);
205 boost::signals2::signal<void(
typename VectorType::value_type,
206 typename VectorType::value_type)>
231 boost::signals2::signal<void(const std::vector<double> &)>
241 template <
typename VectorType>
245 :
Solver<VectorType>(cn, mem)
251 template <
typename VectorType>
259 template <
typename VectorType>
264 const VectorType &)
const 269 template <
typename VectorType>
272 const std::vector<typename VectorType::value_type> &diagonal,
273 const std::vector<typename VectorType::value_type> &offdiagonal,
274 const boost::signals2::signal<
void(
const std::vector<double> &)>
276 const boost::signals2::signal<
void(
double)> &cond_signal)
283 for (
size_type i = 0; i < diagonal.size(); ++i)
285 T(i, i) = diagonal[i];
286 if (i < diagonal.size() - 1)
287 T(i, i + 1) = offdiagonal[i];
289 T.compute_eigenvalues();
291 if (diagonal.size() > 1)
293 auto condition_number = T.eigenvalue(T.n() - 1) / T.eigenvalue(0);
296 cond_signal(std::abs(condition_number));
303 for (
unsigned int j = 0; j < T.n(); ++j)
316 template <
typename VectorType>
317 template <
typename MatrixType,
typename PreconditionerType>
321 const VectorType & b,
322 const PreconditionerType &preconditioner)
324 using number =
typename VectorType::value_type;
336 VectorType &g = *g_pointer;
337 VectorType &d = *d_pointer;
338 VectorType &h = *h_pointer;
341 const bool do_eigenvalues =
347 std::vector<typename VectorType::value_type> diagonal;
348 std::vector<typename VectorType::value_type> offdiagonal;
351 double res = -std::numeric_limits<double>::max();
353 typename VectorType::value_type eigen_beta_alpha = 0;
380 if (std::is_same<PreconditionerType, PreconditionIdentity>::value ==
false)
382 preconditioner.vmult(h, g);
399 number alpha = d * h;
404 res = std::sqrt(std::abs(g.add_and_dot(alpha, h, g)));
412 if (std::is_same<PreconditionerType, PreconditionIdentity>::value ==
415 preconditioner.vmult(h, g);
421 d.sadd(beta, -1., h);
428 d.sadd(beta, -1., g);
438 diagonal.push_back(number(1.) / alpha + eigen_beta_alpha);
439 eigen_beta_alpha = beta / alpha;
440 offdiagonal.push_back(std::sqrt(beta) / alpha);
461 template <
typename VectorType>
462 boost::signals2::connection
464 const std::function<
void(
typename VectorType::value_type,
465 typename VectorType::value_type)> &slot)
472 template <
typename VectorType>
473 boost::signals2::connection
475 const std::function<
void(
double)> &slot,
476 const bool every_iteration)
490 template <
typename VectorType>
491 boost::signals2::connection
493 const std::function<
void(
const std::vector<double> &)> &slot,
494 const bool every_iteration)
510 DEAL_II_NAMESPACE_CLOSE
virtual ~SolverCG() override=default
types::global_dof_index size_type
SolverCG(SolverControl &cn, VectorMemory< VectorType > &mem, const AdditionalData &data=AdditionalData())
boost::signals2::signal< void(const std::vector< double > &)> eigenvalues_signal
std::array< Number, 1 > eigenvalues(const SymmetricTensor< 2, 1, Number > &T)
boost::signals2::signal< void(double)> all_condition_numbers_signal
#define AssertThrow(cond, exc)
static::ExceptionBase & ExcDivideByZero()
unsigned long long int global_dof_index
boost::signals2::signal< void(double)> condition_number_signal
boost::signals2::connection connect_condition_number_slot(const std::function< void(double)> &slot, const bool every_iteration=false)
virtual void print_vectors(const unsigned int step, const VectorType &x, const VectorType &r, const VectorType &d) const
Stop iteration, goal reached.
#define Assert(cond, exc)
AdditionalData additional_data
boost::signals2::signal< SolverControl::State(const unsigned int iteration, const double check_value, const VectorType ¤t_iterate), StateCombiner > iteration_status
boost::signals2::signal< void(typename VectorType::value_type, typename VectorType::value_type)> coefficients_signal
VectorMemory< VectorType > & memory
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const PreconditionerType &preconditioner)
boost::signals2::signal< void(const std::vector< double > &)> all_eigenvalues_signal
boost::signals2::connection connect_eigenvalues_slot(const std::function< void(const std::vector< double > &)> &slot, const bool every_iteration=false)
static void compute_eigs_and_cond(const std::vector< typename VectorType::value_type > &diagonal, const std::vector< typename VectorType::value_type > &offdiagonal, const boost::signals2::signal< void(const std::vector< double > &)> &eigenvalues_signal, const boost::signals2::signal< void(double)> &cond_signal)
boost::signals2::connection connect_coefficients_slot(const std::function< void(typename VectorType::value_type, typename VectorType::value_type)> &slot)