const std::function<void (const internal::SolverGMRES::TmpVectors<VectorType> &)> &slot);
+ /**
+ * Connect a slot to retrieve a notification when the vectors are
+ * re-orthogonalized.
+ */
+ boost::signals2::connection
+ connect_re_orthogonalization_slot(const std::function<void (int)> &slot);
+
+
DeclException1 (ExcTooFewTmpVectors,
int,
<< "The number of temporary vectors you gave ("
*/
boost::signals2::signal<void (const internal::SolverGMRES::TmpVectors<VectorType> &)> krylov_space_signal;
+ /**
+ * Signal used to retrieve a notification
+ * when the vectors are re-orthogonalized.
+ */
+ boost::signals2::signal<void (int)> re_orthogonalize_signal;
+
/**
* Implementation of the computation of the norm of the residual.
*/
* 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<VectorType> &orthogonal_vectors,
- const unsigned int dim,
- const unsigned int accumulated_iterations,
- VectorType &vv,
- Vector<double> &h,
- bool &re_orthogonalize);
+ modified_gram_schmidt
+ (const internal::SolverGMRES::TmpVectors<VectorType> &orthogonal_vectors,
+ const unsigned int dim,
+ const unsigned int accumulated_iterations,
+ VectorType &vv,
+ Vector<double> &h,
+ bool &re_orthogonalize,
+ const boost::signals2::signal<void(int)> &re_orthogonalize_signal = boost::signals2::signal<void(int)>());
/**
* Estimates the eigenvalues from the Hessenberg matrix, H_orig, generated
const unsigned int accumulated_iterations,
VectorType &vv,
Vector<double> &h,
- bool &re_orthogonalize)
+ bool &reorthogonalize,
+ const boost::signals2::signal<void(int)> &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
// 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<typename VectorType::value_type>::epsilon()))
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;
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,
+template <class VectorType>
+boost::signals2::connection
+SolverGMRES<VectorType>::connect_re_orthogonalization_slot
+(const std::function<void (int)> &slot)
+{
+ return re_orthogonalize_signal.connect(slot);
+}
+
+
+
template <class VectorType>
double
SolverGMRES<VectorType>::criterion ()