]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove MGCoarseGridLACIteration 10524/head
authorDaniel Arndt <arndtd@ornl.gov>
Sun, 7 Jun 2020 03:15:52 +0000 (23:15 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Fri, 12 Jun 2020 17:54:06 +0000 (13:54 -0400)
14 files changed:
doc/news/changes/incompatibilities/20200612DanielArndt [new file with mode: 0644]
include/deal.II/multigrid/mg_coarse.h
tests/mpi/multigrid_adaptive.cc
tests/mpi/multigrid_uniform.cc
tests/mpi/step-39-block.cc
tests/mpi/step-39.cc
tests/multigrid/step-16-04.cc
tests/multigrid/step-16-05.cc
tests/multigrid/step-16-06.cc
tests/multigrid/step-16-50-mpi-linear-operator.cc
tests/multigrid/step-16-50-mpi-smoother.cc
tests/multigrid/step-16-50-mpi.cc
tests/multigrid/step-16-50-serial.cc
tests/multigrid/step-50_01.cc

diff --git a/doc/news/changes/incompatibilities/20200612DanielArndt b/doc/news/changes/incompatibilities/20200612DanielArndt
new file mode 100644 (file)
index 0000000..84eb883
--- /dev/null
@@ -0,0 +1,3 @@
+Removed: The deprecated class MGCoarseGridLACIteration has been removed.
+<br>
+(Daniel Arndt, 2020/06/12)
index 5b9c0ac3135c7a520f95c9c8b5d0bcdcf6f32ae8..8095f24bc89f7c192b44f9a46adab10639dbc8a4 100644 (file)
@@ -78,90 +78,6 @@ private:
 };
 
 
-/**
- * Coarse grid solver using LAC iterative methods. This is a little wrapper,
- * transforming a triplet of iterative solver, matrix and preconditioner into
- * a coarse grid solver.
- *
- * The type of the matrix (i.e. the template parameter @p MatrixType) should
- * be derived from @p Subscriptor to allow for the use of a smart pointer to
- * it.
- *
- * @deprecated Use MGCoarseGridIterativeSolver instead.
- */
-template <typename SolverType, class VectorType = Vector<double>>
-class DEAL_II_DEPRECATED MGCoarseGridLACIteration
-  : public MGCoarseGridBase<VectorType>
-{
-public:
-  /**
-   * Default constructor.
-   */
-  MGCoarseGridLACIteration();
-
-  /**
-   * Constructor. Store solver, matrix and preconditioning method for later
-   * use.
-   */
-  template <typename MatrixType, typename PreconditionerType>
-  MGCoarseGridLACIteration(SolverType &,
-                           const MatrixType &,
-                           const PreconditionerType &);
-
-  /**
-   * Destructor freeing the pointers.
-   */
-  ~MGCoarseGridLACIteration();
-
-  /**
-   * Initialize new data.
-   */
-  template <typename MatrixType, typename PreconditionerType>
-  void
-  initialize(SolverType &, const MatrixType &, const PreconditionerType &);
-
-  /**
-   * Clear all pointers.
-   */
-  void
-  clear();
-
-  /**
-   * Implementation of the abstract function. Calls the solver method with
-   * matrix, vectors and preconditioner.
-   */
-  void
-  operator()(const unsigned int level,
-             VectorType &       dst,
-             const VectorType & src) const;
-
-  /**
-   * Set the matrix. This gives the possibility to replace the matrix that
-   * was given to the constructor by a new matrix.
-   */
-  template <typename MatrixType>
-  void
-  set_matrix(const MatrixType &);
-
-private:
-  /**
-   * Reference to the solver.
-   */
-  SmartPointer<SolverType, MGCoarseGridLACIteration<SolverType, VectorType>>
-    solver;
-
-  /**
-   * LinearOperator wrapping a reference to the matrix.
-   */
-  LinearOperator<VectorType> matrix;
-
-  /**
-   * LinearOperator wrapping a reference to the preconditioner.
-   */
-  LinearOperator<VectorType> precondition;
-};
-
-
 
 /**
  * Coarse grid multigrid operator for an iterative solver.
@@ -368,93 +284,6 @@ MGCoarseGridApplySmoother<VectorType>::operator()(const unsigned int level,
   coarse_smooth->smooth(level, dst, src);
 }
 
-/* ------------------ Functions for MGCoarseGridLACIteration ------------ */
-
-
-template <typename SolverType, class VectorType>
-MGCoarseGridLACIteration<SolverType, VectorType>::MGCoarseGridLACIteration()
-  : solver(0, typeid(*this).name())
-  , matrix(0)
-  , precondition(0)
-{}
-
-
-template <typename SolverType, class VectorType>
-template <typename MatrixType, typename PreconditionerType>
-MGCoarseGridLACIteration<SolverType, VectorType>::MGCoarseGridLACIteration(
-  SolverType &              s,
-  const MatrixType &        m,
-  const PreconditionerType &p)
-  : solver(&s, typeid(*this).name())
-{
-  // Workaround: Unfortunately, not every "m" object has a rich enough
-  // interface to populate reinit_(domain|range)_vector. Thus, supply an
-  // empty LinearOperator exemplar.
-  matrix       = linear_operator<VectorType>(LinearOperator<VectorType>(), m);
-  precondition = linear_operator<VectorType>(matrix, p);
-}
-
-
-template <typename SolverType, class VectorType>
-MGCoarseGridLACIteration<SolverType, VectorType>::~MGCoarseGridLACIteration()
-{
-  clear();
-}
-
-
-template <typename SolverType, class VectorType>
-template <typename MatrixType, typename PreconditionerType>
-void
-MGCoarseGridLACIteration<SolverType, VectorType>::initialize(
-  SolverType &              s,
-  const MatrixType &        m,
-  const PreconditionerType &p)
-{
-  solver = &s;
-  // Workaround: Unfortunately, not every "m" object has a rich enough
-  // interface to populate reinit_(domain|range)_vector. Thus, supply an
-  // empty LinearOperator exemplar.
-  matrix       = linear_operator<VectorType>(LinearOperator<VectorType>(), m);
-  precondition = linear_operator<VectorType>(matrix, p);
-}
-
-
-template <typename SolverType, class VectorType>
-void
-MGCoarseGridLACIteration<SolverType, VectorType>::clear()
-{
-  solver       = nullptr;
-  matrix       = LinearOperator<VectorType>();
-  precondition = LinearOperator<VectorType>();
-}
-
-
-template <typename SolverType, class VectorType>
-void
-MGCoarseGridLACIteration<SolverType, VectorType>::
-operator()(const unsigned int /* level */,
-           VectorType &      dst,
-           const VectorType &src) const
-{
-  Assert(solver != nullptr, ExcNotInitialized());
-  solver->solve(matrix, dst, src, precondition);
-}
-
-
-template <typename SolverType, class VectorType>
-template <typename MatrixType>
-void
-MGCoarseGridLACIteration<SolverType, VectorType>::set_matrix(
-  const MatrixType &m)
-{
-  // Workaround: Unfortunately, not every "m" object has a rich enough
-  // interface to populate reinit_(domain|range)_vector. Thus, supply an
-  // empty LinearOperator exemplar.
-  matrix = linear_operator<VectorType>(LinearOperator<VectorType>(), m);
-}
-
-
-
 /* ------------------ Functions for MGCoarseGridIterativeSolver ------------ */
 
 template <class VectorType,
index 092e4df458ce120cc2d9bd318bd11baf551259f5..8f688b401327b82d6bd024fb3812fead0959dc06 100644 (file)
@@ -409,7 +409,10 @@ namespace Step50
     SolverControl         coarse_solver_control(1000, 1e-10, false, false);
     SolverGMRES<vector_t> coarse_solver(coarse_solver_control);
     PreconditionIdentity  id;
-    MGCoarseGridLACIteration<SolverGMRES<vector_t>, vector_t>
+    MGCoarseGridIterativeSolver<vector_t,
+                                SolverGMRES<vector_t>,
+                                matrix_t,
+                                PreconditionIdentity>
       coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
     typedef TrilinosWrappers::PreconditionJacobi         Smoother;
index 85930faed00259f1ed618299843b206baed51e13..ccd1e44190aa5b7305a40141f4b5c01296ffb052 100644 (file)
@@ -396,7 +396,10 @@ namespace Step50
     SolverControl         coarse_solver_control(1000, 1e-10, false, false);
     SolverGMRES<vector_t> coarse_solver(coarse_solver_control);
     PreconditionIdentity  id;
-    MGCoarseGridLACIteration<SolverGMRES<vector_t>, vector_t>
+    MGCoarseGridIterativeSolver<vector_t,
+                                SolverGMRES<vector_t>,
+                                matrix_t,
+                                PreconditionIdentity>
       coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
     typedef TrilinosWrappers::PreconditionJacobi         Smoother;
index c48f5cf093a8cbb8aafd7c51e94aaa9791575ada..5f1773aec619e0b727d9af6b49816225a604d7ce 100644 (file)
@@ -635,8 +635,10 @@ namespace Step39
       coarse_solver_control);
     PreconditionIdentity            identity;
     TrilinosWrappers::SparseMatrix &coarse_matrix = mg_matrix[0];
-    MGCoarseGridLACIteration<SolverCG<TrilinosWrappers::MPI::Vector>,
-                             TrilinosWrappers::MPI::Vector>
+    MGCoarseGridIterativeSolver<TrilinosWrappers::MPI::Vector,
+                                SolverCG<TrilinosWrappers::MPI::Vector>,
+                                TrilinosWrappers::SparseMatrix,
+                                PreconditionIdentity>
       coarse_grid_solver(coarse_solver, coarse_matrix, identity);
 
     typedef RelaxationBlockJacobi<TrilinosWrappers::SparseMatrix,
index 7c028dfeb19a5a93cbcfe4c133f20b2ea9c30544..f7a279e7987f9fad66eb3315f637fdc121739712 100644 (file)
@@ -636,8 +636,10 @@ namespace Step39
       coarse_solver_control);
     PreconditionIdentity            identity;
     TrilinosWrappers::SparseMatrix &coarse_matrix = mg_matrix[0];
-    MGCoarseGridLACIteration<SolverCG<TrilinosWrappers::MPI::Vector>,
-                             TrilinosWrappers::MPI::Vector>
+    MGCoarseGridIterativeSolver<TrilinosWrappers::MPI::Vector,
+                                SolverCG<TrilinosWrappers::MPI::Vector>,
+                                TrilinosWrappers::SparseMatrix,
+                                PreconditionIdentity>
       coarse_grid_solver(coarse_solver, coarse_matrix, identity);
 
     typedef TrilinosWrappers::PreconditionJacobi Smoother;
index 22df24bf8d9eeb3eed03130c5c0d1520f1ac8c04..9a2f1dad7b0227771a266c8aac78b2ff647fe81a 100644 (file)
@@ -380,8 +380,11 @@ LaplaceProblem<dim>::solve()
   SolverCG<LinearAlgebra::distributed::Vector<double>> coarse_solver(
     coarse_solver_control);
   PreconditionIdentity id;
-  MGCoarseGridLACIteration<SolverCG<LinearAlgebra::distributed::Vector<double>>,
-                           LinearAlgebra::distributed::Vector<double>>
+  MGCoarseGridIterativeSolver<
+    LinearAlgebra::distributed::Vector<double>,
+    SolverCG<LinearAlgebra::distributed::Vector<double>>,
+    SparseMatrix<double>,
+    PreconditionIdentity>
     coarse_grid_solver(coarse_solver, mg_matrices[min_level], id);
   deallog << "   Size of coarse grid matrix: " << mg_matrices[min_level].m()
           << std::endl;
index e3b81fa9a05618baba6d9b4b3f3fe60ffa565041..33d64c13c43a8768c2e1164655b7f5bfd87ae3c6 100644 (file)
@@ -381,8 +381,11 @@ LaplaceProblem<dim>::solve()
   SolverCG<LinearAlgebra::distributed::Vector<double>> coarse_solver(
     coarse_solver_control);
   PreconditionIdentity id;
-  MGCoarseGridLACIteration<SolverCG<LinearAlgebra::distributed::Vector<double>>,
-                           LinearAlgebra::distributed::Vector<double>>
+  MGCoarseGridIterativeSolver<
+    LinearAlgebra::distributed::Vector<double>,
+    SolverCG<LinearAlgebra::distributed::Vector<double>>,
+    SparseMatrix<double>,
+    PreconditionIdentity>
     coarse_grid_solver(coarse_solver, mg_matrices[min_level], id);
   deallog << "   Size of coarse grid matrix: " << mg_matrices[min_level].m()
           << std::endl;
index 3b485bcc4600a0ab6028d3aed740388ff9de6d82..a9813fdf517c7616dfe4e5bdcaa070c48f7fc974 100644 (file)
@@ -381,8 +381,11 @@ LaplaceProblem<dim>::solve()
   SolverCG<LinearAlgebra::distributed::Vector<double>> coarse_solver(
     coarse_solver_control);
   PreconditionIdentity id;
-  MGCoarseGridLACIteration<SolverCG<LinearAlgebra::distributed::Vector<double>>,
-                           LinearAlgebra::distributed::Vector<double>>
+  MGCoarseGridIterativeSolver<
+    LinearAlgebra::distributed::Vector<double>,
+    SolverCG<LinearAlgebra::distributed::Vector<double>>,
+    SparseMatrix<double>,
+    PreconditionIdentity>
     coarse_grid_solver(coarse_solver, mg_matrices[min_level], id);
   deallog << "   Size of coarse grid matrix: " << mg_matrices[min_level].m()
           << std::endl;
index 9c78b5efff41289c3feca314c85c65a847265b3c..8b06521f6080ea0eb0b47b817a5d427f06b58e1a 100644 (file)
@@ -442,8 +442,11 @@ namespace Step50
     SolverControl        coarse_solver_control(1000, 1e-10, false, false);
     SolverCG<vector_t>   coarse_solver(coarse_solver_control);
     PreconditionIdentity id;
-    MGCoarseGridLACIteration<SolverCG<vector_t>, vector_t> coarse_grid_solver(
-      coarse_solver, coarse_matrix, id);
+    MGCoarseGridIterativeSolver<vector_t,
+                                SolverCG<vector_t>,
+                                matrix_t,
+                                PreconditionIdentity>
+      coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
     typedef LA::MPI::PreconditionJacobi                  Smoother;
     MGSmootherPrecondition<matrix_t, Smoother, vector_t> mg_smoother;
index d13e0240e977e922665728680fd99b8bd54416d0..4e6fe6a3c3a85de2115e114d5715fbf27072642d 100644 (file)
@@ -16,7 +16,7 @@
 
 // Same as step-16-50, but use Jacobi smoother at the coarsest grid via
 // MGCoarseGridApplySmoother. In this particular case, the number of iterations
-// until convergence is exactly the same as for MGCoarseGridLACIteration.
+// until convergence is exactly the same as for MGCoarseGridIterativeSolver.
 
 #include <deal.II/base/conditional_ostream.h>
 #include <deal.II/base/function.h>
index 14c0d02d775409cbf6c835db02e365bb4783ebd4..d6ed95de0f01afd37086f2cac901c216d72704e2 100644 (file)
@@ -442,8 +442,11 @@ namespace Step50
     SolverControl        coarse_solver_control(1000, 1e-10, false, false);
     SolverCG<vector_t>   coarse_solver(coarse_solver_control);
     PreconditionIdentity id;
-    MGCoarseGridLACIteration<SolverCG<vector_t>, vector_t> coarse_grid_solver(
-      coarse_solver, coarse_matrix, id);
+    MGCoarseGridIterativeSolver<vector_t,
+                                SolverCG<vector_t>,
+                                matrix_t,
+                                PreconditionIdentity>
+      coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
     typedef LA::MPI::PreconditionJacobi                  Smoother;
     MGSmootherPrecondition<matrix_t, Smoother, vector_t> mg_smoother;
index 5fc76f023cdad91c0dc747502bf0326e3b725c19..6e7b0a4f3b3afaabceaec6120c403a5fbde352c2 100644 (file)
@@ -375,8 +375,11 @@ LaplaceProblem<dim>::solve()
   SolverControl        coarse_solver_control(1000, 1e-10, false, false);
   SolverCG<vector_t>   coarse_solver(coarse_solver_control);
   PreconditionIdentity id;
-  MGCoarseGridLACIteration<SolverCG<vector_t>, vector_t> coarse_grid_solver(
-    coarse_solver, coarse_matrix, id);
+  MGCoarseGridIterativeSolver<vector_t,
+                              SolverCG<vector_t>,
+                              matrix_t,
+                              PreconditionIdentity>
+    coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
   typedef PreconditionJacobi<matrix_t>                 Smoother;
   MGSmootherPrecondition<matrix_t, Smoother, vector_t> mg_smoother;
index df241c12deac9cde408be7fd69b72c1d7c098ad5..a7493242b5d729e111c2e8a83deec1094236f665 100644 (file)
@@ -443,8 +443,11 @@ namespace Step50
     SolverControl        coarse_solver_control(1000, 1e-10, false, false);
     SolverCG<vector_t>   coarse_solver(coarse_solver_control);
     PreconditionIdentity id;
-    MGCoarseGridLACIteration<SolverCG<vector_t>, vector_t> coarse_grid_solver(
-      coarse_solver, coarse_matrix, id);
+    MGCoarseGridIterativeSolver<vector_t,
+                                SolverCG<vector_t>,
+                                matrix_t,
+                                PreconditionIdentity>
+      coarse_grid_solver(coarse_solver, coarse_matrix, id);
 
     typedef LA::MPI::PreconditionJacobi                  Smoother;
     MGSmootherPrecondition<matrix_t, Smoother, vector_t> mg_smoother;

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.