]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Updated a few functions in Trilinos sparse matrix and vector classes to make more...
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 7 Nov 2008 12:31:41 +0000 (12:31 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 7 Nov 2008 12:31:41 +0000 (12:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@17511 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/trilinos_block_vector.h
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 daee99f3a3c4c5a2bd3356e68e21f5e139708221..be98f073dade4b81daa2068eb1775e65f55cde8f 100644 (file)
@@ -504,7 +504,7 @@ namespace TrilinosWrappers
                                        * elements.
                                        */
       BlockVector (const BlockVector  &V);
-    
+
                                        /**
                                        * Creates a block vector
                                        * consisting of
@@ -516,7 +516,40 @@ namespace TrilinosWrappers
                                        * reinit of the blocks.
                                        */
       BlockVector (const unsigned int num_blocks);
-    
+
+                                       /**
+                                       * Constructor. Set the number of
+                                       * blocks to <tt>n.size()</tt> and
+                                       * initialize each block with
+                                       * <tt>n[i]</tt> zero elements.
+                                       *
+                                       * References BlockVector.reinit().
+                                       */
+      BlockVector (const std::vector<unsigned int> &N);
+
+                                       /**
+                                        * Constructor. Set the number of
+                                        * blocks to
+                                        * <tt>n.size()</tt>. Initialize the
+                                        * vector with the elements
+                                        * pointed to by the range of
+                                        * iterators given as second and
+                                        * third argument. Apart from the
+                                        * first argument, this
+                                        * constructor is in complete
+                                        * analogy to the respective
+                                        * constructor of the
+                                        * <tt>std::vector</tt> class, but the
+                                        * first argument is needed in
+                                        * order to know how to subdivide
+                                        * the block vector into
+                                        * different blocks.
+                                        */
+      template <typename InputIterator>
+      BlockVector (const std::vector<unsigned int> &n,
+                   const InputIterator              first,
+                   const InputIterator              end);
+
                                        /**
                                        * Destructor. Clears memory
                                        */
@@ -732,6 +765,37 @@ namespace TrilinosWrappers
 
 
 
+  inline
+  BlockVector::BlockVector (const std::vector<unsigned int> &N)
+  {
+    reinit (N);
+  }
+
+
+
+  template <typename InputIterator>
+  BlockVector::BlockVector (const std::vector<unsigned int> &n,
+                            const InputIterator              first,
+                            const InputIterator              end)
+  {
+                                     // first set sizes of blocks, but
+                                     // don't initialize them as we will
+                                     // copy elements soon
+    reinit (n, true);
+    InputIterator start = first;
+    for (unsigned int b=0; b<n.size(); ++b)
+      {
+        InputIterator end = start;
+        std::advance (end, static_cast<signed int>(n[b]));
+
+       for (unsigned int i=0; i<n[b]; ++i, ++start)
+         this->block(b)(i) = *start;
+      }
+    Assert (start == end, ExcIteratorRangeDoesNotMatchVectorSize());
+  }
+
+
+
   inline
   BlockVector::BlockVector (const unsigned int num_blocks)
   {
index b557e10e7666234142a08c9066e694e339b8054a..98129b4a006eebe6bb1850db41ee314647893c9c 100755 (executable)
@@ -416,30 +416,26 @@ namespace TrilinosWrappers
                    const unsigned int  n_max_entries_per_row);
 
                                        /**
-                                        * This constructor is similar to
-                                       * the one above, but it now
-                                       * takes two different Epetra
-                                       * maps for rows and
-                                       * columns. This interface is
-                                       * meant to be used for
-                                       * generating rectangular
-                                       * matrices, where one map
-                                       * specifies the parallel
+                                        * This constructor is similar to the
+                                       * one above, but it now takes two
+                                       * different Epetra maps for rows and
+                                       * columns. This interface is meant
+                                       * to be used for generating
+                                       * rectangular matrices, where one
+                                       * map specifies the parallel
                                        * distribution of rows and the
-                                       * second one specifies the
-                                       * number of columns in the
-                                       * total matrix. It also provides
-                                       * information for the internal
-                                       * arrangement in matrix vector
-                                       * products, but is not used for
-                                       * the distribution of the
-                                       * columns &ndash; rather, all
-                                       * column elements of a row are
-                                       * stored on the same
-                                       * processor. The vector
-                                       * n_entries_per_row specifies
-                                       * the number of entries in each
-                                       * row of the newly generated
+                                       * second one specifies the number of
+                                       * columns in the total matrix. It
+                                       * also provides information for the
+                                       * internal arrangement in matrix
+                                       * vector products, but is not used
+                                       * for the distribution of the
+                                       * columns &ndash; rather, all column
+                                       * elements of a row are stored on
+                                       * the same processor. The vector
+                                       * <tt>n_entries_per_row</tt>
+                                       * specifies the number of entries in
+                                       * each row of the newly generated
                                        * matrix.
                                         */
       SparseMatrix (const Epetra_Map                &InputRowMap,
@@ -452,15 +448,29 @@ namespace TrilinosWrappers
                                         * #n columns. The resulting matrix
                                         * will be completely stored locally.
                                        *
-                                       * The number of columns entries
-                                       * per row is specified as the
-                                       * maximum number of entries
-                                       * argument.
+                                       * The number of columns entries per
+                                       * row is specified as the maximum
+                                       * number of entries argument.
                                         */
       SparseMatrix (const unsigned int  m,
                    const unsigned int  n,
                    const unsigned int  n_max_entries_per_row);
 
+                                       /**
+                                        * Generate a matrix that is completely
+                                        * stored locally, having #m rows and
+                                        * #n columns. The resulting matrix
+                                        * will be completely stored locally.
+                                       *
+                                       * The vector
+                                       * <tt>n_entries_per_row</tt>
+                                       * specifies the number of entries in
+                                       * each row.
+                                        */
+      SparseMatrix (const unsigned int               m,
+                   const unsigned int               n,
+                   const std::vector<unsigned int> &n_entries_per_row);
+
                                        /**
                                         * Copy constructor. Sets the
                                         * calling matrix to be the same
index 2f03830f7ee3ddcabca95a8270fdb8238e4bed7d..a0babddce0a4c72229833a730f4dff2fb75dace0 100755 (executable)
@@ -418,21 +418,28 @@ namespace TrilinosWrappers
                     :
                     map (InputMap)
     {
-      vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map, false));
+      vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map));
       
       const int min_my_id = map.MinMyGID();
-      const unsigned int size = map.NumMyElements();
+      const int size = map.NumMyElements();
 
       Assert (map.MaxLID() == size-1,
              ExcDimensionMismatch(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.
       std::vector<int> indices (size);
+      std::vector<double> values (size);
       for (unsigned int i=0; i<size; ++i)
-       indices[i] = map.GID(i);
+       {
+         indices[i] = map.GID(i);
+         values[i] = v(i);
+       }
 
+      const int ierr = vector->ReplaceGlobalValues (size, &indices[0], 
+                                                   &values[0]);
 
-      const int ierr = vector->ReplaceGlobalValues(size, v.begin()+min_my_id, 
-                                                  &indices[0]);
       AssertThrow (ierr == 0, ExcTrilinosError());
     }
 
@@ -637,15 +644,20 @@ namespace TrilinosWrappers
                  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));
 
     std::vector<int> indices (v.size());
+    std::vector<double> values (v.size());
     for (unsigned int i=0; i<v.size(); ++i)
-      indices[i] = map.GID(i);
+      {
+       indices[i] = map.GID(i);
+       values[i] = v(i);
+      }
+
+    const int ierr = vector->ReplaceGlobalValues((int)v.size(),
+                                                &indices[0], &values[0]);
 
-    const int ierr = vector->ReplaceGlobalValues(v.size(), v.begin(), 
-                                                &indices[0]);
-    AssertThrow (ierr == 0, ExcTrilinosError());
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
   }
 
   
index dfc1be2f26b0ba0021f4596e3e0dc3ddfccd20a3..ea603df8b9b75906fdc349330e7f0433d70fddbd 100755 (executable)
@@ -167,6 +167,25 @@ namespace TrilinosWrappers
                                        int(n_max_entries_per_row), false)))
   {}
 
+  SparseMatrix::SparseMatrix (const unsigned int               m,
+                             const unsigned int               n,
+                             const std::vector<unsigned int> &n_entries_per_row)
+                 :
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+                  row_map (m, 0, Epetra_MpiComm(MPI_COMM_WORLD)),
+                  col_map (n, 0, Epetra_MpiComm(MPI_COMM_WORLD)),
+#else
+                  row_map (m, 0, Epetra_SerialComm()),
+                  col_map (n, 0, Epetra_SerialComm()),
+#endif
+                 last_action (Zero),
+                 compressed (true),
+                 matrix (std::auto_ptr<Epetra_FECrsMatrix>
+                    (new Epetra_FECrsMatrix(Copy, row_map, 
+                       (int*)const_cast<unsigned int*>(&(n_entries_per_row[0])), 
+                                            false)))
+  {}
+
   SparseMatrix::SparseMatrix (const SparseMatrix &InputMatrix)
                  :
                   Subscriptor(),
@@ -609,8 +628,11 @@ namespace TrilinosWrappers
 
     Assert (numbers::is_finite(value),
            ExcMessage("The given value is not finite but either "
-                     "infinite or Not A Number (NaN)"));
+                      "infinite or Not A Number (NaN)"));
 
+                                       // If the last thing we did was an
+                                       // addition, we need to do all the
+                                       // data exchange.
     if (last_action == Add)
       {
        int ierr;
@@ -619,6 +641,9 @@ namespace TrilinosWrappers
        AssertThrow (ierr == 0, ExcTrilinosError(ierr));
       }
 
+                                       // If our last action was not an
+                                       // addition, we tell that at this
+                                       // point here.
     if (last_action != Insert)
       last_action = Insert;
 
@@ -638,13 +663,19 @@ namespace TrilinosWrappers
                                        // command. Otherwise, we're just
                                        // able to replace existing entries.
     if (matrix->Filled() == false)
-      ierr = matrix->InsertGlobalValues (trilinos_i, 1,
-                                        const_cast<double*>(&value), 
-                                        &trilinos_j);
+      {
+       ierr = matrix->InsertGlobalValues (trilinos_i, 1,
+                                          const_cast<double*>(&value), 
+                                          &trilinos_j);
+       if (ierr > 0)
+         ierr = 0;
+      }
     else
-      ierr = matrix->ReplaceGlobalValues (trilinos_i, 1,
-                                         const_cast<double*>(&value), 
-                                         &trilinos_j);
+      {
+       ierr = matrix->ReplaceGlobalValues (trilinos_i, 1,
+                                           const_cast<double*>(&value), 
+                                           &trilinos_j);
+      }
 
     AssertThrow (ierr <= 0, ExcAccessToNonPresentElement(i,j));
     AssertThrow (ierr == 0, ExcTrilinosError(ierr));

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.