]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Allow assignment from block vectors to vectors and back.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 31 Mar 2006 15:19:12 +0000 (15:19 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 31 Mar 2006 15:19:12 +0000 (15:19 +0000)
git-svn-id: https://svn.dealii.org/trunk@12800 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.html
deal.II/lac/include/lac/block_vector.h
deal.II/lac/include/lac/block_vector_base.h
deal.II/lac/include/lac/vector.h
deal.II/lac/include/lac/vector.templates.h
tests/lac/Makefile
tests/lac/block_vector_vector_assign.cc [new file with mode: 0644]
tests/lac/block_vector_vector_assign/.cvsignore [new file with mode: 0644]
tests/lac/block_vector_vector_assign/cmp/generic [new file with mode: 0644]

index 471ad5ee67866a11099fcc3c9cdc9b330b4b4c15..5fe39c8aef6adb3d3161ae00f61dc1030a490293 100644 (file)
@@ -434,13 +434,20 @@ inconvenience this causes.
 
 <ol>
   
+  <li> <p> New: There are now assignment operators from <code
+       class="class">BlockVector</code> to <code
+       class="class">Vector</code> and back.
+       <br>
+       (WB 2006/03/30)
+       </p>
+
   <li> <p> Improved: <code class="class">BlockSparsityPattern</code> can be
-  initialized directly using the vector generated by
-  <code class="class">DoFTools</code>::<code
-  class="member">compute_row_length_vector</code>.
-  <br>
-  (GK 2006/03/30)
-  </p>
+       initialized directly using the vector generated by
+       <code class="class">DoFTools</code>::<code
+       class="member">compute_row_length_vector</code>.
+       <br>
+       (GK 2006/03/30)
+       </p>
 
   <li> <p>
        Improved: All matrix (and some vector) classes now check whether
index 495e1273378fbae8811f9afd679709f7ae70fbc1..c90d7675ba8140e461d0513f6c8c817b8321cf62 100644 (file)
@@ -187,7 +187,14 @@ class BlockVector : public BlockVectorBase<Vector<Number> >
     BlockVector &
     operator= (const BlockVector<Number2> &V);
 
-                                    /**
+                                    /**
+                                     * Copy a regular vector into a
+                                     * block vector.
+                                     */
+    BlockVector &
+    operator= (const Vector<Number> &V);
+
+                                     /**
                                      * Reinitialize the BlockVector to
                                      * contain <tt>num_blocks</tt> blocks of
                                      * size <tt>block_size</tt> each.
@@ -443,6 +450,17 @@ BlockVector<Number>::operator = (const BlockVector &v)
 
 
 
+template <typename Number>
+inline
+BlockVector<Number> &
+BlockVector<Number>::operator = (const Vector<Number> &v)
+{
+  BaseClass::operator = (v);
+  return *this;
+}
+
+
+
 template <typename Number>
 template <typename Number2>
 inline
index eb97e4bad41427d5a13b1caa5cd79c37f2cf616c..7e46c5c3231d655a2502d46b238ea017151fa4d6 100644 (file)
@@ -806,6 +806,13 @@ class BlockVectorBase
     BlockVectorBase &
     operator= (const BlockVectorBase<VectorType2> &V);
 
+                                     /**
+                                      * Copy operator from non-block
+                                      * vectors to block vectors.
+                                      */
+    BlockVectorBase &
+    operator = (const VectorType &v);
+    
                                      /**
                                       * Check for equality of two block vector
                                       * types. This operation is only allowed
@@ -2067,6 +2074,7 @@ void BlockVectorBase<VectorType>::equ (const value_type    a,
 }
 
 
+
 template <class VectorType>
 BlockVectorBase<VectorType>&
 BlockVectorBase<VectorType>::operator = (const value_type s)
@@ -2090,9 +2098,8 @@ BlockVectorBase<VectorType>::operator = (const BlockVectorBase<VectorType>& v)
          ExcDimensionMismatch(n_blocks(), v.n_blocks()));
   
   for (unsigned int i=0;i<n_blocks();++i)
-    {
-      components[i] = v.components[i];
-    }
+    components[i] = v.components[i];
+
   return *this;
 }
 
@@ -2113,6 +2120,23 @@ BlockVectorBase<VectorType>::operator = (const BlockVectorBase<VectorType2> &v)
 
 
 
+template <class VectorType>
+BlockVectorBase<VectorType>&
+BlockVectorBase<VectorType>::operator = (const VectorType &v)
+{
+  Assert (size() == v.size(),
+         ExcDimensionMismatch(size(), v.size()));
+
+  unsigned int index_v = 0;
+  for (unsigned int b=0;b<n_blocks();++b)
+    for (unsigned int i=0; i<block(b).size(); ++i, ++index_v)
+      block(b)(i) = v(index_v);
+
+  return *this;
+}
+
+
+
 template <class VectorType>
 template <class VectorType2>
 inline
index 4f2fcf79127b06add3b3115110811e5e57aebf0c..c02041b365ced579b404d939a26b51ea25111d02 100644 (file)
@@ -34,6 +34,10 @@ namespace PETScWrappers
 
 template<typename number> class LAPACKFullMatrix;
 
+template <typename> class BlockVector;
+
+
+
 /*! @addtogroup Vectors
  *@{
  */
@@ -307,6 +311,13 @@ class Vector :
     template <typename Number2>
     Vector<Number> & operator= (const Vector<Number2> &v);
 
+                                     /**
+                                      * Copy operator for assigning a
+                                      * block vector to a regular
+                                      * vector.
+                                      */
+    Vector<Number> & operator= (const BlockVector<Number> &v);
+    
 #ifdef DEAL_II_USE_PETSC
                                      /**
                                       * Another copy operator: copy the values
index dc188728aee2ebf109c2842bd6f70ec9e315b1c7..887c02bafb045390ea8012d8e0a82ef36cb001d7 100644 (file)
@@ -15,6 +15,7 @@
 
 
 #include <lac/vector.h>
+#include <lac/block_vector.h>
 
 #ifdef DEAL_II_USE_PETSC
 #  include <lac/petsc_vector.h>
@@ -672,7 +673,7 @@ void Vector<Number>::ratio (const Vector<Number> &a, const Vector<Number> &b)
 
 
 template <typename Number>
-Vector<Number>&
+Vector<Number> &
 Vector<Number>::operator = (const Vector<Number>& v)
 {
   if (v.vec_size != vec_size)
@@ -687,7 +688,7 @@ Vector<Number>::operator = (const Vector<Number>& v)
 
 template <typename Number>
 template <typename Number2>
-Vector<Number>&
+Vector<Number> &
 Vector<Number>::operator = (const Vector<Number2>& v)
 {
   if (v.size() != vec_size)
@@ -700,6 +701,23 @@ Vector<Number>::operator = (const Vector<Number2>& v)
 
 
 
+template <typename Number>
+Vector<Number> &
+Vector<Number>::operator = (const BlockVector<Number>& v)
+{
+  if (v.size() != vec_size)
+    reinit (v.size(), true);
+
+  unsigned int this_index = 0;
+  for (unsigned int b=0; b<v.n_blocks(); ++b)
+    for (unsigned int i=0; i<v.block(b).size(); ++i, ++this_index)
+      val[this_index] = v.block(b)(i);
+  
+  return *this;
+}
+
+
+
 #ifdef DEAL_II_USE_PETSC
 
 template <typename Number>
index beaca1e4f3ea62394e59df268c3aa64e8bd78cb0..91fe46527d90a39857486a3f2c821a8f40359556 100644 (file)
@@ -24,10 +24,11 @@ tests_x = vector-vector \
        full_matrix householder tridiagonal_matrix \
        sparsity_pattern sparse_matrices matrices \
        block_sparsity_pattern_01 \
-       block_vector block_vector_iterator block_matrices \
+       block_vector block_vector_* block_matrices \
        matrix_lib matrix_out \
        solver eigen \
-       sparse_ilu sparse_ilu_t lapack lapack_fill block_minres
+       sparse_ilu sparse_ilu_t lapack lapack_fill block_minres \
+       amg_*
 
 
 # from above list of regular expressions, generate the real set of
diff --git a/tests/lac/block_vector_vector_assign.cc b/tests/lac/block_vector_vector_assign.cc
new file mode 100644 (file)
index 0000000..032de9f
--- /dev/null
@@ -0,0 +1,110 @@
+//--------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2006 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//--------------------------------------------------------------------
+
+// check assignment between block vectors and regular vectors
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <lac/block_vector.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+#include <algorithm>
+#include <numeric>
+#include <utility>
+#include <cmath>
+
+template <typename Vector1, typename Vector2>
+bool operator == (const Vector1 &v1,
+                 const Vector2 &v2)
+{
+  if (v1.size() != v2.size())
+    return false;
+  for (unsigned int i=0; i<v1.size(); ++i)
+    if (v1(i) != v2(i))
+      return false;
+  return true;
+}
+
+
+void test ()
+{
+  std::vector<unsigned int> ivector(4);
+  ivector[0] = 2;
+  ivector[1] = 4;
+  ivector[2] = 3;
+  ivector[3] = 5;
+  
+  BlockVector<double> v1(ivector);
+  Vector<double> v2(v1.size());
+
+  for (unsigned int i=0; i<v1.size(); ++i)
+    v1(i) = 1+i*i;
+
+  v2 = v1;
+  Assert (v1==v2, ExcInternalError());
+
+  BlockVector<double> v3 (ivector);
+  v3 = v2;
+  Assert (v3==v2, ExcInternalError());
+  Assert (v3==v1, ExcInternalError());
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+
+
+int main ()
+{
+  std::ofstream logfile("block_vector_vector_assign/output");
+  logfile.setf(std::ios::fixed);
+  logfile.precision(3);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+  
+  try
+    {
+      test ();
+    }
+  catch (std::exception &e)
+    {
+      std::cerr << std::endl << std::endl
+          << "----------------------------------------------------"
+          << std::endl;
+      std::cerr << "Exception on processing: " << e.what() << std::endl
+          << "Aborting!" << std::endl
+          << "----------------------------------------------------"
+          << std::endl;
+                                      // abort
+      return 2;
+    }
+  catch (...) 
+    {
+      std::cerr << std::endl << std::endl
+          << "----------------------------------------------------"
+          << std::endl;
+      std::cerr << "Unknown exception!" << std::endl
+          << "Aborting!" << std::endl
+          << "----------------------------------------------------"
+          << std::endl;
+                                      // abort
+      return 3;
+    };
+  
+  
+  return 0;
+}
+
diff --git a/tests/lac/block_vector_vector_assign/.cvsignore b/tests/lac/block_vector_vector_assign/.cvsignore
new file mode 100644 (file)
index 0000000..d196dfb
--- /dev/null
@@ -0,0 +1 @@
+obj.* exe OK output
diff --git a/tests/lac/block_vector_vector_assign/cmp/generic b/tests/lac/block_vector_vector_assign/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.