From: guido Date: Fri, 24 May 2002 13:46:59 +0000 (+0000) Subject: dynamic vector allocation X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=11f4db1a588c161b952f24f0ca8acc4ca936e67b;p=dealii-svn.git dynamic vector allocation git-svn-id: https://svn.dealii.org/trunk@5875 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/solver_gmres.h b/deal.II/lac/include/lac/solver_gmres.h index e20129f2bc..ad7932bce0 100644 --- a/deal.II/lac/include/lac/solver_gmres.h +++ b/deal.II/lac/include/lac/solver_gmres.h @@ -158,25 +158,109 @@ class SolverGMRES : public Solver void givens_rotation (Vector &h, Vector &b, Vector &ci, Vector &si, int col) const; - /** - * Projected system matrix - */ - FullMatrix H; - /** - * Auxiliary matrix for inverting @p{H} - */ - FullMatrix H1; - + /** + * Projected system matrix + */ + FullMatrix H; + /** + * Auxiliary matrix for inverting @p{H} + */ + FullMatrix H1; + private: - /** - * No copy constructor. - */ - SolverGMRES (const SolverGMRES&); + /** + * Class to hold temporary + * vectors. This class + * automatically allocates a new + * vector, once it is needed. + * + * A future version should also + * be able to shift through + * vectors automatically, + * avoiding restart. + */ + + class TmpVectors + { + public: + /** + * Constructor. Prepares an + * array of @p{VECTOR} of + * length @p{max_size}. + */ + TmpVectors(unsigned int max_size, + VectorMemory& vmem); + /** + * Delete all allocated vectors. + */ + ~TmpVectors(); + /** + * Get vector number + * @p{i}. Allocate it if + * necessary. + */ + VECTOR& operator[] (unsigned int i); + + private: + /** + * Pool were vectors are obtained from. + */ + VectorMemory& mem; + /** + * Field for storing the vectors. + */ + std::vector data; + /** + * Offset of the first + * vector. This is for later + * when vector rotation will + * be implemented. + */ + unsigned int offset; + }; + /** + * No copy constructor. + */ + SolverGMRES (const SolverGMRES&); }; /* --------------------- Inline and template functions ------------------- */ +template +inline +SolverGMRES::TmpVectors::TmpVectors (unsigned int max_size, + VectorMemory& vmem) + : + mem(vmem), + data (max_size, 0), + offset(0) +{} + + +template +inline +SolverGMRES::TmpVectors::~TmpVectors () +{ + for (typename std::vector::iterator v = data.begin(); + v != data.end(); ++v) + if (*v != 0) + mem.free(*v); +} + + +template +inline VECTOR& +SolverGMRES::TmpVectors::operator[] (unsigned int i) +{ + Assert (i+offset SolverGMRES::SolverGMRES (SolverControl &cn, @@ -238,22 +322,10 @@ SolverGMRES::solve (const MATRIX &A, deallog.push("GMRES"); const unsigned int n_tmp_vectors = additional_data.max_n_tmp_vectors; - // allocate an array of n_tmp_vectors - // temporary vectors from the VectorMemory - // object; resize them but do not set their - // values since they will be overwritten soon - // anyway. - - // This is really bad since vectors - // should only be allocated if - // really needed. (GK) - std::vector tmp_vectors (n_tmp_vectors, 0); - for (unsigned int tmp=0; tmpreinit (x, true); - }; - + // Generate an object where basis + // vectors are stored. + TmpVectors tmp_vectors (n_tmp_vectors, memory); + // number of the present iteration; this // number is not reset to zero upon a // restart @@ -282,8 +354,8 @@ unsigned int dim = 0; const bool left_precondition = ! additional_data.right_preconditioning; // define two aliases - VECTOR &v = *tmp_vectors[0]; - VECTOR &p = *tmp_vectors[n_tmp_vectors-1]; + VECTOR &v = tmp_vectors[0]; + VECTOR &p = tmp_vectors[n_tmp_vectors-1]; /////////////////////////////////// @@ -339,14 +411,14 @@ unsigned int dim = 0; { ++accumulated_iterations; // yet another alias - VECTOR& vv = *tmp_vectors[inner_iteration+1]; + VECTOR& vv = tmp_vectors[inner_iteration+1]; if (left_precondition) { - A.vmult(p, *tmp_vectors[inner_iteration]); + A.vmult(p, tmp_vectors[inner_iteration]); precondition.vmult(vv,p); } else { - precondition.vmult(p,*tmp_vectors[inner_iteration]); + precondition.vmult(p, tmp_vectors[inner_iteration]); A.vmult(vv,p); }; @@ -355,8 +427,8 @@ unsigned int dim = 0; /* Orthogonalization */ for (unsigned int i=0 ; i