]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement SolverSelector constructor mentioned in the documentation
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 16 Sep 2018 21:49:33 +0000 (23:49 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sun, 16 Sep 2018 22:11:44 +0000 (00:11 +0200)
doc/news/changes/minor/20180917DanielArndt [new file with mode: 0644]
include/deal.II/lac/solver_selector.h
tests/lac/solver_selector_03.cc [new file with mode: 0644]
tests/lac/solver_selector_03.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20180917DanielArndt b/doc/news/changes/minor/20180917DanielArndt
new file mode 100644 (file)
index 0000000..f61b1bf
--- /dev/null
@@ -0,0 +1,3 @@
+New: SolverSelector gained a new constructor initializing the object completely.
+<br>
+(Daniel Arndt, 2018/09/17)
index e5a1f57833d80096c8f9728506a0d3085fdff715..387ed3c920c8980ade7c6335ad8dd91b7f7da4ac 100644 (file)
@@ -30,7 +30,6 @@
 #include <deal.II/lac/solver_minres.h>
 #include <deal.II/lac/solver_richardson.h>
 #include <deal.II/lac/vector.h>
-#include <deal.II/lac/vector_memory.h>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -47,14 +46,9 @@ DEAL_II_NAMESPACE_OPEN
  *
  * <h3>Usage</h3> The simplest use of this class is the following:
  * @code
- * // generate a @p SolverControl and a @p VectorMemory
- * SolverControl control;
- * VectorMemory<Vector<double> > memory;
- *
- * // Line 3:
  * // generate a @p SolverSelector that calls the @p SolverCG
- * SolverSelector<Vector<double> >
- *   solver_selector("cg", control, memory);
+ * SolverConrtol control;
+ * SolverSelector<Vector<double> > solver_selector ("cg", control);
  *
  * // generate e.g. a @p PreconditionRelaxation
  * PreconditionRelaxation<SparseMatrix<double>, Vector<double> >
@@ -75,13 +69,13 @@ DEAL_II_NAMESPACE_OPEN
  * ...
  * @endcode
  * Assuming that in the users parameter file there exists the line
- * @verbatim
+ * @code
  * set solver = cg
- * @endverbatim
- * then `Line 3' of the above example reads
+ * @endcode
+ * then the constructor call in the above example can be written as
  * @code
  * SolverSelector<SparseMatrix<double>, Vector<double> >
- *   solver_selector(prm.get("solver"), control, memory);
+ *   solver_selector(prm.get("solver"), control);
  * @endcode
  *
  *
@@ -89,9 +83,7 @@ DEAL_II_NAMESPACE_OPEN
  * to change his program. Only in the implementation of the @p SolverSelector
  * the calling of this solver has to be added and each user with program lines
  * quoted above only needs to 'set solver = xyz' in his parameter file to get
- * access to that new solver.  :-)
- *
- * (By the way, thanks to Wolfgang for implementing the @p ParameterHandler.)
+ * access to that new solver.
  *
  * @author Ralf Hartmann, 1999
  */
@@ -109,6 +101,12 @@ public:
    */
   SolverSelector() = default;
 
+  /**
+   * Constructor, selecting the solver @p name
+   * and the SolverControl object @p control already.
+   */
+  SolverSelector(const std::string &name, SolverControl &control);
+
   /**
    * Destructor
    */
@@ -249,11 +247,21 @@ private:
 /* --------------------- Inline and template functions ------------------- */
 
 
+template <typename VectorType>
+SolverSelector<VectorType>::SolverSelector(const std::string &name,
+                                           SolverControl &    solver_control)
+  : solver_name(name)
+  , control(&solver_control)
+{}
+
+
+
 template <typename VectorType>
 SolverSelector<VectorType>::~SolverSelector()
 {}
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::select(const std::string &name)
@@ -262,6 +270,7 @@ SolverSelector<VectorType>::select(const std::string &name)
 }
 
 
+
 template <typename VectorType>
 template <class Matrix, class Preconditioner>
 void
@@ -305,6 +314,7 @@ SolverSelector<VectorType>::solve(const Matrix &        A,
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_control(SolverControl &ctrl)
@@ -313,6 +323,7 @@ SolverSelector<VectorType>::set_control(SolverControl &ctrl)
 }
 
 
+
 template <typename VectorType>
 std::string
 SolverSelector<VectorType>::get_solver_names()
@@ -321,6 +332,7 @@ SolverSelector<VectorType>::get_solver_names()
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -330,6 +342,7 @@ SolverSelector<VectorType>::set_data(
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -339,6 +352,7 @@ SolverSelector<VectorType>::set_data(
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -348,6 +362,7 @@ SolverSelector<VectorType>::set_data(
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -357,6 +372,7 @@ SolverSelector<VectorType>::set_data(
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -366,6 +382,7 @@ SolverSelector<VectorType>::set_data(
 }
 
 
+
 template <typename VectorType>
 void
 SolverSelector<VectorType>::set_data(
@@ -374,7 +391,6 @@ SolverSelector<VectorType>::set_data(
   bicgstab_data = data;
 }
 
-
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
diff --git a/tests/lac/solver_selector_03.cc b/tests/lac/solver_selector_03.cc
new file mode 100644 (file)
index 0000000..2952aff
--- /dev/null
@@ -0,0 +1,76 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Test the SolverSelector class using the constructor that initializes the
+// object already.
+
+#include <deal.II/lac/solver_selector.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/vector_memory.h>
+
+#include "../testmatrix.h"
+#include "../tests.h"
+
+
+template <typename MatrixType, typename VectorType>
+void
+check(const MatrixType &A, const VectorType &f)
+{
+  std::vector<std::string> names;
+  names.push_back("cg");
+  names.push_back("gmres");
+
+  SolverControl                          cont(100, 1.e-7);
+  SolverSelector<VectorType>             solver;
+  PreconditionSSOR<SparseMatrix<double>> pre;
+  pre.initialize(A);
+
+  VectorType u;
+  u.reinit(f);
+
+  for (const auto &name : names)
+    {
+      SolverSelector<VectorType> solver(name, cont);
+      u = 0.;
+      solver.solve(A, u, f, pre);
+    }
+}
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(4);
+
+  unsigned int size = 37;
+  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);
+  f = 1.;
+
+  check(A, f);
+}
diff --git a/tests/lac/solver_selector_03.output b/tests/lac/solver_selector_03.output
new file mode 100644 (file)
index 0000000..a716a5e
--- /dev/null
@@ -0,0 +1,6 @@
+
+DEAL::Size 37 Unknowns 1296
+DEAL:cg::Starting value 36.00
+DEAL:cg::Convergence step 39 value 7.335e-08
+DEAL:GMRES::Starting value 34.45
+DEAL:GMRES::Convergence step 39 value 6.627e-08

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.