]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Added MatrixTools::apply_boundary_values for PETScWrappers::MPI::BlockSparseMatrix...
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 27 Oct 2008 15:39:00 +0000 (15:39 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 27 Oct 2008 15:39:00 +0000 (15:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@17351 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/matrices.h
deal.II/deal.II/source/numerics/matrices.all_dimensions.cc
deal.II/doc/news/changes.h

index 540695fcd065688d390cf8e248631246ec639ad7..6f540a497f1dc599cea3ed0f0906efb869f516d8 100644 (file)
@@ -55,7 +55,9 @@ namespace PETScWrappers
   namespace MPI
   {
     class SparseMatrix;
+    class BlockSparseMatrix;
     class Vector;
+    class BlockVector;
   }
 }
 #endif
@@ -1109,6 +1111,17 @@ class MatrixTools : public MatrixCreator
                           PETScWrappers::Vector       &solution,
                           PETScWrappers::MPI::Vector  &right_hand_side,
                           const bool             eliminate_columns = true);
+    
+                                    /**
+                                      * Same as above but for BlockSparseMatrix.
+                                     */                                        
+    static void
+    apply_boundary_values (const std::map<unsigned int,double>  &boundary_values,
+                       PETScWrappers::MPI::BlockSparseMatrix  &matrix,
+                       PETScWrappers::MPI::BlockVector        &solution,
+                       PETScWrappers::MPI::BlockVector        &right_hand_side,
+                       const bool       eliminate_columns = true);
+    
 #endif
 
 #ifdef DEAL_II_USE_TRILINOS
index 4366991f8e4b3f3c073c4a746130b5921124c420..910c8147f3f793a4b3e40a7e691d29520e4266f4 100644 (file)
@@ -24,6 +24,7 @@
 #  include <lac/petsc_sparse_matrix.h>
 #  include <lac/petsc_parallel_vector.h>
 #  include <lac/petsc_vector.h>
+#  include <lac/petsc_parallel_block_sparse_matrix.h>
 #endif
 
 #ifdef DEAL_II_USE_TRILINOS
@@ -681,6 +682,81 @@ apply_boundary_values (const std::map<unsigned int,double> &boundary_values,
 }
 
 
+void
+MatrixTools::
+apply_boundary_values (const std::map<unsigned int,double>  &boundary_values,
+                       PETScWrappers::MPI::BlockSparseMatrix  &matrix,
+                       PETScWrappers::MPI::BlockVector        &solution,
+                       PETScWrappers::MPI::BlockVector        &right_hand_side,
+                       const bool                            eliminate_columns)
+{
+  Assert (matrix.n() == right_hand_side.size(),
+          ExcDimensionMismatch(matrix.n(), right_hand_side.size()));
+  Assert (matrix.n() == solution.size(),
+          ExcDimensionMismatch(matrix.n(), solution.size()));
+  Assert (matrix.n_block_rows() == matrix.n_block_cols(),
+          ExcNotQuadratic());
+
+  const unsigned int n_blocks = matrix.n_block_rows();
+
+  matrix.compress();
+
+                                   // We need to find the subdivision
+                                   // into blocks for the boundary values.
+                                   // To this end, generate a vector of
+                                   // maps with the respective indices.
+  std::vector<std::map<unsigned int,double> > block_boundary_values(n_blocks);
+  {
+    int offset = 0, block=0;
+    for (std::map<unsigned int,double>::const_iterator
+            dof  = boundary_values.begin();
+          dof != boundary_values.end();
+          ++dof)
+      {
+        if (dof->first >= matrix.block(block,0).m() + offset)
+          {
+            offset += matrix.block(block,0).m();
+            block++;
+          }
+        const unsigned int index = dof->first - offset;
+        block_boundary_values[block].insert(std::pair<unsigned int, double> (index,dof->second));
+      }
+  }
+
+                                   // Now call the non-block variants on
+                                   // the diagonal subblocks and the
+                                   // solution/rhs.
+  for (unsigned int block=0; block<n_blocks; ++block)
+    PETScWrappers::apply_boundary_values(block_boundary_values[block],
+                                         matrix.block(block,block),
+                                         solution.block(block),
+                                         right_hand_side.block(block),
+                                         eliminate_columns);
+
+                                  // Finally, we need to do something
+                                   // about the off-diagonal matrices. This
+                                   // is luckily not difficult. Just clear
+                                   // the whole row.
+  for (unsigned int block_m=0; block_m<n_blocks; ++block_m)
+    {
+      const std::pair<unsigned int, unsigned int> local_range
+        = matrix.block(block_m,0).local_range();
+
+      std::vector<unsigned int> constrained_rows;
+      for (std::map<unsigned int,double>::const_iterator
+            dof  = block_boundary_values[block_m].begin();
+          dof != block_boundary_values[block_m].end();
+          ++dof)
+        if ((dof->first >= local_range.first) &&
+            (dof->first < local_range.second))
+          constrained_rows.push_back (dof->first);
+
+      for (unsigned int block_n=0; block_n<n_blocks; ++block_n)
+        if (block_m != block_n)
+          matrix.block(block_m,block_n).clear_rows(constrained_rows);
+    }
+}
+
 #endif
 
 
index 6e67d2d13768d85a0c331e2d15dcc97913dc9fde..353ff07428f160353f45a6c6f27ccd2a242e88d8 100644 (file)
@@ -319,6 +319,16 @@ inconvenience this causes.
   (WB 2008/09/23)
   </p>
 
+  <li>
+  <p>
+  New: Added the CompressedSimpleSparsityPattern as an alternative to 
+  CompressedSparsityPattern and CompressedSetSparsityPattern, which
+  should be faster in most cases but its memory usage is somewhere
+  inbetween.
+  <br>
+  (Timo Heister 2008/09/03)
+  </p>
+
   <li>
   <p>
   Improved: The CompressedSparsityPattern can now elide some operations
@@ -353,6 +363,13 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: MatrixTools::apply_boundary_values() also for PETScWrappers::MPI::BlockSparseMatrix.
+  <br>
+  (Timo Heister 2008/10/27)
+  </p>
+
   <li>
   <p>
   New: When calling function DoFTools::make_sparsity_pattern with a ConstraintMatrix, it is now possible to set a bool argument keep_constrained_dofs. When this flag is set to false, constrained rows and columns will not be part of the sparsity pattern, which increases the performance of matrix operations and decrease memory consumption in case there are many constraints.

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.