]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce IterationNumberControl to perform a fixed number of solver steps without...
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 20 Nov 2013 09:54:38 +0000 (09:54 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 20 Nov 2013 09:54:38 +0000 (09:54 +0000)
git-svn-id: https://svn.dealii.org/trunk@31722 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/solver_control.h
deal.II/source/lac/solver_control.cc
tests/lac/iteration_number_control.cc [new file with mode: 0644]
tests/lac/iteration_number_control.output [new file with mode: 0644]

index f303cddcf82e5f1321146903ab94158b1d35dd48..80d998fd86700e83846f8415d103795311a6d432 100644 (file)
@@ -462,6 +462,58 @@ protected:
   double reduced_tol;
 };
 
+/**
+ * Specialization of @p SolverControl which returns @p success if a given
+ * number of iteration was performed, irrespective of the actual
+ * residual. This is useful in cases where you don't want to solve exactly,
+ * but rather want to perform a fixed number of iterations, e.g. in an inner
+ * solver. The arguments given to this class are exactly the same as for the
+ * SolverControl class and the solver terminates similarly when one of the
+ * given tolerance or the maximum iteration count were reached. The only
+ * difference to SolverControl is that the solver returns success in the
+ * latter case.
+ *
+ * @author Martin Kronbichler
+ */
+class IterationNumberControl : public SolverControl
+{
+public:
+  /**
+   * Constructor.  Provide exactly the same arguments as the constructor of
+   * the SolverControl class.
+   */
+  IterationNumberControl (const unsigned int maxiter = 100,
+                          const double       tolerance = 1e-12,
+                          const bool     log_history = false,
+                          const bool     log_result  = true);
+
+  /**
+   * Initialize with a SolverControl object. The result will emulate
+   * SolverControl by setting #reduce to zero.
+   */
+  IterationNumberControl (const SolverControl &c);
+
+  /**
+   * Assign a SolverControl object to ReductionControl. The result of the
+   * assignment will emulate SolverControl by setting #reduce to zero.
+   */
+  IterationNumberControl &operator= (const SolverControl &c);
+
+  /**
+   * Virtual destructor is needed as there are virtual functions in this
+   * class.
+   */
+  virtual ~IterationNumberControl();
+
+  /**
+   * Decide about success or failure of an iteration. This function bases
+   * success solely on the fact if a given number of iterations was reached or
+   * the check value reached exactly zero.
+   */
+  virtual State check (const unsigned int step,
+                       const double   check_value);
+};
+
 /*@}*/
 //---------------------------------------------------------------------------
 
index 8f40d154a40dd4e21b5cc957e520921b17069d18..6eab7e4d4f608b61f24a841928c5ef8a3fa5fc9b 100644 (file)
@@ -303,4 +303,41 @@ ReductionControl::parse_parameters (ParameterHandler &param)
   set_reduction (param.get_double("Reduction"));
 }
 
+
+/*---------------------- IterationNumberControl -----------------------------*/
+
+
+IterationNumberControl::IterationNumberControl(const unsigned int n,
+                                               const double       tolerance,
+                                               const bool m_log_history,
+                                               const bool m_log_result)
+  :
+  SolverControl (n, tolerance, m_log_history, m_log_result){}
+
+
+IterationNumberControl::~IterationNumberControl()
+{}
+
+
+SolverControl::State
+IterationNumberControl::check (const unsigned int step,
+                               const double check_value)
+{
+  // check whether the given number of iterations was reached, and return
+  // success in that case. Otherwise, go on to the check of the base class.
+  if (step >= this->maxsteps)
+    {
+      if (m_log_result)
+        deallog << "Convergence step " << step
+                << " value " << check_value << std::endl;
+      lstep  = step;
+      lvalue = check_value;
+
+      lcheck = success;
+      return success;
+    }
+  else
+    return SolverControl::check(step, check_value);
+}
+
 DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/lac/iteration_number_control.cc b/tests/lac/iteration_number_control.cc
new file mode 100644 (file)
index 0000000..26ff41d
--- /dev/null
@@ -0,0 +1,85 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 1998 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// check that IterationNumber control works as expected, i.e., it terminates
+// with success when reaching the maximum number of iterations
+
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iomanip>
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/solver_control.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/precondition.h>
+
+
+void check()
+{
+  // Check 1: should still terminate after one iteration when using a matrix
+  // with only one eigenvalue
+  const unsigned int size = 20;
+  FullMatrix<double> A(size, size);
+  for (unsigned int i=0; i<size; ++i)
+    A(i,i) = 2;
+
+  Vector<double> in(size), out(size);
+  in = 1.;
+  
+  {
+    IterationNumberControl control(5);
+    SolverCG<> solver(control);
+    solver.solve(A, out, in, PreconditionIdentity());
+    Assert (control.last_step() == 1, ExcInternalError());
+  }
+  for (unsigned int i=0; i<size; ++i)
+    Assert(std::abs(out(i)-0.5) < 1e-12, ExcInternalError());
+
+  // Check 2: should only do 5 iterations but the solution should not be exact
+  for (unsigned int i=0; i<size; ++i)
+    A(i,i) = 1+i;
+
+  out = 0;
+  {
+    IterationNumberControl control(5);
+    SolverCG<> solver(control);
+    solver.solve(A, out, in, PreconditionIdentity());
+    Assert (control.last_step() == 5, ExcInternalError());
+  }
+  bool solved_exactly = true;
+  for (unsigned int i=0; i<size; ++i)
+    if (std::abs(out(i)-1./(1+size)) > 1e-8)
+      solved_exactly = false;
+  Assert(solved_exactly == false, ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
+
+int main()
+{
+  std::ofstream logfile("output");
+//  logfile.setf(std::ios::fixed);
+  deallog << std::setprecision(4);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check();
+}
+
diff --git a/tests/lac/iteration_number_control.output b/tests/lac/iteration_number_control.output
new file mode 100644 (file)
index 0000000..e2c05a8
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL:cg::Starting value 4.472
+DEAL:cg::Convergence step 1 value 0
+DEAL:cg::Starting value 4.472
+DEAL:cg::Convergence step 5 value 0.6308
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.