From: Wolfgang Bangerth Date: Tue, 26 Sep 2017 18:59:47 +0000 (-0600) Subject: Avoid raw pointers in favor of VectorMemory::Pointer in SolverQMRS. X-Git-Tag: v9.0.0-rc1~1032^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9b32ed139383977aceeee37e1a2bf96363b1df72;p=dealii.git Avoid raw pointers in favor of VectorMemory::Pointer in SolverQMRS. --- diff --git a/include/deal.II/lac/solver_qmrs.h b/include/deal.II/lac/solver_qmrs.h index a9d7979e88..df6fef22f9 100644 --- a/include/deal.II/lac/solver_qmrs.h +++ b/include/deal.II/lac/solver_qmrs.h @@ -147,24 +147,6 @@ protected: */ virtual double criterion(); - /** - * Temporary vectors, allocated through the @p VectorMemory object at the - * start of the actual solution process and deallocated at the end. - */ - VectorType *Vv; - VectorType *Vp; - VectorType *Vq; - VectorType *Vt; - VectorType *Vd; - /** - * Iteration vector. - */ - VectorType *Vx; - /** - * RHS vector. - */ - const VectorType *Vb; - /** * Within the iteration loop, the square of the residual vector is stored in * this variable. The function @p criterion uses this variable to compute @@ -201,7 +183,14 @@ private: template IterationResult iterate (const MatrixType &A, - const PreconditionerType &preconditioner); + VectorType &x, + const VectorType &b, + const PreconditionerType &preconditioner, + VectorType &v, + VectorType &p, + VectorType &q, + VectorType &t, + VectorType &d); /** * Number of the current iteration (accumulated over restarts) @@ -273,15 +262,15 @@ SolverQMRS::solve (const MatrixType &A, { deallog.push("QMRS"); - // Memory allocation - Vv = this->memory.alloc(); - Vp = this->memory.alloc(); - Vq = this->memory.alloc(); - Vt = this->memory.alloc(); - Vd = this->memory.alloc(); + // temporary vectors, allocated through the @p VectorMemory object at the + // start of the actual solution process and deallocated at the end. + typename VectorMemory::Pointer Vv(this->memory); + typename VectorMemory::Pointer Vp(this->memory); + typename VectorMemory::Pointer Vq(this->memory); + typename VectorMemory::Pointer Vt(this->memory); + typename VectorMemory::Pointer Vd(this->memory); + - Vx = &x; - Vb = &b; // resize the vectors, but do not set // the values since they'd be overwritten // soon anyway. @@ -298,17 +287,10 @@ SolverQMRS::solve (const MatrixType &A, { if (step > 0) deallog << "Restart step " << step << std::endl; - state = iterate(A, preconditioner); + state = iterate(A, x, b, preconditioner, *Vv, *Vp, *Vq, *Vt, *Vd); } while (state.state == SolverControl::iterate); - // Deallocate Memory - this->memory.free(Vv); - this->memory.free(Vp); - this->memory.free(Vq); - this->memory.free(Vt); - this->memory.free(Vd); - // Output deallog.pop(); @@ -325,7 +307,14 @@ template template typename SolverQMRS::IterationResult SolverQMRS::iterate(const MatrixType &A, - const PreconditionerType &preconditioner) + VectorType &x, + const VectorType &b, + const PreconditionerType &preconditioner, + VectorType &v, + VectorType &p, + VectorType &q, + VectorType &t, + VectorType &d) { /* Remark: the matrix A in the article is the preconditioned matrix. * Therefore, we have to precondition x before we compute the first residual. @@ -335,15 +324,6 @@ SolverQMRS::iterate(const MatrixType &A, SolverControl::State state = SolverControl::iterate; - // define some aliases for simpler access - VectorType &v = *Vv; - VectorType &p = *Vp; - VectorType &q = *Vq; - VectorType &t = *Vt; - VectorType &d = *Vd; - VectorType &x = *Vx; - const VectorType &b = *Vb; - int it=0; double tau, rho, theta=0;