From 9df3783e3f8ee0171ab79ce38a91506a1966bf64 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 28 Aug 2017 13:23:42 -0600 Subject: [PATCH] Convert more places to VectorMemory::Pointer. --- include/deal.II/lac/eigen.h | 7 ++-- include/deal.II/lac/filtered_matrix.h | 12 +++---- include/deal.II/lac/householder.h | 10 +++--- include/deal.II/lac/solver_bicgstab.h | 46 +++++++++------------------ source/lac/trilinos_sparse_matrix.cc | 32 +++++-------------- 5 files changed, 33 insertions(+), 74 deletions(-) diff --git a/include/deal.II/lac/eigen.h b/include/deal.II/lac/eigen.h index ddd9e9a2dd..e9da6fa6c2 100644 --- a/include/deal.II/lac/eigen.h +++ b/include/deal.II/lac/eigen.h @@ -345,10 +345,10 @@ EigenInverse::solve (double &value, unsigned int goal = additional_data.start_adaption; // Auxiliary vector - VectorType *Vy = this->memory.alloc (); + typename VectorMemory::Pointer Vy (this->memory); VectorType &y = *Vy; y.reinit (x); - VectorType *Vr = this->memory.alloc (); + typename VectorMemory::Pointer Vr (this->memory); VectorType &r = *Vr; r.reinit (x); @@ -416,9 +416,6 @@ EigenInverse::solve (double &value, old_value = value; } - this->memory.free(Vy); - this->memory.free(Vr); - deallog.pop(); // in case of failure: throw diff --git a/include/deal.II/lac/filtered_matrix.h b/include/deal.II/lac/filtered_matrix.h index 2280233363..bb946d6509 100644 --- a/include/deal.II/lac/filtered_matrix.h +++ b/include/deal.II/lac/filtered_matrix.h @@ -900,7 +900,7 @@ FilteredMatrix::vmult (VectorType &dst, const VectorType &src) const if (!expect_constrained_source) { GrowingVectorMemory mem; - VectorType *tmp_vector = mem.alloc(); + typename VectorMemory::Pointer tmp_vector (mem); // first copy over src vector and // pre-filter tmp_vector->reinit(src, true); @@ -908,7 +908,6 @@ FilteredMatrix::vmult (VectorType &dst, const VectorType &src) const pre_filter (*tmp_vector); // then let matrix do its work matrix->vmult (dst, *tmp_vector); - mem.free(tmp_vector); } else { @@ -929,7 +928,7 @@ FilteredMatrix::Tvmult (VectorType &dst, const VectorType &src) cons if (!expect_constrained_source) { GrowingVectorMemory mem; - VectorType *tmp_vector = mem.alloc(); + typename VectorMemory::Pointer tmp_vector (mem); // first copy over src vector and // pre-filter tmp_vector->reinit(src, true); @@ -937,7 +936,6 @@ FilteredMatrix::Tvmult (VectorType &dst, const VectorType &src) cons pre_filter (*tmp_vector); // then let matrix do its work matrix->Tvmult (dst, *tmp_vector); - mem.free(tmp_vector); } else { @@ -958,7 +956,7 @@ FilteredMatrix::vmult_add (VectorType &dst, const VectorType &src) c if (!expect_constrained_source) { GrowingVectorMemory mem; - VectorType *tmp_vector = mem.alloc(); + typename VectorMemory::Pointer tmp_vector (mem); // first copy over src vector and // pre-filter tmp_vector->reinit(src, true); @@ -966,7 +964,6 @@ FilteredMatrix::vmult_add (VectorType &dst, const VectorType &src) c pre_filter (*tmp_vector); // then let matrix do its work matrix->vmult_add (dst, *tmp_vector); - mem.free(tmp_vector); } else { @@ -987,7 +984,7 @@ FilteredMatrix::Tvmult_add (VectorType &dst, const VectorType &src) if (!expect_constrained_source) { GrowingVectorMemory mem; - VectorType *tmp_vector = mem.alloc(); + typename VectorMemory::Pointer tmp_vector (mem); // first copy over src vector and // pre-filter tmp_vector->reinit(src, true); @@ -995,7 +992,6 @@ FilteredMatrix::Tvmult_add (VectorType &dst, const VectorType &src) pre_filter (*tmp_vector); // then let matrix do its work matrix->Tvmult_add (dst, *tmp_vector); - mem.free(tmp_vector); } else { diff --git a/include/deal.II/lac/householder.h b/include/deal.II/lac/householder.h index cec7b58eb5..d1272117fc 100644 --- a/include/deal.II/lac/householder.h +++ b/include/deal.II/lac/householder.h @@ -208,7 +208,7 @@ Householder::least_squares (Vector &dst, const size_type m = this->m(), n = this->n(); GrowingVectorMemory > mem; - Vector *aux = mem.alloc(); + typename VectorMemory >::Pointer aux (mem); aux->reinit(src, true); *aux = src; // m > n, m = src.n, n = dst.n @@ -235,11 +235,11 @@ Householder::least_squares (Vector &dst, // Compute solution this->backward(dst, *aux); - mem.free(aux); - return std::sqrt(sum); } + + template template double @@ -253,7 +253,7 @@ Householder::least_squares (BlockVector &dst, const size_type m = this->m(), n = this->n(); GrowingVectorMemory > mem; - BlockVector *aux = mem.alloc(); + typename VectorMemory >::Pointer aux (mem); aux->reinit(src, true); *aux = src; // m > n, m = src.n, n = dst.n @@ -289,8 +289,6 @@ Householder::least_squares (BlockVector &dst, //to the BlockVector dst = v_dst; - mem.free(aux); - return std::sqrt(sum); } diff --git a/include/deal.II/lac/solver_bicgstab.h b/include/deal.II/lac/solver_bicgstab.h index be4aaa2467..ede0758a68 100644 --- a/include/deal.II/lac/solver_bicgstab.h +++ b/include/deal.II/lac/solver_bicgstab.h @@ -160,31 +160,31 @@ protected: /** * Auxiliary vector. */ - VectorType *Vr; + typename VectorMemory::Pointer Vr; /** * Auxiliary vector. */ - VectorType *Vrbar; + typename VectorMemory::Pointer Vrbar; /** * Auxiliary vector. */ - VectorType *Vp; + typename VectorMemory::Pointer Vp; /** * Auxiliary vector. */ - VectorType *Vy; + typename VectorMemory::Pointer Vy; /** * Auxiliary vector. */ - VectorType *Vz; + typename VectorMemory::Pointer Vz; /** * Auxiliary vector. */ - VectorType *Vt; + typename VectorMemory::Pointer Vt; /** * Auxiliary vector. */ - VectorType *Vv; + typename VectorMemory::Pointer Vv; /** * Right hand side vector. */ @@ -296,15 +296,6 @@ SolverBicgstab::SolverBicgstab (SolverControl &cn, const AdditionalData &data) : Solver(cn), - Vx(nullptr), - Vr(nullptr), - Vrbar(nullptr), - Vp(nullptr), - Vy(nullptr), - Vz(nullptr), - Vt(nullptr), - Vv(nullptr), - Vb(nullptr), alpha(0.), beta(0.), omega(0.), @@ -459,19 +450,20 @@ SolverBicgstab::solve(const MatrixType &A, const PreconditionerType &preconditioner) { deallog.push("Bicgstab"); - Vr = this->memory.alloc(); + Vr = typename VectorMemory::Pointer(this->memory); + Vrbar = typename VectorMemory::Pointer(this->memory); + Vp = typename VectorMemory::Pointer(this->memory); + Vy = typename VectorMemory::Pointer(this->memory); + Vz = typename VectorMemory::Pointer(this->memory); + Vt = typename VectorMemory::Pointer(this->memory); + Vv = typename VectorMemory::Pointer(this->memory); + Vr->reinit(x, true); - Vrbar = this->memory.alloc(); Vrbar->reinit(x, true); - Vp = this->memory.alloc(); Vp->reinit(x, true); - Vy = this->memory.alloc(); Vy->reinit(x, true); - Vz = this->memory.alloc(); Vz->reinit(x, true); - Vt = this->memory.alloc(); Vt->reinit(x, true); - Vv = this->memory.alloc(); Vv->reinit(x, true); Vx = &x; @@ -496,14 +488,6 @@ SolverBicgstab::solve(const MatrixType &A, } while (state.breakdown == true); - this->memory.free(Vr); - this->memory.free(Vrbar); - this->memory.free(Vp); - this->memory.free(Vy); - this->memory.free(Vz); - this->memory.free(Vt); - this->memory.free(Vv); - deallog.pop(); // in case of failure: throw exception diff --git a/source/lac/trilinos_sparse_matrix.cc b/source/lac/trilinos_sparse_matrix.cc index 3a4e008675..c233063843 100644 --- a/source/lac/trilinos_sparse_matrix.cc +++ b/source/lac/trilinos_sparse_matrix.cc @@ -2973,7 +2973,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // Initialize intermediate vector const Epetra_Map &first_op_init_map = first_op.OperatorDomainMap(); @@ -2997,8 +2997,6 @@ namespace TrilinosWrappers first_op.Apply(tril_src, tril_dst); const int ierr = tril_dst.Update (1.0, tril_int, 1.0); AssertThrow (ierr == 0, ExcTrilinosError(ierr)); - - vector_memory.free(i); }; return_op.Tvmult = [first_op, second_op](Domain &tril_dst, @@ -3007,7 +3005,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // These operators may themselves be transposed or not, so we let them // decide what the intended outcome is @@ -3042,8 +3040,6 @@ namespace TrilinosWrappers // Reset transpose flag const_cast(first_op).transpose(); const_cast(second_op).transpose(); - - vector_memory.free(i); }; return_op.inv_vmult = [first_op, second_op](Domain &tril_dst, @@ -3052,7 +3048,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // Initialize intermediate vector const Epetra_Map &first_op_init_map = first_op.OperatorRangeMap(); @@ -3076,8 +3072,6 @@ namespace TrilinosWrappers first_op.ApplyInverse(tril_src, tril_dst); const int ierr = tril_dst.Update (1.0, tril_int, 1.0); AssertThrow (ierr == 0, ExcTrilinosError(ierr)); - - vector_memory.free(i); }; return_op.inv_Tvmult = [first_op, second_op](Range &tril_dst, @@ -3086,7 +3080,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // These operators may themselves be transposed or not, so we let them // decide what the intended outcome is @@ -3121,8 +3115,6 @@ namespace TrilinosWrappers // Reset transpose flag const_cast(first_op).transpose(); const_cast(second_op).transpose(); - - vector_memory.free(i); }; return return_op; @@ -3150,7 +3142,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // Initialize intermediate vector const Epetra_Map &first_op_init_map = first_op.OperatorDomainMap(); @@ -3172,8 +3164,6 @@ namespace TrilinosWrappers // decide what the intended outcome is second_op.Apply(tril_src, tril_int); first_op.Apply(tril_int, tril_dst); - - vector_memory.free(i); }; return_op.Tvmult = [first_op, second_op](Domain &tril_dst, @@ -3182,7 +3172,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // These operators may themselves be transposed or not, so we let them // decide what the intended outcome is @@ -3214,8 +3204,6 @@ namespace TrilinosWrappers // Reset transpose flag const_cast(first_op).transpose(); const_cast(second_op).transpose(); - - vector_memory.free(i); }; return_op.inv_vmult = [first_op, second_op](Domain &tril_dst, @@ -3224,7 +3212,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // Initialize intermediate vector const Epetra_Map &first_op_init_map = first_op.OperatorRangeMap(); @@ -3246,8 +3234,6 @@ namespace TrilinosWrappers // and the same order as Tvmult first_op.ApplyInverse(tril_src, tril_int); second_op.ApplyInverse(tril_int, tril_dst); - - vector_memory.free(i); }; return_op.inv_Tvmult = [first_op, second_op](Range &tril_dst, @@ -3256,7 +3242,7 @@ namespace TrilinosWrappers // Duplicated from LinearOperator::operator* // TODO: Template the constructor on GrowingVectorMemory vector type? GrowingVectorMemory vector_memory; - GVMVectorType *i = vector_memory.alloc(); + VectorMemory::Pointer i (vector_memory); // These operators may themselves be transposed or not, so we let them // decide what the intended outcome is @@ -3291,8 +3277,6 @@ namespace TrilinosWrappers // Reset transpose flag const_cast(first_op).transpose(); const_cast(second_op).transpose(); - - vector_memory.free(i); }; return return_op; -- 2.39.5