]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove, without deprecation, the Eisenstat preconditioner. 12446/head
authorDavid Wells <drwells@email.unc.edu>
Fri, 11 Jun 2021 16:59:28 +0000 (12:59 -0400)
committerDavid Wells <drwells@email.unc.edu>
Fri, 11 Jun 2021 16:59:28 +0000 (12:59 -0400)
This has never worked in parallel and has suffered from a variety of other
issues over the years.

doc/news/changes/incompatibilities/20210611DavidWells [new file with mode: 0644]
include/deal.II/lac/petsc_precondition.h
source/lac/petsc_precondition.cc
tests/petsc/reinit_preconditioner_01.cc
tests/petsc/reinit_preconditioner_02.cc
tests/petsc/solver_03_precondition_eisenstat.cc [deleted file]
tests/petsc/solver_03_precondition_eisenstat.output [deleted file]

diff --git a/doc/news/changes/incompatibilities/20210611DavidWells b/doc/news/changes/incompatibilities/20210611DavidWells
new file mode 100644 (file)
index 0000000..5078277
--- /dev/null
@@ -0,0 +1,4 @@
+Removed: The PETSc Eisenstat preconditioner wrapper, which has never worked
+correctly in parallel, has been removed.
+<br>
+(David Wells, 2021/06/11)
index 9e5585ecd1c89e582f0975f541599d26dbd6f378..b17c5e2a9f1d4421ff45790f9d6239c6d10c9b1f 100644 (file)
@@ -397,72 +397,6 @@ namespace PETScWrappers
     AdditionalData additional_data;
   };
 
-
-
-  /**
-   * A class that implements the interface to use the PETSc Eisenstat
-   * preconditioner, which implements SSOR on the diagonal block owned by
-   * each processor.
-   *
-   * See the comment in the base class
-   * @ref PreconditionBase
-   * for when this preconditioner may or may not work.
-   *
-   * @ingroup PETScWrappers
-   */
-  class PreconditionEisenstat : public PreconditionBase
-  {
-  public:
-    /**
-     * Standardized data struct to pipe additional flags to the
-     * preconditioner.
-     */
-    struct AdditionalData
-    {
-      /**
-       * Constructor. By default, set the damping parameter to one.
-       */
-      AdditionalData(const double omega = 1);
-
-      /**
-       * Relaxation parameter.
-       */
-      double omega;
-    };
-
-    /**
-     * Empty Constructor. You need to call initialize() before using this
-     * object.
-     */
-    PreconditionEisenstat() = default;
-
-    /**
-     * Constructor. Take the matrix which is used to form the preconditioner,
-     * and additional flags if there are any.
-     */
-    PreconditionEisenstat(
-      const MatrixBase &    matrix,
-      const AdditionalData &additional_data = AdditionalData());
-
-    /**
-     * Initialize the preconditioner object and calculate all data that is
-     * necessary for applying it in a solver. This function is automatically
-     * called when calling the constructor with the same arguments and is only
-     * used if you create the preconditioner without arguments.
-     */
-    void
-    initialize(const MatrixBase &    matrix,
-               const AdditionalData &additional_data = AdditionalData());
-
-  protected:
-    /**
-     * Store a copy of the flags for this particular preconditioner.
-     */
-    AdditionalData additional_data;
-  };
-
-
-
   /**
    * A class that implements the interface to use the PETSc Incomplete
    * Cholesky preconditioner.
index 267d2c63c3bec1e4a422832777f88f6ad710c807..cdcc9caac63ef507717c233ec04365217b708ab1 100644 (file)
@@ -307,48 +307,6 @@ namespace PETScWrappers
   }
 
 
-  /* ----------------- PreconditionEisenstat -------------------- */
-
-  PreconditionEisenstat::AdditionalData::AdditionalData(const double omega)
-    : omega(omega)
-  {}
-
-
-
-  PreconditionEisenstat::PreconditionEisenstat(
-    const MatrixBase &    matrix,
-    const AdditionalData &additional_data)
-  {
-    initialize(matrix, additional_data);
-  }
-
-
-  void
-  PreconditionEisenstat::initialize(const MatrixBase &    matrix_,
-                                    const AdditionalData &additional_data_)
-  {
-    clear();
-
-    matrix          = static_cast<Mat>(matrix_);
-    additional_data = additional_data_;
-
-    create_pc();
-
-    PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCEISENSTAT));
-    AssertThrow(ierr == 0, ExcPETScError(ierr));
-
-    // then set flags as given
-    ierr = PCEisenstatSetOmega(pc, additional_data.omega);
-    AssertThrow(ierr == 0, ExcPETScError(ierr));
-
-    ierr = PCSetFromOptions(pc);
-    AssertThrow(ierr == 0, ExcPETScError(ierr));
-
-    ierr = PCSetUp(pc);
-    AssertThrow(ierr == 0, ExcPETScError(ierr));
-  }
-
-
   /* ----------------- PreconditionICC -------------------- */
 
 
index ae4252f5617d333f3a8ad3eb072817de281f8731..11ba3928af276d1746ae95371da1a934381bb7c3 100644 (file)
@@ -98,7 +98,6 @@ main(int argc, char **argv)
 
   test<PETScWrappers::PreconditionJacobi>();
   test<PETScWrappers::PreconditionBlockJacobi>();
-  test<PETScWrappers::PreconditionEisenstat>();
   test<PETScWrappers::PreconditionBoomerAMG>();
   test<PETScWrappers::PreconditionParaSails>();
   test<PETScWrappers::PreconditionNone>();
index 6be18763173b3c0dc2dad2b84451a606122e2485..9c063b4ba99b5163e0193cf31c6e49852807700a 100644 (file)
@@ -89,7 +89,6 @@ main(int argc, char **argv)
   test<PETScWrappers::PreconditionBlockJacobi>();
   test<PETScWrappers::PreconditionSOR>();
   test<PETScWrappers::PreconditionSSOR>();
-  // todo: this crashes test<PETScWrappers::PreconditionEisenstat> ();
   test<PETScWrappers::PreconditionICC>();
   test<PETScWrappers::PreconditionILU>();
   test<PETScWrappers::PreconditionLU>();
diff --git a/tests/petsc/solver_03_precondition_eisenstat.cc b/tests/petsc/solver_03_precondition_eisenstat.cc
deleted file mode 100644 (file)
index 2865a21..0000000
+++ /dev/null
@@ -1,70 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 2004 - 2020 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 PETSc CG solver
-
-
-#include <deal.II/lac/petsc_precondition.h>
-#include <deal.II/lac/petsc_solver.h>
-#include <deal.II/lac/petsc_sparse_matrix.h>
-#include <deal.II/lac/petsc_vector.h>
-#include <deal.II/lac/vector_memory.h>
-
-#include <iostream>
-#include <typeinfo>
-
-#include "../tests.h"
-
-#include "../testmatrix.h"
-
-
-int
-main(int argc, char **argv)
-{
-  initlog();
-  deallog << std::setprecision(4);
-
-  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
-  {
-    SolverControl control(100, 1.e-3);
-
-    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);
-    PETScWrappers::SparseMatrix A(dim, dim, 5);
-    testproblem.five_point(A);
-
-    IndexSet indices(dim);
-    indices.add_range(0, dim);
-    PETScWrappers::MPI::Vector f(indices, MPI_COMM_WORLD);
-    PETScWrappers::MPI::Vector u(indices, MPI_COMM_WORLD);
-    f = 1.;
-    A.compress(VectorOperation::insert);
-
-    PETScWrappers::SolverCG              solver(control);
-    PETScWrappers::PreconditionEisenstat preconditioner(A);
-    deallog << "Solver type: " << typeid(solver).name() << std::endl;
-
-    check_solver_within_range(solver.solve(A, u, f, preconditioner),
-                              control.last_step(),
-                              19,
-                              21);
-  }
-}
diff --git a/tests/petsc/solver_03_precondition_eisenstat.output b/tests/petsc/solver_03_precondition_eisenstat.output
deleted file mode 100644 (file)
index 1587b4f..0000000
+++ /dev/null
@@ -1,4 +0,0 @@
-
-DEAL::Size 32 Unknowns 961
-DEAL::Solver type: N6dealii13PETScWrappers8SolverCGE
-DEAL::Solver stopped within 19 - 21 iterations

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.