]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix some compiling problems. Fix sparse matrix iterator.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 23 Nov 2009 02:35:34 +0000 (02:35 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 23 Nov 2009 02:35:34 +0000 (02:35 +0000)
git-svn-id: https://svn.dealii.org/trunk@20147 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/trilinos_sparse_matrix.h
deal.II/lac/include/lac/trilinos_vector.h
deal.II/lac/source/trilinos_sparse_matrix.cc

index e40e8caacee66f0d9c1ba46dc5327c3165c9e69e..2b64416e3a816219084ed0eb72817842a9b1915e 100644 (file)
@@ -167,13 +167,13 @@ namespace TrilinosWrappers
                                        * than one accessor can access
                                        * this data if necessary.
                                        */
-            std_cxx1x::shared_ptr<const std::vector<unsigned int> > colnum_cache;
+            std_cxx1x::shared_ptr<std::vector<unsigned int> > colnum_cache;
 
                                        /**
                                        * Similar cache for the values
                                        * of this row.
                                        */
-            std_cxx1x::shared_ptr<const std::vector<TrilinosScalar> > value_cache;
+            std_cxx1x::shared_ptr<std::vector<TrilinosScalar> > value_cache;
 
                                       /**
                                        * Discard the old row caches
@@ -3235,7 +3235,7 @@ namespace TrilinosWrappers
 
     temp_vector.reinit(dst, true);
 
-    vmult (temp_vector, src);
+    Tvmult (temp_vector, src);
     dst += temp_vector;
   }
 
index 5bff7053ebfc0d5fb3940e3c4b62ee840a85c9f9..ec0a4e82e79602ac46e82ed65757192bf5ec75b4 100644 (file)
@@ -367,6 +367,15 @@ namespace TrilinosWrappers
        explicit Vector (const Epetra_Map &parallel_partitioning,
                         const VectorBase &v);
 
+                                      /**
+                                       * Reinitialize from a deal.II
+                                       * vector. The Epetra_Map specifies the
+                                       * %parallel partitioning.
+                                       */
+       template <typename number>
+       void reinit (const Epetra_Map             &parallel_partitioner,
+                    const dealii::Vector<number> &v);
+
                                       /**
                                        * Reinit functionality. This
                                        * function destroys the old
@@ -469,19 +478,7 @@ namespace TrilinosWrappers
     Vector::Vector (const Epetra_Map             &input_map,
                     const dealii::Vector<number> &v)
     {
-      vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(input_map));
-
-      const int min_my_id = input_map.MinMyGID();
-      const int size = input_map.NumMyElements();
-
-      Assert (input_map.MaxLID() == size-1,
-             ExcDimensionMismatch(input_map.MaxLID(), size-1));
-
-                                  // Need to copy out values, since the
-                                  // deal.II might not use doubles, so
-                                  // that a direct access is not possible.
-      for (int i=0; i<size; ++i)
-       (*vector)[0][i] = v(i);
+      reinit (input_map, v);
     }
 
 
@@ -497,6 +494,25 @@ namespace TrilinosWrappers
 
 
 
+
+    template <typename number>
+    void Vector::reinit (const Epetra_Map             &parallel_partitioner,
+                        const dealii::Vector<number> &v)
+    {
+      if (&*vector != 0 && vector->Map().SameAs(parallel_partitioner))
+       vector = std::auto_ptr<Epetra_FEVector>
+         (new Epetra_FEVector(parallel_partitioner));
+
+      const int size = parallel_partitioner.NumMyElements();
+
+                                  // Need to copy out values, since the
+                                  // deal.II might not use doubles, so
+                                  // that a direct access is not possible.
+      for (int i=0; i<size; ++i)
+       (*vector)[0][i] = v(parallel_partitioner.GID(i));
+    }
+
+
     inline
     Vector &
     Vector::operator = (const TrilinosScalar s)
@@ -513,7 +529,7 @@ namespace TrilinosWrappers
     {
       if (size() != v.size())
        {
-         *vector = std::auto_ptr<Epetra_FEVector>
+         vector = std::auto_ptr<Epetra_FEVector>
          (new Epetra_FEVector(Epetra_Map (v.size(), 0,
 #ifdef DEAL_II_COMPILER_SUPPORTS_MPI
                                           Epetra_MpiComm(MPI_COMM_SELF)
@@ -763,8 +779,7 @@ namespace TrilinosWrappers
        vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map));
       }
 
-    Epetra_Map & map = vector_partitioner();
-    const int min_my_id = map.MinMyGID();
+    const Epetra_Map & map = vector_partitioner();
     const int size = map.NumMyElements();
 
     Assert (map.MaxLID() == size-1,
index 94de3a0f70833e8e0389fe852fb21a4f8568b9e5..363c8a1788084cff7277e80b6a86df97d2dc5ccd 100644 (file)
@@ -57,12 +57,24 @@ namespace TrilinosWrappers
                                  // row
       int ncols;
       int colnums = matrix->n();
-      TrilinosScalar *values = new TrilinosScalar(colnums);
+      if (value_cache.get() == 0)
+       {
+         value_cache.reset (new std::vector<TrilinosScalar> (matrix->n()));
+         colnum_cache.reset (new std::vector<unsigned int> (matrix->n()));
+       }
+      else
+       {
+         value_cache->resize (matrix->n());
+         colnum_cache->resize (matrix->n());
+       }
 
-      int ierr;
-      ierr = matrix->trilinos_matrix().ExtractGlobalRowCopy((int)this->a_row,
-                                                           colnums,
-                                                           ncols, &(values[0]));
+      int ierr = matrix->trilinos_matrix().
+       ExtractGlobalRowCopy((int)this->a_row,
+                            colnums,
+                            ncols, &((*value_cache)[0]),
+                            reinterpret_cast<int*>(&((*colnum_cache)[0])));
+      value_cache->resize (ncols);
+      colnum_cache->resize (ncols);
       AssertThrow (ierr == 0, ExcTrilinosError(ierr));
 
                                  // copy it into our caches if the
@@ -71,10 +83,6 @@ namespace TrilinosWrappers
                                  // we shouldn't have initialized an
                                  // iterator for an empty line (what
                                  // would it point to?)
-      Assert (ncols != 0, ExcInternalError());
-      colnum_cache.reset (new std::vector<unsigned int> (colnums,
-                                                        colnums+ncols));
-      value_cache.reset (new std::vector<TrilinosScalar> (values, values+ncols));
     }
   }
 
@@ -1397,7 +1405,7 @@ namespace TrilinosWrappers
          {
            matrix->ExtractMyRowView (i, num_entries, values, indices);
            for (int j=0; j<num_entries; ++j)
-             out << "(" << i << "," << indices[matrix->GRID(j)] << ") "
+             out << "(" << matrix->GRID(i) << "," << indices[matrix->GCID(j)] << ") "
                  << values[j] << std::endl;
          }
       }

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.