From f67261612d80dd71a4746a0f6d12fe38295feb34 Mon Sep 17 00:00:00 2001 From: guido Date: Thu, 22 Nov 2001 17:25:39 +0000 Subject: [PATCH] parameterized reinit functions and some constructor fixes git-svn-id: https://svn.dealii.org/trunk@5254 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/block_vector.h | 9 ++-- .../lac/include/lac/block_vector.templates.h | 5 ++- deal.II/lac/include/lac/solver_richardson.h | 45 ++++++++++++++----- deal.II/lac/include/lac/vector.h | 13 ++++-- deal.II/lac/include/lac/vector.templates.h | 3 +- deal.II/lac/source/block_vector.cc | 9 ++++ deal.II/lac/source/vector.cc | 4 ++ 7 files changed, 69 insertions(+), 19 deletions(-) diff --git a/deal.II/lac/include/lac/block_vector.h b/deal.II/lac/include/lac/block_vector.h index 5aeff86099..a9fd7c994c 100644 --- a/deal.II/lac/include/lac/block_vector.h +++ b/deal.II/lac/include/lac/block_vector.h @@ -639,8 +639,8 @@ class BlockVector * use blocks of different * sizes. */ - BlockVector (unsigned int num_blocks = 0, - unsigned int block_size = 0); + explicit BlockVector (unsigned int num_blocks = 0, + unsigned int block_size = 0); /** * Copy-Constructor. Dimension set to @@ -748,7 +748,8 @@ class BlockVector * this function is the same as calling * @p{reinit (V.size(), fast)}. */ - void reinit (const BlockVector &V, + template + void reinit (const BlockVector &V, const bool fast=false); /** @@ -1129,6 +1130,8 @@ class BlockVector */ friend class iterator; friend class const_iterator; + + template friend class BlockVector; }; diff --git a/deal.II/lac/include/lac/block_vector.templates.h b/deal.II/lac/include/lac/block_vector.templates.h index 59d31a3a93..2b53f56671 100644 --- a/deal.II/lac/include/lac/block_vector.templates.h +++ b/deal.II/lac/include/lac/block_vector.templates.h @@ -96,7 +96,8 @@ void BlockVector::reinit (const std::vector &n, template -void BlockVector::reinit (const BlockVector& v, +template +void BlockVector::reinit (const BlockVector& v, const bool fast) { block_indices = v.block_indices; @@ -493,7 +494,7 @@ BlockVector::operator = (const BlockVector& v) template template BlockVector& -BlockVector::operator = (const BlockVector< Number2>& v) +BlockVector::operator = (const BlockVector& v) { Assert (num_blocks == v.num_blocks, ExcDimensionMismatch(num_blocks, v.num_blocks)); diff --git a/deal.II/lac/include/lac/solver_richardson.h b/deal.II/lac/include/lac/solver_richardson.h index 7a825b8078..a59da8dc8a 100644 --- a/deal.II/lac/include/lac/solver_richardson.h +++ b/deal.II/lac/include/lac/solver_richardson.h @@ -57,14 +57,21 @@ class SolverRichardson : public Solver * set the damping parameter * to one. */ - AdditionalData(double omega=1): - omega(omega) + AdditionalData(double omega=1, + bool use_preconditioned_residual = true): + omega(omega), + use_preconditioned_residual(use_preconditioned_residual) {}; /** * Relaxation parameter. */ double omega; + /** + * Parameter for stopping criterion. + */ + bool use_preconditioned_residual; + }; /** @@ -125,7 +132,7 @@ class SolverRichardson : public Solver * Implementation of the computation of * the norm of the residual. */ - virtual double criterion(); + virtual typename VECTOR::value_type criterion(); /** * Temporary vectors, allocated through @@ -151,7 +158,7 @@ class SolverRichardson : public Solver * norm of the residual vector and thus * the square root of the @p{res2} value. */ - double res; + typename VECTOR::value_type res; }; @@ -196,13 +203,31 @@ SolverRichardson::solve (const MATRIX &A, // but do it in 2 steps A.vmult(r,x); r.sadd(-1.,1.,b); - res=sqrt(r*r); - - conv = control().check (iter, criterion()); - if (conv != SolverControl::iterate) - break; + if (!additional_data.use_preconditioned_residual) + { + res = r*r; +// cout << '[' << res << ' '; + res=sqrt(res); +// cout << res << ' ' << r.l1_norm() << ']'; +// r.print(cout); + conv = control().check (iter, criterion()); + if (conv != SolverControl::iterate) + break; + } + precondition.vmult(d,r); + if (additional_data.use_preconditioned_residual) + { + res = d*d; +// cout << '[' << res << ' '; + res=sqrt(res); +// cout << res << ' ' << r.l1_norm() << ']'; +// r.print(cout); + conv = control().check (iter, criterion()); + if (conv != SolverControl::iterate) + break; + } x.add(additional_data.omega,d); print_vectors(iter,x,r,d); } @@ -281,7 +306,7 @@ SolverRichardson::print_vectors(const unsigned int, template -inline double +inline typename VECTOR::value_type SolverRichardson::criterion() { return res; diff --git a/deal.II/lac/include/lac/vector.h b/deal.II/lac/include/lac/vector.h index e4513eb27e..83e1363e1a 100644 --- a/deal.II/lac/include/lac/vector.h +++ b/deal.II/lac/include/lac/vector.h @@ -76,7 +76,11 @@ class Vector /** * Copy-Constructor. Dimension set to * that of V, all components are copied - * from V + * from V. + * + * We would like to make this + * constructor explicit, but STL + * insists on using it implicitly. */ Vector (const Vector& V); @@ -102,7 +106,7 @@ class Vector * Constructor. Set dimension to @p{n} and * initialize all elements with zero. */ - Vector (const unsigned int n); + explicit Vector (const unsigned int n); /** * Initialize the vector with a @@ -160,7 +164,8 @@ class Vector * this function is the same as calling * @p{reinit (V.size(), fast)}. */ - void reinit (const Vector &V, + template + void reinit (const Vector &V, const bool fast=false); /** @@ -535,6 +540,8 @@ class Vector * Pointer to the array of components. */ Number *val; + + template friend class Vector; }; diff --git a/deal.II/lac/include/lac/vector.templates.h b/deal.II/lac/include/lac/vector.templates.h index 842fdf818f..a1850e4fd7 100644 --- a/deal.II/lac/include/lac/vector.templates.h +++ b/deal.II/lac/include/lac/vector.templates.h @@ -71,7 +71,8 @@ Vector::Vector (const Vector& v) : template -void Vector::reinit (const Vector& v, const bool fast) +template +void Vector::reinit (const Vector& v, const bool fast) { reinit (v.size(), fast); }; diff --git a/deal.II/lac/source/block_vector.cc b/deal.II/lac/source/block_vector.cc index 601deb685b..b4cf91a8e9 100644 --- a/deal.II/lac/source/block_vector.cc +++ b/deal.II/lac/source/block_vector.cc @@ -15,7 +15,16 @@ // explicit instantiations template class BlockVector; +template BlockVector& BlockVector::operator=( + const BlockVector&); +template void BlockVector::reinit(const BlockVector&, bool); +template void BlockVector::reinit(const BlockVector&, bool); + template class BlockVector; +template BlockVector& BlockVector::operator=( + const BlockVector&); +template void BlockVector::reinit(const BlockVector&, bool); +template void BlockVector::reinit(const BlockVector&, bool); namespace BlockVectorIterators { diff --git a/deal.II/lac/source/vector.cc b/deal.II/lac/source/vector.cc index d1e40da33b..73a210a9be 100644 --- a/deal.II/lac/source/vector.cc +++ b/deal.II/lac/source/vector.cc @@ -19,11 +19,15 @@ template class Vector; template Vector& Vector::operator=(const Vector&); template double Vector::operator*(const Vector&) const; template double Vector::operator*(const Vector&) const; +template void Vector::reinit(const Vector&, bool); +template void Vector::reinit(const Vector&, bool); template class Vector; template Vector& Vector::operator=(const Vector&); template float Vector::operator*(const Vector&) const; template float Vector::operator*(const Vector&) const; +template void Vector::reinit(const Vector&, bool); +template void Vector::reinit(const Vector&, bool); // see the .h file for why these functions are disabled. // template Vector::Vector (const Vector& v); -- 2.39.5