</p>
<li> <p>
- New: <code class="class">Solver*</code>: There is a virtual function
+ New: <code class="class">SolverXX</code>: There is a virtual function
<code>print_vectors</code> called in every step. It is void in the
solver itself but can be used to print intermediate iteration
vectors.
</p>
+
+ <li> <p>
+ New: there are now functions <code class="member">Vector::swap</code>
+ and <code class="member">BlockVector::swap</code>, as well as
+ global functions <code class="member">swap(u,v)</code> 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.
+ </p>
</ol>
~BlockVector ();
/**
- * Access to a single block.
- */
- Vector<Number>& 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<unsigned int>& 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<Number>& block(unsigned int i) const;
+ void reinit (const BlockVector<n_blocks,Number> &V,
+ const bool fast=false);
/**
* Set all entries to zero. Equivalent to
*/
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<n_blocks,Number> &v);
+
+ /**
+ * Access to a single block.
+ */
+ Vector<Number>& block(unsigned int i);
+
+ /**
+ * Read-only access to a single block.
+ */
+ const Vector<Number>& block(unsigned int i) const;
+
/**
* $U(0-N) = s$: fill all components.
*/
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<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 BlockVector<n_blocks,Number> &V,
- const bool fast=false);
-
/**
* Return dimension of the vector. This is the
* sum of the dimensions of all components.
}
+/**
+ * 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 <int n_blocks, typename Number>
+inline
+void swap (BlockVector<n_blocks,Number> &u, BlockVector<n_blocks,Number> &v)
+{
+ u.swap (v);
+};
+
+
+
+
#endif
}
+
template <int n_blocks, typename Number>
BlockVector<n_blocks,Number>::~BlockVector ()
{
}
+
+template <int n_blocks, typename Number>
+void BlockVector<n_blocks,Number>::swap (BlockVector<n_blocks,Number> &v)
+{
+ for (unsigned int i=0; i<n_blocks; ++i)
+ {
+ swap (components[i], v.components[i]);
+ swap (start[i], v.start[i]);
+ };
+ swap (start[n_blocks], v.start[n_blocks]);
+};
+
+
+
template <int n_blocks, typename Number>
void BlockVector<n_blocks,Number>::clear ()
{
v.reinit(VS,true);
// some values needed
- vector<double> delta(3);
- vector<double> f(2);
- vector<double> e(2);
+ double delta[3];
+ double f[2];
+ double e[2];
double r_l2 = 0;
double r0 = 0;
double d = 0;
// The iteration step.
- int j = 1;
+ unsigned int j = 1;
// Start of the solution process
// 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];
* 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<Number> &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<Number> &v);
/**
* $U(0-N) = s$: fill all components.
*/
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<Number> &V,
- const bool fast=false);
-
/**
* Return dimension of the vector. This
* function was formerly called #n()#, but
* 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.
*/
}
+
+/**
+ * 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 <typename Number>
+inline
+void swap (Vector<Number> &u, Vector<Number> &v)
+{
+ u.swap (v);
+};
+
+
+
#endif
#include <algorithm>
+
template <typename Number>
static inline Number sqr (const Number x) {
return x*x;
// }
+template <typename Number>
+Vector<Number>::~Vector ()
+{
+ if (val) delete[] val;
+}
+
+
+
template <typename Number>
void Vector<Number>::reinit (const unsigned int n, const bool fast) {
if (n==0)
template <typename Number>
-void Vector<Number>::reinit (const Vector<Number>& v, const bool fast) {
+void Vector<Number>::reinit (const Vector<Number>& v, const bool fast)
+{
reinit (v.size(), fast);
};
template <typename Number>
-Vector<Number>::~Vector ()
+void Vector<Number>::clear ()
{
- if (val) delete[] val;
+ if (dim>0)
+ fill (begin(), end(), 0.);
}
+
template <typename Number>
-void Vector<Number>::clear () {
- if (dim>0)
- fill (begin(), end(), 0.);
-}
+void
+Vector<Number>::swap (Vector<Number> &v)
+{
+ std::swap (dim, v.dim);
+ std::swap (maxdim, v.maxdim);
+ std::swap (val, v.val);
+};
+
template <typename Number>
-bool Vector<Number>::all_zero () const {
+bool Vector<Number>::all_zero () const
+{
Assert (dim!=0, ExcEmptyVector());
const_iterator p = begin(),