]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add the possibility to create block sparse matrices without an Epetra_Map when runnin...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 21 Sep 2008 04:20:06 +0000 (04:20 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 21 Sep 2008 04:20:06 +0000 (04:20 +0000)
git-svn-id: https://svn.dealii.org/trunk@16875 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/trilinos_block_sparse_matrix.h
deal.II/lac/source/trilinos_block_sparse_matrix.cc

index 2167a2617b686c39cc313b6c3eb1984fd41b576d..4b9d3a761a62a52e15ea4120cd03602cbca006e5 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2004, 2005, 2006, 2007 by the deal.II authors
+//    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
@@ -182,6 +182,15 @@ namespace TrilinosWrappers
       void reinit (const std::vector<Epetra_Map> &input_maps,
                   const BlockSparsityPattern    &block_sparsity_pattern);
 
+                                       /**
+                                        * Resize the matrix and initialize it
+                                        * by the given sparsity pattern. Since
+                                        * no distribution map is given, the
+                                        * result is a block matrix for which
+                                        * all elements are stored locally.
+                                        */
+      void reinit (const BlockSparsityPattern &block_sparsity_pattern);
+      
                                        /**
                                         * This function initializes the
                                        * Trilinos matrix using the deal.II
index 1b8d171490305d584a7730cff8377fc82377d291..8d36581430373219bea95ea32429888123c77413 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id: trilinos_block_sparse_matrix.cc 15631 2008-01-17 23:47:31Z bangerth $
 //    Version: $Name$
 //
-//    Copyright (C) 2004, 2005, 2006, 2008 by the deal.II authors
+//    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
@@ -93,6 +93,13 @@ namespace TrilinosWrappers
   reinit (const std::vector<Epetra_Map> &input_maps,
          const BlockSparsityPattern    &block_sparsity_pattern)
   {
+    Assert (input_maps.size() == block_sparsity_pattern.n_block_rows(),
+           ExcDimensionMismatch (input_maps.size(),
+                                 block_sparsity_pattern.n_block_rows()));
+    Assert (input_maps.size() == block_sparsity_pattern.n_block_cols(),
+           ExcDimensionMismatch (input_maps.size(),
+                                 block_sparsity_pattern.n_block_cols()));
+    
     const unsigned int n_block_rows = input_maps.size();
 
     if (input_maps[0].Comm().MyPID()==0)
@@ -122,6 +129,38 @@ namespace TrilinosWrappers
 
 
 
+  void
+  BlockSparseMatrix::
+  reinit (const BlockSparsityPattern    &block_sparsity_pattern)
+  {
+    Assert (block_sparsity_pattern.n_block_rows() ==
+           block_sparsity_pattern.n_block_cols(),
+           ExcDimensionMismatch (block_sparsity_pattern.n_block_rows(),
+                                 block_sparsity_pattern.n_block_cols()));
+    Assert (block_sparsity_pattern.n_rows() ==
+           block_sparsity_pattern.n_cols(),
+           ExcDimensionMismatch (block_sparsity_pattern.n_rows(),
+                                 block_sparsity_pattern.n_cols()));
+    
+                                    // produce a dummy local map and pass it
+                                    // off to the other function
+#ifdef DEAL_II_COMPILER_SUPPORTS_MPI
+    Epetra_MpiComm    trilinos_communicator (MPI_COMM_WORLD);
+#else
+    Epetra_SerialComm trilinos_communicator;
+#endif
+
+    std::vector<Epetra_Map> input_maps;
+    for (unsigned int i=0; i<block_sparsity_pattern.n_block_rows(); ++i)
+      input_maps.push_back (Epetra_Map(block_sparsity_pattern.block(i,0).n_rows(),
+                                      0,
+                                      trilinos_communicator));
+
+    reinit (input_maps, block_sparsity_pattern);
+  }
+  
+
+
   void
   BlockSparseMatrix::
   reinit (const std::vector<Epetra_Map>             &input_maps,

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.