virtual void set_solver_type (KSP &ksp) const;
};
+
+/**
+ * An implementation of the solver interface using the PETSc PREONLY
+ * solver. Actually this is NOT a real solution algorithm. Its only
+ * purpose is to provide a solver object, when the preconditioner
+ * should be used as real solver. It is very useful in conjunction with
+ * the complete LU decomposition preconditioner <tt> PreconditionLU </tt>,
+ * which in conjunction with this solver class becomes a direct solver.
+ *
+ * @ingroup PETScWrappers
+ * @author Wolfgang Bangerth, 2004, Oliver Kayser-Herold, 2004
+ */
+ class SolverPreOnly : public SolverBase
+ {
+ public:
+ /**
+ * Standardized data struct to
+ * pipe additional data to the
+ * solver.
+ */
+ struct AdditionalData
+ {};
+
+ /**
+ * Constructor. In contrast to
+ * deal.II's own solvers, there is no
+ * need to give a vector memory
+ * object. However, PETSc solvers want
+ * to have an MPI communicator context
+ * over which computations are
+ * parallelized. By default,
+ * @p PETSC_COMM_SELF is used here,
+ * but you can change this. Note that
+ * for single processor (non-MPI)
+ * versions, this parameter does not
+ * have any effect.
+ *
+ * The last argument takes a structure
+ * with additional, solver dependent
+ * flags for tuning.
+ *
+ * Note that the communicator used here
+ * must match the communicator used in
+ * the system matrix, solution, and
+ * right hand side object of the solve
+ * to be done with this
+ * solver. Otherwise, PETSc will
+ * generate hard to track down errors,
+ * see the documentation of the
+ * SolverBase class.
+ */
+ SolverPreOnly (SolverControl &cn,
+ const MPI_Comm &mpi_communicator = PETSC_COMM_SELF,
+ const AdditionalData &data = AdditionalData());
+
+ protected:
+ /**
+ * Store a copy of the flags for this
+ * particular solver.
+ */
+ const AdditionalData additional_data;
+
+ /**
+ * Function that takes a Krylov
+ * Subspace Solver context object, and
+ * sets the type of solver that is
+ * appropriate for this class.
+ */
+ virtual void set_solver_type (KSP &ksp) const;
+ };
+
}
#endif // DEAL_II_USE_PETSC
preconditioner.set_preconditioner_type (pc);
- // in the deal.II solvers, we always
- // honor the initial guess in the
- // solution vector. do so here as well:
- KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
-
// then a convergence monitor
// function. that function simply checks
// with the solver_control object we have
preconditioner.set_preconditioner_type (solver_data->pc);
- // in the deal.II solvers, we always
- // honor the initial guess in the
- // solution vector. do so here as
- // well:
- KSPSetInitialGuessNonzero (solver_data->ksp, PETSC_TRUE);
-
// then a convergence monitor
// function. that function simply
// checks with the solver_control
// set the damping factor from the data
ierr = KSPRichardsonSetScale (ksp, additional_data.omega);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPCG));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPBICG));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
ierr = KSPSetPreconditionerSide(ksp, PC_RIGHT);
AssertThrow (ierr == 0, ExcPETScError(ierr));
}
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPBCGS));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPCGS));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPTFQMR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPTCQMR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPCR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
}
int ierr;
ierr = KSPSetType (ksp, const_cast<char *>(KSPLSQR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // in the deal.II solvers, we always
+ // honor the initial guess in the
+ // solution vector. do so here as well:
+ KSPSetInitialGuessNonzero (ksp, PETSC_TRUE);
+ }
+
+
+/* ---------------------- SolverPreOnly ------------------------ */
+
+ SolverPreOnly::SolverPreOnly (SolverControl &cn,
+ const MPI_Comm &mpi_communicator,
+ const AdditionalData &data)
+ :
+ SolverBase (cn, mpi_communicator),
+ additional_data (data)
+ {}
+
+
+ void
+ SolverPreOnly::set_solver_type (KSP &ksp) const
+ {
+ // set the type of solver. work around a
+ // problem in PETSc 2.1.6, where it asks
+ // for a char*, even though KSPXXXX is of
+ // type const char*
+ int ierr;
+ ierr = KSPSetType (ksp, const_cast<char *>(KSPPREONLY));
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ // The KSPPREONLY solver of
+ // PETSc never calls the convergence
+ // monitor, which leads to failure
+ // even when everything was ok.
+ // Therefore the SolverControl status
+ // is set to some nice values, which
+ // guarantee a nice result at the end
+ // of the solution process.
+ solver_control.check (1, 0.0);
+
+ // Using the PREONLY solver with
+ // a nonzero initial guess leads
+ // PETSc to produce some error messages.
+ KSPSetInitialGuessNonzero (ksp, PETSC_FALSE);
}
-
}
#else