]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
The template functions that copied a deal.II vector to a Trilinos vector were there...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 23 Oct 2008 03:06:09 +0000 (03:06 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 23 Oct 2008 03:06:09 +0000 (03:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@17309 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/trilinos_vector.h
deal.II/lac/include/lac/trilinos_vector_base.h
deal.II/lac/source/trilinos_vector.cc
deal.II/lac/source/trilinos_vector_base.cc
deal.II/lac/source/trilinos_vector_base.inst.in [new file with mode: 0644]

index 1386d52fa0e748129d5d9aab9250c03886b9f1b6..12308c0e09871e1810628ea7d68376bcfaa994c8 100755 (executable)
@@ -326,16 +326,10 @@ namespace TrilinosWrappers
                                        * change the map, use the
                                        * reinit(const Epetra_Map
                                        * &input_map) function.
-                                       *
-                                       * Since Trilinos only works on
-                                       * doubles, this function is
-                                       * limited to accept only one
-                                       * possible number type in the
-                                       * deal.II vector.
                                        */
        template <typename Number>
        Vector & 
-         operator = (const ::dealii::Vector<Number> &v);
+       operator = (const ::dealii::Vector<Number> &v);
 
                                       /**
                                        * This reinit function is
@@ -414,7 +408,7 @@ namespace TrilinosWrappers
     {
       u.swap (v);
     }
-
+  
 
 #ifndef DOXYGEN
 
@@ -453,6 +447,16 @@ namespace TrilinosWrappers
       return *this;
     }
 
+
+    template <typename Number>
+    Vector & 
+    Vector::operator = (const ::dealii::Vector<Number> &v)
+    {
+      VectorBase::operator = (v);
+      return *this;
+    }
+    
+    
 #endif
 
   } /* end of namespace MPI */
@@ -594,8 +598,11 @@ namespace TrilinosWrappers
        operator = (const Vector &V); 
 
     private:
+                                      /**
+                                       * A map indicating the size of the
+                                       * vector.
+                                       */
       Epetra_LocalMap map;
-
   };
 
 
@@ -612,46 +619,54 @@ namespace TrilinosWrappers
  * @relates TrilinosWrappers::Vector
  * @author Martin Kronbichler, Wolfgang Bangerth, 2008
  */
-    inline
-    void swap (Vector &u, Vector &v)
-    {
-      u.swap (v);
-    }
+  inline
+  void swap (Vector &u, Vector &v)
+  {
+    u.swap (v);
+  }
 
 
 #ifndef DOXYGEN
 
-    template <typename number>
-    Vector::Vector (const dealii::Vector<number> &v)
-                    :
+  template <typename number>
+  Vector::Vector (const dealii::Vector<number> &v)
+                 :
 #ifdef DEAL_II_COMPILER_SUPPORTS_MPI
-                    map (v.size(), 0, Epetra_MpiComm(MPI_COMM_WORLD))
+                 map (v.size(), 0, Epetra_MpiComm(MPI_COMM_WORLD))
 #else
-                   map (v.size(), 0, Epetra_SerialComm())
+                 map (v.size(), 0, Epetra_SerialComm())
 #endif
-    {
-      vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map, false));
+  {
+    vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map, false));
 
-      std::vector<int> indices (v.size());
-      for (unsigned int i=0; i<v.size(); ++i)
-       indices[i] = map.GID(i);
+    std::vector<int> indices (v.size());
+    for (unsigned int i=0; i<v.size(); ++i)
+      indices[i] = map.GID(i);
 
-      const int ierr = vector->ReplaceGlobalValues(v.size(), 0, v.begin(), 
-                                                  &indices[0]);
-      AssertThrow (ierr == 0, ExcTrilinosError());
-    }
+    const int ierr = vector->ReplaceGlobalValues(v.size(), 0, v.begin(), 
+                                                &indices[0]);
+    AssertThrow (ierr == 0, ExcTrilinosError());
+  }
 
   
   
-    inline
-    Vector &
-    Vector::operator = (const TrilinosScalar s)
-    {
-      VectorBase::operator = (s);
+  inline
+  Vector &
+  Vector::operator = (const TrilinosScalar s)
+  {
+    VectorBase::operator = (s);
 
-      return *this;
-    }
+    return *this;
+  }
 
+  template <typename Number>
+  Vector & 
+  Vector::operator = (const ::dealii::Vector<Number> &v)
+  {
+    VectorBase::operator = (v);
+    return *this;
+  }
+  
 #endif
 
 
index 7bffc8adb8e7cdbff7339060c805eeb817fa966a..0448e3c63d9970db226958d4e09998101bc9b1d6 100644 (file)
@@ -364,6 +364,27 @@ namespace TrilinosWrappers
       VectorBase &
        operator = (const TrilinosScalar s);
 
+                                       /**
+                                       * Another copy function. This
+                                       * one takes a deal.II vector and
+                                       * copies it into a
+                                       * TrilinosWrapper vector. Note
+                                       * that since we do not provide
+                                       * any Epetra_map that tells
+                                       * about the partitioning of the
+                                       * vector among the MPI
+                                       * processes, the size of the
+                                       * TrilinosWrapper vector has to
+                                       * be the same as the size of the
+                                       * input vector. In order to
+                                       * change the map, use the
+                                       * reinit(const Epetra_Map
+                                       * &input_map) function.
+                                       */
+      template <typename Number>
+      VectorBase & 
+      operator = (const ::dealii::Vector<Number> &v);
+
                                        /**
                                         * Test for equality. This
                                         * function assumes that the
index caa5d8d56f2cde15925d82ae62f96ba96c14ee77..24eed161befbc597e3d291dd42fc1bc67c93d5a8 100755 (executable)
@@ -191,35 +191,6 @@ namespace TrilinosWrappers
 
 
 
-    template <typename Number>
-    Vector &
-    Vector::operator = (const ::dealii::Vector<Number> &v)
-    {
-      Assert (size() == v.size(),
-             ExcDimensionMismatch(size(), v.size()));
-
-      vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map));
-      
-      const int min_my_id = map.MinMyGID();
-      const unsigned int size = map.NumMyElements();
-
-      Assert (map.MaxLID() == size-1,
-             ExcDimensionMismatch(map.MaxLID(), size-1));
-
-      std::vector<int> indices (size);
-      for (unsigned int i=0; i<size; ++i)
-       indices[i] = map.GID(i);
-
-
-      const int ierr = vector->ReplaceGlobalValues(size, 0, v.begin()+min_my_id, 
-                                                  &indices[0]);
-      AssertThrow (ierr == 0, ExcTrilinosError());
-
-      return *this;
-    }
-
-
-
     void
     Vector::do_data_exchange (const TrilinosWrappers::SparseMatrix &m,
                              const Vector                         &v)
@@ -424,30 +395,7 @@ namespace TrilinosWrappers
 
     return *this;
   }
-
-
-
-  template <typename Number>
-  Vector &
-  Vector::operator = (const ::dealii::Vector<Number> &v)
-  {
-    if (size() != v.size())
-      {
-       map = LocalMap (v.size(), 0, vector->Comm());
-       vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map));
-      }
-
-    std::vector<int> indices (v.size());
-    for (unsigned int i=0; i<v.size(); ++i)
-      indices[i] = map.GID(i);
-
-    const int ierr = vector->ReplaceGlobalValues(v.size(), 0, v.begin(), 
-                                                &indices[0]);
-    AssertThrow (ierr == 0, ExcTrilinosError());
-
-    return *this;
-  }
-
+  
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 2522478eab81da42b795d06a661ef4ae4ad323c9..7d74519947904812f9cd24ce5d674cc42c967487 100644 (file)
@@ -146,6 +146,40 @@ namespace TrilinosWrappers
   }
 
 
+  template <typename number>
+  VectorBase &
+  VectorBase::operator = (const ::dealii::Vector<number> &v)
+  {
+    Assert (size() == v.size(),
+           ExcDimensionMismatch(size(), v.size()));
+    
+                                    // this is probably not very efficient
+                                    // but works. in particular, we could do
+                                    // better if we know that
+                                    // number==TrilinosScalar because then we
+                                    // could elide the copying of elements
+                                    //
+                                    // let's hope this isn't a
+                                    // particularly frequent operation
+    std::pair<unsigned int, unsigned int>
+      local_range = this->local_range ();
+    std::vector<unsigned int>   indices (local_range.second -
+                                        local_range.first);
+    std::vector<TrilinosScalar> values (local_range.second -
+                                       local_range.first);
+
+    for (unsigned int i=0; i<local_range.second-local_range.first; ++i)
+      {
+       indices[i] = i + local_range.first;
+       values[i] = v(i + local_range.first);
+      }
+    
+    set (indices.size(), &indices[0], &values[0]);
+      
+    return *this;
+  }
+  
+
 
   bool
   VectorBase::operator == (const VectorBase &v) const
@@ -949,6 +983,11 @@ namespace TrilinosWrappers
 } /* end of namespace TrilinosWrappers */
 
 
+namespace TrilinosWrappers
+{
+#include "trilinos_vector_base.inst"
+}
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif // DEAL_II_USE_TRILINOS
diff --git a/deal.II/lac/source/trilinos_vector_base.inst.in b/deal.II/lac/source/trilinos_vector_base.inst.in
new file mode 100644 (file)
index 0000000..c8c4b12
--- /dev/null
@@ -0,0 +1,20 @@
+//---------------------------------------------------------------------------
+//    $Id: vector.cc 15443 2007-11-05 03:43:52Z bangerth $
+//    Version: $Name$
+//
+//    Copyright (C) 2008 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.
+//
+//---------------------------------------------------------------------------
+
+
+for (SCALAR : REAL_SCALARS)
+  {
+    template
+      VectorBase &
+      VectorBase::operator = (const ::dealii::Vector<SCALAR> &v);
+  }

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.