]> https://gitweb.dealii.org/ - dealii.git/commitdiff
introduce get_history_data function 4253/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Thu, 13 Apr 2017 02:25:20 +0000 (20:25 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Thu, 13 Apr 2017 18:59:55 +0000 (12:59 -0600)
include/deal.II/lac/solver_control.h
source/lac/solver_control.cc
tests/lac/solver_control_04.cc [new file with mode: 0644]
tests/lac/solver_control_04.output [new file with mode: 0644]

index d3ca328c4283b92d709de237721ef36ce8354c2e..ec6d0435b2391d0c1446bc8bd69037a8d14c9c8d 100644 (file)
@@ -250,6 +250,11 @@ public:
    */
   void enable_history_data();
 
+  /**
+   * Provide read access to the collected residual data.
+   */
+  const std::vector<double> &get_history_data() const;
+
   /**
    * Average error reduction over all steps.
    *
index 9fa6490dfe234fe242878cc1f8d4da57dbb582f9..82d3a8452a6250fde28087d6fdf54383cfae95ae 100644 (file)
@@ -63,8 +63,6 @@ SolverControl::check (const unsigned int step,
   if (step==0)
     {
       initial_val = check_value;
-      if (history_data_enabled)
-        history_data.resize(maxsteps);
     }
 
   if (m_log_history && ((step % m_log_frequency) == 0))
@@ -83,7 +81,7 @@ SolverControl::check (const unsigned int step,
     }
 
   if (history_data_enabled)
-    history_data[step] = check_value;
+    history_data.push_back(check_value);
 
   if (check_value <= tol)
     {
@@ -158,6 +156,20 @@ SolverControl::enable_history_data ()
 }
 
 
+
+const std::vector<double> &SolverControl::get_history_data() const
+{
+  Assert (history_data_enabled, ExcHistoryDataRequired());
+  Assert (history_data.size() > 0,
+          ExcMessage("The SolverControl object was asked for the solver history "
+                     "data, but there is no data. Possibly you requested the data before the "
+                     "solver was run."));
+
+  return history_data;
+}
+
+
+
 double
 SolverControl::average_reduction() const
 {
diff --git a/tests/lac/solver_control_04.cc b/tests/lac/solver_control_04.cc
new file mode 100644 (file)
index 0000000..e710158
--- /dev/null
@@ -0,0 +1,114 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test SolverControl's get_history_data function
+// This test is adapted from tests/trilinos/solver_03.cc
+
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include "../testmatrix.h"
+#include <cmath>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/vector_memory.h>
+#include <deal.II/lac/solver_control.h>
+#include <deal.II/lac/solver_richardson.h>
+#include <deal.II/lac/solver_relaxation.h>
+#include <deal.II/lac/precondition.h>
+
+template<typename MatrixType, typename VectorType, class PRECONDITION>
+void
+check_solve (SolverControl       &solver_control,
+             const MatrixType    &A,
+             VectorType          &u,
+             VectorType          &f,
+             const PRECONDITION  &P,
+             const bool           expected_result)
+{
+  SolverCG<VectorType> solver(solver_control);
+
+  u = 0.;
+  f = 1.;
+  bool success = false;
+  try
+    {
+      solver.solve(A,u,f,P);
+      deallog << "Success. " << std::endl;
+      success = true;
+    }
+  catch (std::exception &e)
+    {
+      deallog << "Failure. " << std::endl;
+    }
+
+  deallog << "Solver history data: " << solver_control.get_history_data()
+          << std::endl;
+  Assert(success == expected_result, ExcMessage("Incorrect result."));
+}
+
+
+int main(int argc, char **argv)
+{
+  std::ofstream logfile("output");
+  logfile.precision(4);
+  deallog.attach(logfile);
+  deallog.threshold_double(1.e-10);
+
+  {
+    const unsigned int size = 32;
+    unsigned int dim = (size-1)*(size-1);
+
+    deallog << "Size " << size << " Unknowns " << dim << std::endl;
+
+    // Make matrix
+    FDMatrix testproblem(size, size);
+    SparsityPattern structure(dim, dim, 5);
+    testproblem.five_point_structure(structure);
+    structure.compress();
+    SparseMatrix<double>  A(structure);
+    testproblem.five_point(A);
+
+    Vector<double>  f(dim);
+    Vector<double>  u(dim);
+    f = 1.;
+
+    PreconditionJacobi<> preconditioner;
+    preconditioner.initialize(A);
+
+    deallog.push("Abs tol");
+    {
+      // Expects success
+      SolverControl solver_control(2000, 1.e-3);
+      solver_control.enable_history_data();
+
+      check_solve (solver_control, A,u,f, preconditioner, true);
+    }
+    deallog.pop();
+    deallog.push("Iterations");
+    {
+      // Expects failure
+      SolverControl solver_control(20, 1.e-3);
+      solver_control.enable_history_data();
+
+      check_solve (solver_control, A,u,f, preconditioner, false);
+    }
+    deallog.pop();
+  }
+}
diff --git a/tests/lac/solver_control_04.output b/tests/lac/solver_control_04.output
new file mode 100644 (file)
index 0000000..e7a3eb4
--- /dev/null
@@ -0,0 +1,10 @@
+
+DEAL::Size 32 Unknowns 961
+DEAL:Abs tol:cg::Starting value 31.0000
+DEAL:Abs tol:cg::Convergence step 43 value 0.000818932
+DEAL:Abs tol::Success. 
+DEAL:Abs tol::Solver history data: 31.0000 83.4701 77.7927 79.1157 77.5543 68.9713 66.7436 59.1467 56.1131 49.6461 46.1902 40.5905 36.9687 32.0317 28.3918 23.9735 20.3921 16.3806 12.8817 9.19955 6.00198 3.44349 2.20296 2.04617 2.24528 1.46106 0.914875 0.832960 0.599735 0.417881 0.381758 0.254475 0.196013 0.129650 0.0908879 0.0642667 0.0399487 0.0267517 0.0157501 0.00919834 0.00527074 0.00269181 0.00156656 0.000818932
+DEAL:Iterations:cg::Starting value 31.0000
+DEAL:Iterations:cg::Failure step 20 value 6.00198
+DEAL:Iterations::Failure. 
+DEAL:Iterations::Solver history data: 31.0000 83.4701 77.7927 79.1157 77.5543 68.9713 66.7436 59.1467 56.1131 49.6461 46.1902 40.5905 36.9687 32.0317 28.3918 23.9735 20.3921 16.3806 12.8817 9.19955 6.00198

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.