* Base class for solver classes using the SLEPc solvers which are
* selected based on flags passed to the eigenvalue problem solver
* context. Derived classes set the right flags to set the right
- * solver. On the other hand, note that: the
- * <code>AdditionalData</code> structure is a dummy structure that
- * currently exists for backward/forward compatibility only.
+ * solver. Note that: the <code>AdditionalData</code> structure is a
+ * dummy structure that currently exists for backward/forward
+ * compatibility only.
*
* The SLEPc solvers are intended to be used for solving the
- * generalized eigenspectrum problem $(A-\lambda M)x=0$, for $x\neq0$;
- * where $A$ is a system matrix, $M$ is a mass matrix, and $\lambda,
+ * generalized eigenspectrum problem $(A-\lambda B)x=0$, for $x\neq0$;
+ * where $A$ is a system matrix, $B$ is a mass matrix, and $\lambda,
* x$ are a set of eigenvalues and eigenvectors respectively. The
* emphasis is on methods and techniques appropriate for problems in
* which the associated matrices are sparse. Most of the methods
* following way:
@verbatim
SolverControl solver_control (1000, 1e-9);
- SolverArnoldi system (solver_control,
- mpi_communicator);
- system.solve (A, M, lambda, x, size_of_spectrum);
+ SolverArnoldi system (solver_control, mpi_communicator);
+ system.solve (A, B, lambda, x, size_of_spectrum);
@endverbatim
- * for the generalized eigenvalue problem $Ax=M\lambda x$, where the
+ * for the generalized eigenvalue problem $Ax=B\lambda x$, where the
* variable <code>const unsigned int size_of_spectrum</code> tells
* SLEPc the number of eigenvector/eigenvalue pairs to solve for: See
* also <code>step-36</code> for a hands-on example.
* @author Toby D. Young 2008, 2009, 2010; and Rickard Armiento 2008.
*
* @note Various tweaks and enhancments contributed by Eloy Romero and
- * Jose E. Roman 2009, 2010.
+ * Jose E. Roman 2009, 2010.
*/
namespace SLEPcWrappers
*/
virtual void set_solver_type (EPS &eps) const = 0;
- /**
+ /*
* Attributes that store the
* relevant information for the
- * eigenproblem solver context.
+ * eigenproblem solver
+ * context. These are: which
+ * portion of the spectrum to solve
+ * from; and the matrices and
+ * vectors that discribe the
+ * problem.
*/
EPSWhich set_which;
const PETScWrappers::MatrixBase *opB;
const PETScWrappers::VectorBase *initial_vector;
+ /**
+ * Pointer to an an object that
+ * describes transformations that
+ * can be applied to the eigenvalue
+ * problem.
+ */
SLEPcWrappers::TransformationBase *transformation;
private:
};
+/**
+ * An implementation of the solver interface using the SLEPc Lanczos
+ * solver. Usage: Largest values of spectrum only, all problem types,
+ * complex.
+ *
+ * @ingroup SLEPcWrappers
+ * @author Toby D. Young 2010
+ */
+ class SolverPower : public SolverBase
+ {
+ public:
+ /**
+ * Standardized data struct to pipe
+ * additional data to the solver,
+ * should it be needed.
+ */
+ struct AdditionalData
+ {};
+
+ /**
+ * SLEPc solvers will want to have
+ * an MPI communicator context over
+ * which computations are
+ * parallelized. By default, this
+ * carries the same behaviour has
+ * the PETScWrappers, but you can
+ * change that.
+ */
+ SolverPower (SolverControl &cn,
+ const MPI_Comm &mpi_communicator = PETSC_COMM_WORLD,
+ const AdditionalData &data = AdditionalData());
+
+ protected:
+
+ /**
+ * Store a copy of the flags for
+ * this particular solver.
+ */
+ const AdditionalData additional_data;
+
+ /**
+ * Function that takes a Eigenvalue
+ * Problem Solver context object,
+ * and sets the type of solver that
+ * is appropriate for this class.
+ */
+ virtual void set_solver_type (EPS &eps) const;
+
+ };
+
// --------------------------- inline and template functions -----------
/**
void
SolverBase::set_matrices (const PETScWrappers::MatrixBase &A)
{
+ // standard eigenspectrum problem
opA = &A;
opB = NULL;
}
SolverBase::set_matrices (const PETScWrappers::MatrixBase &A,
const PETScWrappers::MatrixBase &B)
{
+ // generalized eigenspectrum problem
opA = &A;
opB = &B;
}
solver_data.reset (new SolverData());
// create eigensolver context and
- // set operators.
+ // set operators
ierr = EPSCreate (mpi_communicator, &solver_data->eps);
AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+ // set eigenspectrum problem type
+ // (general/standard)
AssertThrow (opA, ExcSLEPcWrappersUsageError());
if (opB)
ierr = EPSSetOperators (solver_data->eps, *opA, *opB);
ierr = EPSSetOperators (solver_data->eps, *opA, PETSC_NULL);
AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+ // set the initial vector(s) if any
if (initial_vector && initial_vector->size() != 0)
{
AssertThrow (ierr == 0, ExcSLEPcError(ierr));
}
+ // set transformation type if any
if (transformation)
transformation->set_context(solver_data->eps);
- // set runtime options.
+ // set runtime options
set_solver_type (solver_data->eps);
+ // set which portion of the
+ // eigenspectrum to solve for
ierr = EPSSetWhichEigenpairs (solver_data->eps, set_which);
AssertThrow (ierr == 0, ExcSLEPcError(ierr));
AssertThrow (ierr == 0, ExcSLEPcError(ierr));
}
+ /* ----------------------- Power ------------------------- */
+
+ SolverPower::SolverPower (SolverControl &cn,
+ const MPI_Comm &mpi_communicator,
+ const AdditionalData &data)
+ :
+ SolverBase (cn, mpi_communicator),
+ additional_data (data)
+ {}
+
+ void
+ SolverPower::set_solver_type (EPS &eps) const
+ {
+ int ierr;
+ ierr = EPSSetType (eps, const_cast<char *>(EPSPOWER));
+ AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+
+ // hand over the absolute
+ // tolerance and the maximum
+ // number of iteration steps to
+ // the SLEPc convergence
+ // criterion.
+ ierr = EPSSetTolerances(eps, this->solver_control.tolerance(),
+ this->solver_control.max_steps());
+ AssertThrow (ierr == 0, ExcSLEPcError(ierr));
+ }
+
}
DEAL_II_NAMESPACE_CLOSE