// $Id$
// Version: $Name$
//
-// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009 by the deal.II authors
+// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2009, 2010 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
/**
* Change the relaxaton parameter.
*/
- void initialize (const AdditionalData parameters);
+ void initialize (const AdditionalData ¶meters);
/**
* Change the relaxaton parameter
*/
template <class MATRIX>
void initialize (const MATRIX&,
- const AdditionalData parameters);
+ const AdditionalData ¶meters);
/**
* Apply preconditioner.
*/
template<class VECTOR>
void Tvmult_add (VECTOR&, const VECTOR&) const;
-
+
private:
/**
* The relaxation parameter
* function of the matrix.
*/
typedef void ( MATRIX::* function_ptr)(VECTOR&, const VECTOR&) const;
-
+
/**
* Constructor.
* This constructor stores a
*/
PreconditionUseMatrix(const MATRIX &M,
const function_ptr method);
-
+
/**
* Execute preconditioning. Calls the
* function passed to the constructor
* Pointer to the matrix in use.
*/
const MATRIX &matrix;
-
+
/**
* Pointer to the preconditioning
* function.
/**
* Relaxation parameter.
*/
- double relaxation;
+ double relaxation;
};
-
+
/**
* Initialize matrix and
* relaxation parameter. The
* reasons. It defaults to 1.
*/
void initialize (const MATRIX &A,
- AdditionalData parameters = AdditionalData());
-
+ const AdditionalData & parameters = AdditionalData());
+
/**
* Release the matrix and reset
* its pointer.
*/
void clear();
-
+
protected:
/**
* Pointer to the matrix object.
*/
typedef PreconditionRelaxation<MATRIX> BaseClass;
-
+
/**
* Initialize matrix and
* relaxation parameter. The
* reasons. It defaults to 1.
*/
void initialize (const MATRIX &A,
- typename BaseClass::AdditionalData parameters = typename BaseClass::AdditionalData());
+ const typename BaseClass::AdditionalData ¶meters = typename BaseClass::AdditionalData());
/**
* Apply preconditioner.
void initialize (const MATRIX &A,
const std::vector<unsigned int> &permutation,
const std::vector<unsigned int> &inverse_permutation,
- typename PreconditionRelaxation<MATRIX>::AdditionalData
+ const typename PreconditionRelaxation<MATRIX>::AdditionalData &
parameters = typename PreconditionRelaxation<MATRIX>::AdditionalData());
-
+
/**
* Apply preconditioner.
*/
void initialize (SOLVER&,
const MATRIX&,
const PRECONDITION&);
-
+
/**
* Execute preconditioning.
*/
* The matrix in use.
*/
SmartPointer<const MATRIX,PreconditionLACSolver<SOLVER,MATRIX,PRECONDITION> > matrix;
-
+
/**
* The preconditioner to use.
*/
template <class MATRIX=SparseMatrix<double>, class VECTOR=Vector<double> >
class PreconditionChebyshev : public Subscriptor
{
-public:
- /**
- * Standardized data struct to
- * pipe additional parameters
- * to the preconditioner.
- */
- struct AdditionalData
- {
- /**
- * Constructor.
- */
- AdditionalData (const unsigned int degree = 0,
- const double smoothing_range = 0.,
- const bool nonzero_starting = false,
- const unsigned int eig_cg_n_iterations = 8,
- const double eig_cg_residual = 1e-2);
-
- /**
- * This determines the degree of the
- * Chebyshev polynomial. The degree
- * of the polynomial gives the number
- * of matrix-vector products to be
- * performed for one application of
- * the vmult() operation. Degree zero
- * corresponds to a damped Jacobi
- * method.
- */
- unsigned int degree;
-
- /**
- * This sets the range between the
- * largest eigenvalue in the matrix
- * and the smallest eigenvalue to be
- * treated. If the parameter is zero,
- * an estimate for the largest and
- * for the smallest eigenvalue will
- * be calculated
- * internally. Otherwise, the
- * Chebyshev polynomial will act in
- * the interval
- * $[\lambda_\mathrm{max}/
- * \tt{smoothing\_range},
- * \lambda_\mathrm{max}]$, where
- * $\lambda_\mathrm{max}$ is an
- * estimate of the maximum eigenvalue
- * of the matrix. A choice of
- * <tt>smoothing_range</tt> between 5
- * and 20 is useful in case the
- * preconditioner is used as a
- * smoother in multigrid.
- */
- double smoothing_range;
-
- /**
- * When this flag is set to
- * <tt>true</tt>, it enables the
- * method <tt>vmult(dst, src)</tt> to
- * use non-zero data in the vector
- * <tt>dst</tt>, appending to it the
- * Chebyshev corrections. This can be
- * useful in some situations
- * (e.g. when used for high-frequency
- * error smoothing in a multigrid
- * algorithm), but not the way the
- * solver classes expect a
- * preconditioner to work (where one
- * ignores the content in
- * <tt>dst</tt> for the
- * preconditioner application).
- */
- bool nonzero_starting;
-
- /**
- * Maximum number of CG iterations
- * performed for finding the maximum
- * eigenvalue.
- */
- unsigned int eig_cg_n_iterations;
-
- /**
- * Tolerance for CG iterations
- * performed for finding the maximum
- * eigenvalue.
- */
- double eig_cg_residual;
-
- /**
- * Stores the inverse of the diagonal
- * of the underlying matrix.
- */
- VECTOR matrix_diagonal_inverse;
- };
-
- PreconditionChebyshev ();
-
- /**
- * Initialize function. Takes the
- * matrix which is used to form the
- * preconditioner, and additional
- * flags if there are any. This
- * function works only if the input
- * matrix has an operator
- * <tt>el(i,i)</tt> for accessing all
- * the elements in the
- * diagonal. Alternatively, the
- * diagonal can be supplied with the
- * help of the AdditionalData field.
- *
- * This function calculates an
- * estimate of the eigenvalue range
- * of the matrix weighted by its
- * diagonal using a modified CG
- * iteration.
- */
- void initialize (const MATRIX &matrix,
- const AdditionalData &additional_data = AdditionalData());
-
- /**
- * Computes the action of the
- * preconditioner on <tt>src</tt>,
- * storing the result in
- * <tt>dst</tt>.
- */
- void vmult (VECTOR &dst,
- const VECTOR &src) const;
-
- /**
- * Computes the action of the
- * transposed preconditioner on
- * <tt>src</tt>, storing the result
- * in <tt>dst</tt>.
- */
- void Tvmult (VECTOR &dst,
- const VECTOR &src) const;
-
- /**
- * Resets the preconditioner.
- */
- void clear ();
-
-private:
-
- /**
- * A pointer to the underlying
- * matrix.
- */
+ public:
+ /**
+ * Standardized data struct to
+ * pipe additional parameters
+ * to the preconditioner.
+ */
+ struct AdditionalData
+ {
+ /**
+ * Constructor.
+ */
+ AdditionalData (const unsigned int degree = 0,
+ const double smoothing_range = 0.,
+ const bool nonzero_starting = false,
+ const unsigned int eig_cg_n_iterations = 8,
+ const double eig_cg_residual = 1e-2);
+
+ /**
+ * This determines the degree of the
+ * Chebyshev polynomial. The degree
+ * of the polynomial gives the number
+ * of matrix-vector products to be
+ * performed for one application of
+ * the vmult() operation. Degree zero
+ * corresponds to a damped Jacobi
+ * method.
+ */
+ unsigned int degree;
+
+ /**
+ * This sets the range between the
+ * largest eigenvalue in the matrix
+ * and the smallest eigenvalue to be
+ * treated. If the parameter is zero,
+ * an estimate for the largest and
+ * for the smallest eigenvalue will
+ * be calculated
+ * internally. Otherwise, the
+ * Chebyshev polynomial will act in
+ * the interval
+ * $[\lambda_\mathrm{max}/
+ * \tt{smoothing\_range},
+ * \lambda_\mathrm{max}]$, where
+ * $\lambda_\mathrm{max}$ is an
+ * estimate of the maximum eigenvalue
+ * of the matrix. A choice of
+ * <tt>smoothing_range</tt> between 5
+ * and 20 is useful in case the
+ * preconditioner is used as a
+ * smoother in multigrid.
+ */
+ double smoothing_range;
+
+ /**
+ * When this flag is set to
+ * <tt>true</tt>, it enables the
+ * method <tt>vmult(dst, src)</tt> to
+ * use non-zero data in the vector
+ * <tt>dst</tt>, appending to it the
+ * Chebyshev corrections. This can be
+ * useful in some situations
+ * (e.g. when used for high-frequency
+ * error smoothing in a multigrid
+ * algorithm), but not the way the
+ * solver classes expect a
+ * preconditioner to work (where one
+ * ignores the content in
+ * <tt>dst</tt> for the
+ * preconditioner application).
+ */
+ bool nonzero_starting;
+
+ /**
+ * Maximum number of CG iterations
+ * performed for finding the maximum
+ * eigenvalue.
+ */
+ unsigned int eig_cg_n_iterations;
+
+ /**
+ * Tolerance for CG iterations
+ * performed for finding the maximum
+ * eigenvalue.
+ */
+ double eig_cg_residual;
+
+ /**
+ * Stores the inverse of the diagonal
+ * of the underlying matrix.
+ */
+ VECTOR matrix_diagonal_inverse;
+ };
+
+ PreconditionChebyshev ();
+
+ /**
+ * Initialize function. Takes the
+ * matrix which is used to form the
+ * preconditioner, and additional
+ * flags if there are any. This
+ * function works only if the input
+ * matrix has an operator
+ * <tt>el(i,i)</tt> for accessing all
+ * the elements in the
+ * diagonal. Alternatively, the
+ * diagonal can be supplied with the
+ * help of the AdditionalData field.
+ *
+ * This function calculates an
+ * estimate of the eigenvalue range
+ * of the matrix weighted by its
+ * diagonal using a modified CG
+ * iteration.
+ */
+ void initialize (const MATRIX &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
+ /**
+ * Computes the action of the
+ * preconditioner on <tt>src</tt>,
+ * storing the result in
+ * <tt>dst</tt>.
+ */
+ void vmult (VECTOR &dst,
+ const VECTOR &src) const;
+
+ /**
+ * Computes the action of the
+ * transposed preconditioner on
+ * <tt>src</tt>, storing the result
+ * in <tt>dst</tt>.
+ */
+ void Tvmult (VECTOR &dst,
+ const VECTOR &src) const;
+
+ /**
+ * Resets the preconditioner.
+ */
+ void clear ();
+
+ private:
+
+ /**
+ * A pointer to the underlying
+ * matrix.
+ */
SmartPointer<const MATRIX,PreconditionChebyshev<MATRIX,VECTOR> > matrix_ptr;
- /**
- * Internal vector used for the
- * <tt>vmult</tt> operation.
- */
- mutable VECTOR update1;
-
- /**
- * Internal vector used for the
- * <tt>vmult</tt> operation.
- */
- mutable VECTOR update2;
-
- /**
- * Stores the additional data
- * provided to the initialize
- * function.
- */
- AdditionalData data;
-
- /**
- * Average of the largest and
- * smallest eigenvalue under
- * consideration.
- */
- double theta;
-
- /**
- * Half the interval length between
- * the largest and smallest
- * eigenvalue under consideration.
- */
- double delta;
-
- /**
- * Stores whether the preconditioner
- * has been set up.
- */
- bool is_initialized;
+ /**
+ * Internal vector used for the
+ * <tt>vmult</tt> operation.
+ */
+ mutable VECTOR update1;
+
+ /**
+ * Internal vector used for the
+ * <tt>vmult</tt> operation.
+ */
+ mutable VECTOR update2;
+
+ /**
+ * Stores the additional data
+ * provided to the initialize
+ * function.
+ */
+ AdditionalData data;
+
+ /**
+ * Average of the largest and
+ * smallest eigenvalue under
+ * consideration.
+ */
+ double theta;
+
+ /**
+ * Half the interval length between
+ * the largest and smallest
+ * eigenvalue under consideration.
+ */
+ double delta;
+
+ /**
+ * Stores whether the preconditioner
+ * has been set up.
+ */
+ bool is_initialized;
};
inline void
PreconditionRichardson::initialize (
- const PreconditionRichardson::AdditionalData parameters)
+ const PreconditionRichardson::AdditionalData ¶meters)
{
relaxation = parameters.relaxation;
}
inline void
PreconditionRichardson::initialize (
const MATRIX&,
- const PreconditionRichardson::AdditionalData parameters)
+ const PreconditionRichardson::AdditionalData ¶meters)
{
relaxation = parameters.relaxation;
}
template <class MATRIX>
inline void
PreconditionRelaxation<MATRIX>::initialize (const MATRIX &rA,
- AdditionalData parameters)
+ const AdditionalData ¶meters)
{
A = &rA;
relaxation = parameters.relaxation;
template <class MATRIX>
inline void
PreconditionSSOR<MATRIX>::initialize (const MATRIX &rA,
- typename BaseClass::AdditionalData parameters)
+ const typename BaseClass::AdditionalData ¶meters)
{
this->PreconditionRelaxation<MATRIX>::initialize (rA, parameters);
// in case we have a SparseMatrix class,
// we can extract information about the
// diagonal.
- const SparseMatrix<typename MATRIX::value_type> * mat =
- dynamic_cast<const SparseMatrix<typename MATRIX::value_type> *>(&*this->A);
+ const SparseMatrix<typename MATRIX::value_type> * mat =
+ dynamic_cast<const SparseMatrix<typename MATRIX::value_type> *>(&*this->A);
// calculate the positions first after
// the diagonal.
if (mat != 0)
{
- const std::size_t * rowstart_ptr =
+ const std::size_t * rowstart_ptr =
mat->get_sparsity_pattern().get_rowstart_indices();
- const unsigned int * const colnums =
+ const unsigned int * const colnums =
mat->get_sparsity_pattern().get_column_numbers();
const unsigned int n = this->A->n();
pos_right_of_diagonal.resize(n);
// note: the first entry in each
// line denotes the diagonal element,
// which we need not check.
- pos_right_of_diagonal[row] =
+ pos_right_of_diagonal[row] =
std::lower_bound (&colnums[*rowstart_ptr+1],
&colnums[*(rowstart_ptr+1)],
row)
const MATRIX &rA,
const std::vector<unsigned int> &p,
const std::vector<unsigned int> &ip,
- typename PreconditionRelaxation<MATRIX>::AdditionalData parameters)
+ const typename PreconditionRelaxation<MATRIX>::AdditionalData ¶meters)
{
permutation = &p;
inverse_permutation = &ip;
{
Assert (solver !=0 && matrix != 0 && precondition != 0,
ExcNotInitialized());
-
+
solver->solve(*matrix, dst, src, *precondition);
}
template <class MATRIX, class VECTOR>
inline
-void
+void
PreconditionChebyshev<MATRIX,VECTOR>::initialize (const MATRIX &matrix,
const AdditionalData &additional_data)
{
// calculate largest eigenvalue using a
// hand-tuned CG iteration on the matrix
- // weighted by its diagonal. we start
- // with a vector that consists of ones
+ // weighted by its diagonal. we start
+ // with a vector that consists of ones
// only, weighted by the length.
- //
+ //
// TODO: can we obtain this with the
// regular CG implementation? we would need
// to read the logfile in that case, which
// include a safety factor since the CG
// method will in general not be converged
const double beta = 1.2 * max_eigenvalue;
- const double alpha = (data.smoothing_range > 0 ?
+ const double alpha = (data.smoothing_range > 0 ?
max_eigenvalue / data.smoothing_range :
max_eigenvalue / min_eigenvalue);
delta = (beta-alpha)*0.5;
template <class MATRIX, class VECTOR>
inline
-void
+void
PreconditionChebyshev<MATRIX,VECTOR>::vmult (VECTOR &dst,
const VECTOR &src) const
{
template <class MATRIX, class VECTOR>
inline
-void
+void
PreconditionChebyshev<MATRIX,VECTOR>::Tvmult (VECTOR &dst,
const VECTOR &src) const
{
projects / dealii.git / commitdiff