]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
implement MGMatrixBlockVector
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 25 Sep 2010 20:29:59 +0000 (20:29 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 25 Sep 2010 20:29:59 +0000 (20:29 +0000)
git-svn-id: https://svn.dealii.org/trunk@22157 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 7b75ed41547023499c94f4127a6c5cc324415d01..6cebc0c8014308e5475c6316174f7f6aefafc44e 100644 (file)
@@ -91,10 +91,20 @@ class MGLevelObject : public Subscriptor
                                     /**
                                      * Coarsest level for multigrid.
                                      */
+    unsigned int min_level () const;
+    
+                                    /**
+                                     * Finest level for multigrid.
+                                     */
+    unsigned int max_level () const;
+    
+                                    /**
+                                     * @deprecated Replaced by min_level()
+                                     */
     unsigned int get_minlevel () const;
     
                                     /**
-                                     * Ignored
+                                     * @deprecated Replaced by max_level()
                                      */
     unsigned int get_maxlevel () const;
     
@@ -204,6 +214,22 @@ MGLevelObject<Object>::get_maxlevel () const
 }
 
 
+template<class Object>
+unsigned int
+MGLevelObject<Object>::min_level () const
+{
+  return minlevel;
+}
+
+
+template<class Object>
+unsigned int
+MGLevelObject<Object>::max_level () const
+{
+  return minlevel + objects.size() - 1;
+}
+
+
 template<class Object>
 unsigned int
 MGLevelObject<Object>::memory_consumption () const
index 00597ecd4486749ec5b73fd9ed21ea1b1f802b18..e44ddc07eb16418f36fae0a128240aaba750967a 100644 (file)
@@ -18,6 +18,7 @@
 #include <base/smartpointer.h>
 #include <base/std_cxx1x/shared_ptr.h>
 #include <base/memory_consumption.h>
+#include <base/mg_level_object.h>
 #include <lac/block_indices.h>
 #include <lac/block_sparsity_pattern.h>
 #include <lac/sparse_matrix.h>
@@ -408,7 +409,7 @@ class MatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<MAT
  * @author Baerbel Janssen, Guido Kanschat, 2010
  */
 template <class MATRIX>
-class MGMatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<MATRIX> > >
+class MGMatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MGLevelObject<MatrixBlock<MATRIX> > > >
 {
   public:
                                     /**
@@ -426,18 +427,26 @@ class MGMatrixBlockVector : public NamedData<std_cxx1x::shared_ptr<MatrixBlock<M
     void add(unsigned int row, unsigned int column,
             const std::string& name,
             const BlockIndices* block_indices = 0);
-    
+
+                                    /**
+                                     * For matrices using a
+                                     * SparsityPattern, this function
+                                     * reinitializes each matrix in
+                                     * the vector with the correct
+                                     * pattern from the block system.
+                                     */
+    void reinit(const MGLevelObject<BlockSparsityPattern>& sparsity);
                                     /**
                                      * Access a constant reference to
                                      * the block at position <i>i</i>.
                                      */
-    const MatrixBlock<MATRIX>& block(unsigned int i) const;
+    const MGLevelObject<MatrixBlock<MATRIX> >& block(unsigned int i) const;
     
                                     /**
                                      * Access a reference to
                                      * the block at position <i>i</i>.
                                      */
-    MatrixBlock<MATRIX>& block(unsigned int i);
+    MGLevelObject<MatrixBlock<MATRIX> >& block(unsigned int i);
     
                                     /**
                                      * The memory used by this object.
@@ -461,10 +470,11 @@ namespace internal
 {
   template <class MATRIX>
   void
-  reinit(MatrixBlockVector<MATRIX>& v, const BlockSparsityPattern& p)
+  reinit(MatrixBlockVector<MATRIX>&, const BlockSparsityPattern&)
   {
     Assert(false, ExcNotImplemented());
   }
+
   
   template <typename number>
   void
@@ -473,6 +483,25 @@ namespace internal
     for (unsigned int i=0;i<v.size();++i)
       v(i)->matrix.reinit(p.block(v(i)->row, v(i)->column));
   }
+
+
+  template <class MATRIX>
+  void
+  reinit(MGMatrixBlockVector<MATRIX>&, const MGLevelObject<BlockSparsityPattern>&)
+  {
+    Assert(false, ExcNotImplemented());
+  }
+  
+  
+  template <typename number>
+  void
+  reinit(MGMatrixBlockVector<dealii::SparseMatrix<number> >& v,
+        const MGLevelObject<BlockSparsityPattern>& sparsity)
+  {
+    for (unsigned int i=0;i<v.size();++i)
+      for (unsigned int level=sparsity.min_level(); level <= sparsity.max_level();++level)
+       v(i)->matrix[level].reinit(sparsity[level].block(v(i)->row, v(i)->column));
+  }
   
 }
 
@@ -719,6 +748,55 @@ MatrixBlockVector<MATRIX>::matrix(unsigned int i)
 
 
 
+//----------------------------------------------------------------------//
+
+template <class MATRIX>
+inline void
+MGMatrixBlockVector<MATRIX>::add(
+  unsigned int row, unsigned int column,
+  const std::string& name,
+  const BlockIndices* block_indices)
+{
+  std_cxx1x::shared_ptr<MGLevelObject<MatrixBlock<MATRIX> > >
+    p(new MGLevelObject<MatrixBlock<MATRIX> >(row, column, block_indices));
+  NamedData<std_cxx1x::shared_ptr<MatrixBlock<MATRIX> > >::add(p, name);
+}
+
+
+template <class MATRIX>
+inline void
+MGMatrixBlockVector<MATRIX>::reinit(const MGLevelObject<BlockSparsityPattern>& sparsity)
+{
+  internal::reinit(*this, sparsity);
+}
+
+
+
+template <class MATRIX>
+inline const MGLevelObject<MatrixBlock<MATRIX> >&
+MGMatrixBlockVector<MATRIX>::block(unsigned int i) const
+{
+  return *this->read(i);
+}
+
+
+template <class MATRIX>
+inline MGLevelObject<MatrixBlock<MATRIX> >&
+MGMatrixBlockVector<MATRIX>::block(unsigned int i)
+{
+  return *(*this)(i);
+}
+
+
+template <class MATRIX>
+inline MATRIX&
+MGMatrixBlockVector<MATRIX>::matrix(unsigned int i)
+{
+  return (*this)(i)->matrix;
+}
+
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif

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.