From: Daniel Arndt Date: Sat, 2 Sep 2017 13:48:40 +0000 (+0200) Subject: Add a slot for notification about re-orthogonalization in SolverGMRES X-Git-Tag: v9.0.0-rc1~1130^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=846b6762178edc0200e8d80103dc5b18dada4a9c;p=dealii.git Add a slot for notification about re-orthogonalization in SolverGMRES --- diff --git a/include/deal.II/lac/solver_gmres.h b/include/deal.II/lac/solver_gmres.h index 08d51fcdde..4920926e28 100644 --- a/include/deal.II/lac/solver_gmres.h +++ b/include/deal.II/lac/solver_gmres.h @@ -292,6 +292,14 @@ public: const std::function &)> &slot); + /** + * Connect a slot to retrieve a notification when the vectors are + * re-orthogonalized. + */ + boost::signals2::connection + connect_re_orthogonalization_slot(const std::function &slot); + + DeclException1 (ExcTooFewTmpVectors, int, << "The number of temporary vectors you gave (" @@ -347,6 +355,12 @@ protected: */ boost::signals2::signal &)> krylov_space_signal; + /** + * Signal used to retrieve a notification + * when the vectors are re-orthogonalized. + */ + boost::signals2::signal re_orthogonalize_signal; + /** * Implementation of the computation of the norm of the residual. */ @@ -368,14 +382,17 @@ protected: * should be applied twice. The algorithm checks loss of orthogonality in * the procedure every fifth step and sets the flag to true in that case. * All subsequent iterations use re-orthogonalization. + * Calls the signal re_orthogonalize_signal if it is connected. */ static double - modified_gram_schmidt (const internal::SolverGMRES::TmpVectors &orthogonal_vectors, - const unsigned int dim, - const unsigned int accumulated_iterations, - VectorType &vv, - Vector &h, - bool &re_orthogonalize); + modified_gram_schmidt + (const internal::SolverGMRES::TmpVectors &orthogonal_vectors, + const unsigned int dim, + const unsigned int accumulated_iterations, + VectorType &vv, + Vector &h, + bool &re_orthogonalize, + const boost::signals2::signal &re_orthogonalize_signal = boost::signals2::signal()); /** * Estimates the eigenvalues from the Hessenberg matrix, H_orig, generated @@ -642,14 +659,17 @@ SolverGMRES::modified_gram_schmidt const unsigned int accumulated_iterations, VectorType &vv, Vector &h, - bool &re_orthogonalize) + bool &reorthogonalize, + const boost::signals2::signal &reorthogonalize_signal) { Assert(dim > 0, ExcInternalError()); const unsigned int inner_iteration = dim - 1; // need initial norm for detection of re-orthogonalization, see below double norm_vv_start = 0; - if (re_orthogonalize == false && inner_iteration % 5 == 4) + const bool consider_reorthogonalize = (reorthogonalize == false) && + (inner_iteration % 5 == 4); + if (consider_reorthogonalize) norm_vv_start = vv.l2_norm(); // Orthogonalization @@ -665,7 +685,7 @@ SolverGMRES::modified_gram_schmidt // orthogonalization. If vv became very small (here: less than the square // root of the machine precision times 10), it is almost in the span of the // previous vectors, which indicates loss of precision. - if (re_orthogonalize == false && inner_iteration % 5 == 4) + if (consider_reorthogonalize) { if (norm_vv > 10. * norm_vv_start * std::sqrt(std::numeric_limits::epsilon())) @@ -673,13 +693,13 @@ SolverGMRES::modified_gram_schmidt else { - re_orthogonalize = true; - deallog << "Re-orthogonalization enabled at step " - << accumulated_iterations << std::endl; + reorthogonalize = true; + if (!reorthogonalize_signal.empty()) + reorthogonalize_signal(accumulated_iterations); } } - if (re_orthogonalize == true) + if (reorthogonalize == true) { double htmp = vv * orthogonal_vectors[0]; h(0) += htmp; @@ -913,7 +933,8 @@ SolverGMRES::solve (const MatrixType &A, const double s = modified_gram_schmidt(tmp_vectors, dim, accumulated_iterations, - vv, h, re_orthogonalize); + vv, h, re_orthogonalize, + re_orthogonalize_signal); h(inner_iteration+1) = s; //s=0 is a lucky breakdown, the solver will reach convergence, @@ -1099,6 +1120,16 @@ SolverGMRES::connect_krylov_space_slot +template +boost::signals2::connection +SolverGMRES::connect_re_orthogonalization_slot +(const std::function &slot) +{ + return re_orthogonalize_signal.connect(slot); +} + + + template double SolverGMRES::criterion ()