* double b, const Vector& z);
* };
* \end{verbatim}
+ *
+ * \subsection{AdditionalData}
+ * Several solvers need additional data, like the damping parameter #omega# of the
+ * #SolverRichardson# class or the maximum number of tmp vectors of the #SolverGMRES#.
+ * To have a standardized constructor for each solver class the #struct AdditionalData#
+ * has been introduced to each solver class. Some solver needs no additional data, like
+ * the SolverCG or SolverBicgstab. For these solvers the struct #AdditionalData# is empty
+ * and the calling
+ * of the constructor has not change as a default #AdditionalData# is set by default.
+ *
+ * Now the generating of a solver looks like
+ * \begin{verbatim}
+ * // GMRES with 50 tmp vectors
+ * SolverGMRES solver_gmres (solver_control, vector_memory,
+ * SolverGMRES::AdditionalData(50));
+ *
+ * // Richardson with omega=0.8
+ * SolverRichardson solver_richardson (solver_control, vector_memory,
+ * SolverGMRES::AdditionalData(0.8));
+ *
+ * // CG with default AdditionalData
+ * SolverCG solver_cg (solver_control, vector_memory);
+ * \end{verbatim}
*/
template <class Matrix, class Vector>
class Solver
/**
* Bicgstab algorithm by van der Vorst.
+ *
+ * The use of the #AdditionalData# struct is described in the #solver#
+ * base class.
*/
template<class Matrix, class Vector>
class SolverBicgstab //<Matrix,Vector>
: public Solver<Matrix,Vector>
{
public:
+ /**
+ * Standardized data struct to
+ * pipe additional data to the
+ * solver. This solver does not
+ * need additional data yet.
+ */
+ struct AdditionalData {};
+
/**
* Constructor.
*/
- SolverBicgstab(SolverControl &cn, VectorMemory<Vector> &mem);
+ SolverBicgstab(SolverControl &cn,
+ VectorMemory<Vector> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Solve primal problem only.
template<class Matrix, class Vector>
SolverBicgstab<Matrix, Vector>::SolverBicgstab(SolverControl &cn,
- VectorMemory<Vector> &mem)
- :
- Solver<Matrix,Vector>(cn,mem)
-{}
+ VectorMemory<Vector> &mem,
+ const AdditionalData &) :
+ Solver<Matrix,Vector>(cn,mem) {}
template<class Matrix, class Vector>
/**
* Preconditioned cg method.
* @author Original implementation by G. Kanschat, R. Becker and F.-T. Suttmeier, reworking and documentation by Wolfgang Bangerth
+ *
+ * The use of the #AdditionalData# struct is described in the #solver#
+ * base class.
*/
template<class Matrix, class Vector>
class SolverCG : public Solver<Matrix,Vector>
{
public:
+ /**
+ * Standardized data struct to
+ * pipe additional data to the
+ * solver. This solver does not
+ * need additional data yet.
+ */
+ struct AdditionalData {};
/**
* Constructor.
*/
- SolverCG (SolverControl &cn, VectorMemory<Vector> &mem) :
- Solver<Matrix,Vector>(cn,mem) {};
+ SolverCG (SolverControl &cn,
+ VectorMemory<Vector> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Solver method.
/*------------------------- Implementation ----------------------------*/
+template<class Matrix, class Vector>
+SolverCG<Matrix,Vector>::SolverCG(SolverControl &cn,
+ VectorMemory<Vector> &mem,
+ const AdditionalData &) :
+ Solver<Matrix,Vector>(cn,mem) {};
+
+
template<class Matrix, class Vector>
double SolverCG<Matrix,Vector>::criterion()
{
* computational speed against memory. One of the few other
* possibilities is to use a good preconditioner.
*
+ * The use of the #AdditionalData# struct is described in the #solver#
+ * base class.
+ *
* @author Wolfgang Bangerth
*/
template<class Matrix, class Vector>
class SolverGMRES : public Solver<Matrix, Vector>
{
public:
+ /**
+ * Standardized data struct to
+ * pipe additional data to the
+ * solver.
+ */
+ struct AdditionalData
+ {
+ AdditionalData(const unsigned int max_n_tmp_vectors=30):
+ max_n_tmp_vectors(max_n_tmp_vectors) {};
+
+ /**
+ * Maximum number of
+ * tmp vectors.
+ */
+ const unsigned int max_n_tmp_vectors;
+ };
+
/**
* Constructor.
*/
SolverGMRES (SolverControl &cn,
VectorMemory<Vector> &mem,
- const unsigned int max_n_tmp_vectors);
-
+ const AdditionalData &data=AdditionalData());
+
/**
* Solver method.
*/
<< "any results, and much more for reasonable ones.");
protected:
- const unsigned int max_n_tmp_vectors;
+ /**
+ * Includes the maximum number of
+ * tmp vectors.
+ */
+ AdditionalData additional_data;
/**
* Implementation of the computation of
template <class Matrix, class Vector>
SolverGMRES<Matrix,Vector>::SolverGMRES (SolverControl &cn,
VectorMemory<Vector> &mem,
- const unsigned int max_n_tmp_vectors) :
+ const AdditionalData &data) :
Solver<Matrix,Vector> (cn,mem),
- max_n_tmp_vectors (max_n_tmp_vectors)
+ additional_data(data)
{
- Assert (max_n_tmp_vectors >= 10, ExcTooFewTmpVectors (max_n_tmp_vectors));
+ Assert (additional_data.max_n_tmp_vectors >= 10,
+ ExcTooFewTmpVectors (additional_data.max_n_tmp_vectors));
};
// if the size of the matrix is very small,
// then only allocate a number of vectors
// which is needed
- const unsigned int n_tmp_vectors = (A.m()+3 > max_n_tmp_vectors ?
- max_n_tmp_vectors :
+ const unsigned int n_tmp_vectors = (A.m()+3 > additional_data.max_n_tmp_vectors ?
+ additional_data.max_n_tmp_vectors :
A.m()+3);
// allocate an array of n_tmp_vectors
* Implementation of the richardson iteration method. The stopping criterion
* is the norm of the residual.
*
+ * The use of the #AdditionalData# struct is described in the #solver#
+ * base class.
+ *
* @author Ralf Hartmann
*/
template<class Matrix, class Vector>
class SolverRichardson : public Solver<Matrix, Vector>
{
public:
+ /**
+ * Standardized data struct to
+ * pipe additional data to the
+ * solver.
+ */
+ struct AdditionalData
+ {
+ AdditionalData(double omega=1):
+ omega(omega) {};
+
+ /**
+ * Damping parameter.
+ */
+ double omega;
+ };
+
/**
* Constructor.
*/
- SolverRichardson (SolverControl &cn, VectorMemory<Vector> &mem) :
- Solver<Matrix,Vector> (cn,mem), omega(1.)
- {};
+ SolverRichardson (SolverControl &cn,
+ VectorMemory<Vector> &mem,
+ const AdditionalData &data=AdditionalData());
/**
* Solve $Ax=b$ for $x$.
Vector *Vr, *Vd;
/**
- * Damping-coefficient.
+ * Damping parameter.
*/
- double omega; ;
+ AdditionalData additional_data;
/**
* Within the iteration loop, the
/*----------------- Implementation of the Richardson Method ------------------*/
+template<class Matrix, class Vector>
+SolverRichardson<Matrix,Vector>::SolverRichardson(SolverControl &cn,
+ VectorMemory<Vector> &mem,
+ const AdditionalData &data):
+ Solver<Matrix,Vector> (cn,mem),
+ additional_data(data) {};
+
+
template<class Matrix,class Vector>
template<class Preconditioner>
Solver<Matrix,Vector>::ReturnState
break;
precondition(d,r);
- x.add(omega,d);
+ x.add(additional_data.omega,d);
}
// Deallocate Memory
return sqrt(res2);
}
+
template<class Matrix,class Vector>
inline void
SolverRichardson<Matrix,Vector>::set_omega(double om)
{
- omega=om;
+ additional_data.omega=om;
}