]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove some old PETSc workarounds. 14161/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 27 Jul 2022 15:30:19 +0000 (11:30 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 27 Jul 2022 16:01:24 +0000 (12:01 -0400)
include/deal.II/lac/slepc_spectral_transformation.h
source/lac/petsc_parallel_vector.cc
source/lac/petsc_solver.cc
source/lac/slepc_spectral_transformation.cc

index 8a98dd0e457970c2e0d19eac27a18df63b6f5d06..4031ff83c0cdc04fa74023678a614eede8047d08 100644 (file)
@@ -203,9 +203,13 @@ namespace SLEPcWrappers
    * Spectrum Folding. This transformation type has been removed in SLEPc
    * 3.5.0 and thus cannot be used in the newer versions.
    *
+   * @deprecated Since deal.II requires PETSc 3.7 or newer this class no longer
+   * does anything.
+   *
    * @ingroup SLEPcWrappers
    */
-  class TransformationSpectrumFolding : public TransformationBase
+  class DEAL_II_DEPRECATED TransformationSpectrumFolding
+    : public TransformationBase
   {
   public:
     /**
index 97496965b6a2eb43083b34b123ee02b2b02de13a..4e7e29a553191069ef014a5fbd5dc249c5dfd2f6 100644 (file)
@@ -337,18 +337,6 @@ namespace PETScWrappers
                           static_cast<PetscInt>(ghost_indices.n_elements()));
       }
 #  endif
-
-
-      // in PETSc versions up to 3.5, VecCreateGhost zeroed out the locally
-      // owned vector elements but forgot about the ghost elements. we need to
-      // do this ourselves
-      //
-      // see https://code.google.com/p/dealii/issues/detail?id=233
-#  if DEAL_II_PETSC_VERSION_LT(3, 6, 0)
-      PETScWrappers::MPI::Vector zero;
-      zero.reinit(communicator, this->size(), locally_owned_size);
-      *this = zero;
-#  endif
     }
 
 
index b1515e6c66e8a55bb465f11def0341dd93bbebdd..2bb04ffe55ae48b5336419143db5e7e8840df95a 100644 (file)
@@ -94,16 +94,7 @@ namespace PETScWrappers
 
         // setting the preconditioner overwrites the used matrices.
         // hence, we need to set the matrices after the preconditioner.
-#  if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
-        // the last argument is irrelevant here,
-        // since we use the solver only once anyway
-        ierr = KSPSetOperators(solver_data->ksp,
-                               A,
-                               preconditioner,
-                               SAME_PRECONDITIONER);
-#  else
         ierr = KSPSetOperators(solver_data->ksp, A, preconditioner);
-#  endif
         AssertThrow(ierr == 0, ExcPETScError(ierr));
 
         // then a convergence monitor
@@ -398,37 +389,10 @@ namespace PETScWrappers
     PetscErrorCode ierr = KSPSetType(ksp, KSPGMRES);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
 
-    // set the restart parameter from the
-    // data. we would like to use the simple
-    // code that is commented out, but this
-    // leads to nasty warning and error
-    // messages due to some stupidity on
-    // PETSc's side: KSPGMRESSetRestart is
-    // implemented as a macro in which return
-    // statements are hidden. This may work
-    // if people strictly follow the PETSc
-    // coding style of always having
-    // functions return an integer error
-    // code, but the present function isn't
-    // like this.
-    /*
-        ierr = KSPGMRESSetRestart (ksp, additional_data.restart_parameter);
-        AssertThrow (ierr == 0, ExcPETScError(ierr));
-    */
-    // so rather expand their macros by hand,
-    // and do some equally nasty stuff that at
-    // least doesn't yield warnings...
-    int (*fun_ptr)(KSP, int);
-    ierr = PetscObjectQueryFunction(reinterpret_cast<PetscObject>(ksp),
-                                    "KSPGMRESSetRestart_C",
-                                    reinterpret_cast<void (**)()>(&fun_ptr));
+    ierr = KSPGMRESSetRestart(ksp, additional_data.restart_parameter);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
 
-    ierr = (*fun_ptr)(ksp, additional_data.restart_parameter);
-    AssertThrow(ierr == 0, ExcPETScError(ierr));
-
-    // Set preconditioning side to
-    // right
+    // Set preconditioning side to right
     if (additional_data.right_preconditioning)
       {
         ierr = KSPSetPCSide(ksp, PC_RIGHT);
@@ -711,12 +675,7 @@ namespace PETScWrappers
          * set the matrices involved. the last argument is irrelevant here,
          * since we use the solver only once anyway
          */
-#    if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
-        ierr =
-          KSPSetOperators(solver_data->ksp, A, A, DIFFERENT_NONZERO_PATTERN);
-#    else
         ierr = KSPSetOperators(solver_data->ksp, A, A);
-#    endif
         AssertThrow(ierr == 0, ExcPETScError(ierr));
 
         /*
index f1b856be34f5399287c7b5f03e632479fd97f294..dc68c3c8c71de09541f9dd145ad76a3e3007074d 100644 (file)
@@ -111,20 +111,12 @@ namespace SLEPcWrappers
     : TransformationBase(mpi_communicator)
     , additional_data(data)
   {
-#  if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
-    PetscErrorCode ierr = STSetType(st, const_cast<char *>(STFOLD));
-    AssertThrow(ierr == 0, SolverBase::ExcSLEPcError(ierr));
-
-    ierr = STSetShift(st, additional_data.shift_parameter);
-    AssertThrow(ierr == 0, SolverBase::ExcSLEPcError(ierr));
-#  else
-    // PETSc/SLEPc version must be < 3.5.0.
-    (void)st;
-    Assert((false),
+    // This feature is only in PETSc/SLEPc versions 3.4 or older, which we no
+    // longer support.
+    Assert(false,
            ExcMessage(
              "Folding transformation has been removed in SLEPc 3.5.0 and newer."
              " You cannot use this transformation anymore."));
-#  endif
   }
 
   /* ------------------- TransformationCayley --------------------- */

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.