]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add reinit function
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 6 Jul 2010 17:55:58 +0000 (17:55 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 6 Jul 2010 17:55:58 +0000 (17:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@21457 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/matrix_block.h

index c3752db5b299166e03fe160bf1edaf666d6fa9ac..ce414d21e9161e2f6251182403fb48798347fc47 100644 (file)
@@ -18,6 +18,8 @@
 #include <base/smartpointer.h>
 #include <base/std_cxx1x/shared_ptr.h>
 #include <lac/block_indices.h>
+#include <lac/block_sparsity_pattern.h>
+#include <lac/sparse_matrix.h>
 #include <lac/full_matrix.h>
 
 
@@ -352,16 +354,28 @@ class MatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<MAT
     void add(unsigned int row, unsigned int column,
             const std::string& name,
             const BlockIndices* block_indices);
+
+                                    /**
+                                     * For matrices using a
+                                     * SparsityPattern, this function
+                                     * reinitializes each matrix in
+                                     * the vector with the correct
+                                     * pattern from the block system.
+                                     */
+    void reinit(const BlockSparsityPattern& sparsity);
+    
                                     /**
                                      * Access a constant reference to
                                      * the block at position <i>i</i>.
                                      */
     const MatrixBlock<MATRIX>& block(unsigned int i) const;
+    
                                     /**
                                      * Access a reference to
                                      * the block at position <i>i</i>.
                                      */
     MatrixBlock<MATRIX>& block(unsigned int i);
+    
                                     /**
                                      * Access the matrix at position
                                      * <i>i</i> for read and write
@@ -389,21 +403,31 @@ class MGMatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<M
                                     /**
                                      * Add a new matrix block at the
                                      * position <tt>(row,column)</tt>
-                                     * in the block system.
+                                     * in the block system. The third
+                                     * argument allows to give the
+                                     * matrix a name for later
+                                     * identification.
+                                     *
+                                     * @deprecated The final
+                                     * argument is ignored and will
+                                     * be removed in a future release.
                                      */
     void add(unsigned int row, unsigned int column,
             const std::string& name,
-            const BlockIndices* block_indices);
+            const BlockIndices* block_indices = 0);
+    
                                     /**
                                      * Access a constant reference to
                                      * the block at position <i>i</i>.
                                      */
     const MatrixBlock<MATRIX>& block(unsigned int i) const;
+    
                                     /**
                                      * Access a reference to
                                      * the block at position <i>i</i>.
                                      */
     MatrixBlock<MATRIX>& block(unsigned int i);
+    
                                     /**
                                      * Access the matrix at position
                                      * <i>i</i> for read and write
@@ -417,6 +441,25 @@ class MGMatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<M
 
 //----------------------------------------------------------------------//
 
+namespace internal
+{
+  template <class MATRIX>
+  void
+  reinit(MatrixBlockVector<MATRIX>& v, const BlockSparsityPattern& p)
+  {
+    Assert(false, ExcNotImplemented());
+  }
+  
+  template <typename number>
+  void
+  reinit(MatrixBlockVector<dealii::SparseMatrix<number> >& v, const BlockSparsityPattern& p)
+  {
+    for (unsigned int i=0;i<v.size();++i)
+      v(i)->matrix.reinit(p.block(v(i)->row, v(i)->column));
+  }
+  
+}
+
 template <class MATRIX>
 inline
 MatrixBlock<MATRIX>::MatrixBlock()
@@ -616,6 +659,15 @@ MatrixBlockVector<MATRIX>::add(
 }
 
 
+template <class MATRIX>
+inline void
+MatrixBlockVector<MATRIX>::reinit(const BlockSparsityPattern& sparsity)
+{
+  internal::reinit(*this, sparsity);
+}
+
+
+
 template <class MATRIX>
 inline const MatrixBlock<MATRIX>&
 MatrixBlockVector<MATRIX>::block(unsigned int i) const

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.