From 575c24430567a6104cdd9d106bcbab161dc4732c Mon Sep 17 00:00:00 2001
From: Wolfgang Bangerth
Date: Tue, 2 May 2000 08:59:26 +0000
Subject: [PATCH] Implement swap for vectors and block vectors.
git-svn-id: https://svn.dealii.org/trunk@2778 0785d39b-7218-0410-832d-ea1e28bc413d
---
deal.II/doc/news/2000/c-3-0.html | 11 +-
deal.II/lac/include/lac/block_vector.h | 116 ++++++++++++------
.../lac/include/lac/block_vector.templates.h | 15 +++
deal.II/lac/include/lac/solver_minres.h | 31 +++--
deal.II/lac/include/lac/vector.h | 111 +++++++++++------
deal.II/lac/include/lac/vector.templates.h | 33 +++--
6 files changed, 227 insertions(+), 90 deletions(-)
diff --git a/deal.II/doc/news/2000/c-3-0.html b/deal.II/doc/news/2000/c-3-0.html
index 92f71e8668..b38c17af0f 100644
--- a/deal.II/doc/news/2000/c-3-0.html
+++ b/deal.II/doc/news/2000/c-3-0.html
@@ -73,11 +73,20 @@
- New: Solver*
: There is a virtual function
+ New: SolverXX
: There is a virtual function
print_vectors
called in every step. It is void in the
solver itself but can be used to print intermediate iteration
vectors.
+
+
+ New: there are now functions Vector::swap
+ and BlockVector::swap
, as well as
+ global functions swap(u,v)
that
+ exchange the data of two vectors without needing a temporary
+ vector and without copying around data. Their run-time is
+ therefore independent of the length of the vectors.
+
diff --git a/deal.II/lac/include/lac/block_vector.h b/deal.II/lac/include/lac/block_vector.h
index d192382e2d..30d81f8bbf 100644
--- a/deal.II/lac/include/lac/block_vector.h
+++ b/deal.II/lac/include/lac/block_vector.h
@@ -77,14 +77,36 @@ class BlockVector
~BlockVector ();
/**
- * Access to a single block.
- */
- Vector& block(unsigned int i);
+ * Change the dimension of the vector to
+ * #N#. The reserved memory for this vector
+ * remains unchanged if possible, to make
+ * things faster, but this may waste some
+ * memory, so take this in the back of your
+ * head.
+ * However, if #N==0# all memory is freed,
+ * i.e. if you want to resize the vector
+ * and release the memory not needed, you
+ * have to first call #reinit(0)# and then
+ * #reinit(N)#. This cited behaviour is
+ * analogous to that of the STL containers.
+ *
+ * On #fast==false#, the vector is filled by
+ * zeros.
+ */
+ void reinit (const vector& N,
+ const bool fast=false);
/**
- * Read-only access to a single block.
+ * Change the dimension to that of the
+ * vector #V#. The same applies as for
+ * the other #reinit# function.
+ *
+ * The elements of #V# are not copied, i.e.
+ * this function is the same as calling
+ * #reinit (V.size(), fast)#.
*/
- const Vector& block(unsigned int i) const;
+ void reinit (const BlockVector &V,
+ const bool fast=false);
/**
* Set all entries to zero. Equivalent to
@@ -95,6 +117,41 @@ class BlockVector
*/
void clear ();
+ /**
+ * Swap the contents of this
+ * vector and the other vector
+ * #v#. One could do this
+ * operation with a temporary
+ * variable and copying over the
+ * data elements, but this
+ * function is significantly more
+ * efficient since it only swaps
+ * the pointers to the data of
+ * the two vectors and therefore
+ * does not need to allocate
+ * temporary storage and move
+ * data around.
+ *
+ * This function is analog to the
+ * the #swap# function of all C++
+ * standard containers. Also,
+ * there is a global function
+ * #swap(u,v)# that simply calls
+ * #u.swap(v)#, again in analogy
+ * to standard functions.
+ */
+ void swap (BlockVector &v);
+
+ /**
+ * Access to a single block.
+ */
+ Vector& block(unsigned int i);
+
+ /**
+ * Read-only access to a single block.
+ */
+ const Vector& block(unsigned int i) const;
+
/**
* $U(0-N) = s$: fill all components.
*/
@@ -150,38 +207,6 @@ class BlockVector
Number linfty_norm () const;
- /**
- * Change the dimension of the vector to
- * #N#. The reserved memory for this vector
- * remains unchanged if possible, to make
- * things faster, but this may waste some
- * memory, so take this in the back of your
- * head.
- * However, if #N==0# all memory is freed,
- * i.e. if you want to resize the vector
- * and release the memory not needed, you
- * have to first call #reinit(0)# and then
- * #reinit(N)#. This cited behaviour is
- * analogous to that of the STL containers.
- *
- * On #fast==false#, the vector is filled by
- * zeros.
- */
- void reinit (const vector& N,
- const bool fast=false);
-
- /**
- * Change the dimension to that of the
- * vector #V#. The same applies as for
- * the other #reinit# function.
- *
- * The elements of #V# are not copied, i.e.
- * this function is the same as calling
- * #reinit (V.size(), fast)#.
- */
- void reinit (const BlockVector &V,
- const bool fast=false);
-
/**
* Return dimension of the vector. This is the
* sum of the dimensions of all components.
@@ -525,4 +550,21 @@ BlockVector::block(unsigned int i) const
}
+/**
+ * Global function #swap# which overloads the default implementation
+ * of the C++ standard library which uses a temporary object. The
+ * function simply exchanges the data of the two vectors.
+ *
+ * @author Wolfgang Bangerth, 2000
+ */
+template
+inline
+void swap (BlockVector &u, BlockVector &v)
+{
+ u.swap (v);
+};
+
+
+
+
#endif
diff --git a/deal.II/lac/include/lac/block_vector.templates.h b/deal.II/lac/include/lac/block_vector.templates.h
index 6d54d753bf..e7c262d854 100644
--- a/deal.II/lac/include/lac/block_vector.templates.h
+++ b/deal.II/lac/include/lac/block_vector.templates.h
@@ -104,12 +104,27 @@ void BlockVector::reinit (const BlockVector& v
}
+
template
BlockVector::~BlockVector ()
{
}
+
+template
+void BlockVector::swap (BlockVector &v)
+{
+ for (unsigned int i=0; i
void BlockVector::clear ()
{
diff --git a/deal.II/lac/include/lac/solver_minres.h b/deal.II/lac/include/lac/solver_minres.h
index 6bf5e26a35..bff486943c 100644
--- a/deal.II/lac/include/lac/solver_minres.h
+++ b/deal.II/lac/include/lac/solver_minres.h
@@ -208,9 +208,9 @@ SolverMinRes::solve (const MATRIX &A,
v.reinit(VS,true);
// some values needed
- vector delta(3);
- vector f(2);
- vector e(2);
+ double delta[3];
+ double f[2];
+ double e[2];
double r_l2 = 0;
double r0 = 0;
@@ -222,7 +222,7 @@ SolverMinRes::solve (const MATRIX &A,
double d = 0;
// The iteration step.
- int j = 1;
+ unsigned int j = 1;
// Start of the solution process
@@ -314,11 +314,24 @@ SolverMinRes::solve (const MATRIX &A,
// All vectors have to be shifted
// one iteration step.
// This should be changed one time.
- m[2] = m[1];
- m[1] = m[0];
-
- u[0] = u[1];
- u[1] = u[2];
+ //
+ // the previous code was like this:
+ // m[2] = m[1];
+ // m[1] = m[0];
+ // but it can be made more efficient,
+ // since the value of m[0] is no more
+ // needed in the next iteration
+ swap (m[2], m[1]);
+ swap (m[1], m[0]);
+
+ // likewise, but reverse direction:
+ // u[0] = u[1];
+ // u[1] = u[2];
+ swap (u[0], u[1]);
+ swap (u[1], u[2]);
+
+ // these are scalars, so need
+ // to bother
f[0] = f[1];
e[0] = e[1];
delta[0] = delta[1];
diff --git a/deal.II/lac/include/lac/vector.h b/deal.II/lac/include/lac/vector.h
index 6e39eaa79d..76d7c140b0 100644
--- a/deal.II/lac/include/lac/vector.h
+++ b/deal.II/lac/include/lac/vector.h
@@ -104,7 +104,64 @@ class Vector {
* the size of the vector, unlike the
* STL's #vector<>::clear# function.
*/
- void clear ();
+ void clear ();
+
+ /**
+ * Change the dimension of the vector to
+ * #N#. The reserved memory for this vector
+ * remains unchanged if possible, to make
+ * things faster, but this may waste some
+ * memory, so take this in the back of your
+ * head.
+ * However, if #N==0# all memory is freed,
+ * i.e. if you want to resize the vector
+ * and release the memory not needed, you
+ * have to first call #reinit(0)# and then
+ * #reinit(N)#. This cited behaviour is
+ * analogous to that of the STL containers.
+ *
+ * On #fast==false#, the vector is filled by
+ * zeros.
+ */
+ void reinit (const unsigned int N,
+ const bool fast=false);
+
+ /**
+ * Change the dimension to that of the
+ * vector #V#. The same applies as for
+ * the other #reinit# function.
+ *
+ * The elements of #V# are not copied, i.e.
+ * this function is the same as calling
+ * #reinit (V.size(), fast)#.
+ */
+ void reinit (const Vector &V,
+ const bool fast=false);
+
+ /**
+ * Swap the contents of this
+ * vector and the other vector
+ * #v#. One could do this
+ * operation with a temporary
+ * variable and copying over the
+ * data elements, but this
+ * function is significantly more
+ * efficient since it only swaps
+ * the pointers to the data of
+ * the two vectors and therefore
+ * does not need to allocate
+ * temporary storage and move
+ * data around.
+ *
+ * This function is analog to the
+ * the #swap# function of all C++
+ * standard containers. Also,
+ * there is a global function
+ * #swap(u,v)# that simply calls
+ * #u.swap(v)#, again in analogy
+ * to standard functions.
+ */
+ void swap (Vector &v);
/**
* $U(0-N) = s$: fill all components.
@@ -158,39 +215,6 @@ class Vector {
*/
Number linfty_norm () const;
-
- /**
- * Change the dimension of the vector to
- * #N#. The reserved memory for this vector
- * remains unchanged if possible, to make
- * things faster, but this may waste some
- * memory, so take this in the back of your
- * head.
- * However, if #N==0# all memory is freed,
- * i.e. if you want to resize the vector
- * and release the memory not needed, you
- * have to first call #reinit(0)# and then
- * #reinit(N)#. This cited behaviour is
- * analogous to that of the STL containers.
- *
- * On #fast==false#, the vector is filled by
- * zeros.
- */
- void reinit (const unsigned int N,
- const bool fast=false);
-
- /**
- * Change the dimension to that of the
- * vector #V#. The same applies as for
- * the other #reinit# function.
- *
- * The elements of #V# are not copied, i.e.
- * this function is the same as calling
- * #reinit (V.size(), fast)#.
- */
- void reinit (const Vector &V,
- const bool fast=false);
-
/**
* Return dimension of the vector. This
* function was formerly called #n()#, but
@@ -204,7 +228,7 @@ class Vector {
* Return whether the vector contains only
* elements with value zero. This function
* is mainly for internal consistency
- * check and should seldomly be used when
+ * checks and should seldomly be used when
* not in debug mode since it uses quite
* some time.
*/
@@ -524,4 +548,21 @@ Number& Vector::operator() (const unsigned int i)
}
+
+/**
+ * Global function #swap# which overloads the default implementation
+ * of the C++ standard library which uses a temporary object. The
+ * function simply exchanges the data of the two vectors.
+ *
+ * @author Wolfgang Bangerth, 2000
+ */
+template
+inline
+void swap (Vector &u, Vector &v)
+{
+ u.swap (v);
+};
+
+
+
#endif
diff --git a/deal.II/lac/include/lac/vector.templates.h b/deal.II/lac/include/lac/vector.templates.h
index 5d79eb4257..c303a580ba 100644
--- a/deal.II/lac/include/lac/vector.templates.h
+++ b/deal.II/lac/include/lac/vector.templates.h
@@ -19,6 +19,7 @@
#include
+
template
static inline Number sqr (const Number x) {
return x*x;
@@ -76,6 +77,14 @@ Vector::Vector (const Vector& v) :
// }
+template
+Vector::~Vector ()
+{
+ if (val) delete[] val;
+}
+
+
+
template
void Vector::reinit (const unsigned int n, const bool fast) {
if (n==0)
@@ -100,27 +109,35 @@ void Vector::reinit (const unsigned int n, const bool fast) {
template
-void Vector::reinit (const Vector& v, const bool fast) {
+void Vector::reinit (const Vector& v, const bool fast)
+{
reinit (v.size(), fast);
};
template
-Vector::~Vector ()
+void Vector::clear ()
{
- if (val) delete[] val;
+ if (dim>0)
+ fill (begin(), end(), 0.);
}
+
template
-void Vector::clear () {
- if (dim>0)
- fill (begin(), end(), 0.);
-}
+void
+Vector::swap (Vector &v)
+{
+ std::swap (dim, v.dim);
+ std::swap (maxdim, v.maxdim);
+ std::swap (val, v.val);
+};
+
template
-bool Vector::all_zero () const {
+bool Vector::all_zero () const
+{
Assert (dim!=0, ExcEmptyVector());
const_iterator p = begin(),
--
2.39.5