]> https://gitweb.dealii.org/ - dealii.git/commitdiff
introduce MatrixBlockVector
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 18 Mar 2010 04:46:44 +0000 (04:46 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 18 Mar 2010 04:46:44 +0000 (04:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@20842 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/mesh_worker_assembler.h
deal.II/deal.II/include/numerics/mesh_worker_info.h
deal.II/deal.II/include/numerics/mesh_worker_info.templates.h
deal.II/deal.II/source/numerics/mesh_worker_info.cc
deal.II/lac/include/lac/matrix_block.h

index 36524292e162e147d72d2536c1a7981076d1879c..015bfca1e90cf55ede17fc4d9778ace7c33da9ba 100644 (file)
@@ -662,9 +662,6 @@ namespace MeshWorker
     class MatrixLocalBlocksToGlobalBlocks
     {
       public:
-                                        /// The object that is stored
-       typedef boost::shared_ptr<MatrixBlock<MATRIX> > MatrixPtr;
-       
                                         /**
                                          * Constructor, initializing
                                          * the #threshold, which
@@ -681,7 +678,7 @@ namespace MeshWorker
                                          * cell matrix vectors.
                                          */
       void initialize(const BlockInfo* block_info,
-                     std::vector<MatrixPtr>& matrices);
+                     MatrixBlockVector<MATRIX>& matrices);
                                         /**
                                          * Initialize the local data
                                          * in the
@@ -737,7 +734,8 @@ namespace MeshWorker
                                          * stored as a vector of
                                          * pointers.
                                          */
-       std::vector<MatrixPtr> matrices;
+       SmartPointer<MatrixBlockVector<MATRIX>,
+                    MatrixLocalBlocksToGlobalBlocks<MATRIX, number> > matrices;
       
       /**
        * A pointer to the object containing the block structure.
@@ -1400,10 +1398,12 @@ namespace MeshWorker
     
     template <class MATRIX, typename number>
     inline void
-    MatrixLocalBlocksToGlobalBlocks<MATRIX, number>::initialize(const BlockInfo* b, std::vector<MatrixPtr>& m)
+    MatrixLocalBlocksToGlobalBlocks<MATRIX, number>::initialize(
+      const BlockInfo* b,
+      MatrixBlockVector<MATRIX>& m)
     {
       block_info = b;
-      matrices = m;
+      matrices = &m;
     }
     
 
@@ -1415,7 +1415,7 @@ namespace MeshWorker
       DoFInfo<dim>& info,
       bool interior_face) const
     {
-      info.initialize_matrices(matrices, interior_face);
+      info.initialize_matrices(*matrices, interior_face);
     }
 
 
@@ -1467,14 +1467,14 @@ namespace MeshWorker
     MatrixLocalBlocksToGlobalBlocks<MATRIX, number>::assemble(
       const DoFInfo<dim>& info)
     {
-      for (unsigned int i=0;i<matrices.size();++i)
+      for (unsigned int i=0;i<matrices->size();++i)
        {
                                           // Row and column index of
                                           // the block we are dealing with
-         const unsigned int row = matrices[i]->row;
-         const unsigned int col = matrices[i]->column;
+         const unsigned int row = matrices->block(i).row;
+         const unsigned int col = matrices->block(i).column;
 
-         assemble(matrices[i]->matrix, info.matrix(i,false).matrix, row, col, info.indices, info.indices);
+         assemble(matrices->matrix(i), info.matrix(i,false).matrix, row, col, info.indices, info.indices);
        }
     }
 
@@ -1486,17 +1486,17 @@ namespace MeshWorker
       const DoFInfo<dim>& info1,
       const DoFInfo<dim>& info2)
     {
-      for (unsigned int i=0;i<matrices.size();++i)
+      for (unsigned int i=0;i<matrices->size();++i)
        {
                                           // Row and column index of
                                           // the block we are dealing with
-         const unsigned int row = matrices[i]->row;
-         const unsigned int col = matrices[i]->column;
+         const unsigned int row = matrices->block(i).row;
+         const unsigned int col = matrices->block(i).column;
 
-         assemble(matrices[i]->matrix, info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices);
-         assemble(matrices[i]->matrix, info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices);
-         assemble(matrices[i]->matrix, info2.matrix(i,false).matrix, row, col, info2.indices, info2.indices);
-         assemble(matrices[i]->matrix, info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices);
+         assemble(matrices->matrix(i), info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices);
+         assemble(matrices->matrix(i), info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices);
+         assemble(matrices->matrix(i), info2.matrix(i,false).matrix, row, col, info2.indices, info2.indices);
+         assemble(matrices->matrix(i), info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices);
        }
     }
 
index 5e69d90ee39f1affe84bead7cba04ee2c66514a3..9b13df27d8ca0d6481319976da5909dfbf22b7d6 100644 (file)
@@ -84,8 +84,8 @@ namespace MeshWorker
                                        * instantiation in order to
                                        * provide row and column info.
                                        */
-      template <class MatrixPtr>
-      void initialize_matrices(const std::vector<MatrixPtr>& matrices,
+      template <class MATRIX>
+      void initialize_matrices(const MatrixBlockVector<MATRIX>& matrices,
                               bool both);
 
                                       /**
@@ -672,7 +672,8 @@ namespace MeshWorker
                                        * #values, #gradients and
                                        * #hessians.
                                        */
-      void fill_local_data(const DoFInfo<dim, spacedim>& info, const bool split_fevalues);
+      template<typename number>
+      void fill_local_data(const DoFInfo<dim, spacedim, number>& info, bool split_fevalues);
 
                                       /**
                                        * The global data vector
@@ -782,10 +783,10 @@ namespace MeshWorker
 
 
   template <typename number>
-  template <class MatrixPtr>
+  template <class MATRIX>
   inline void
   LocalResults<number>::initialize_matrices(
-    const std::vector<MatrixPtr>& matrices,
+    const MatrixBlockVector<MATRIX>& matrices,
     bool both)
   {
     M1.resize(matrices.size());
@@ -793,8 +794,8 @@ namespace MeshWorker
       M2.resize(matrices.size());
     for (unsigned int i=0;i<matrices.size();++i)
       {
-       const unsigned int row = matrices[i]->row;
-       const unsigned int col = matrices[i]->column;
+       const unsigned int row = matrices.block(i).row;
+       const unsigned int col = matrices.block(i).column;
 
        M1[i].row = row;
        M1[i].column = col;
@@ -1142,6 +1143,35 @@ namespace MeshWorker
 
 //----------------------------------------------------------------------//
 
+  template<int dim, int sdim>
+  template <class FEVALUES>
+  inline void
+  IntegrationInfo<dim,sdim>::initialize(
+    const FiniteElement<dim,sdim>& el,
+    const Mapping<dim,sdim>& mapping,
+    const Quadrature<FEVALUES::integral_dimension>& quadrature,
+    const UpdateFlags flags,
+    const BlockInfo* block_info)
+  {
+    if (block_info == 0 || block_info->local().size() == 0)
+      {
+       fevalv.resize(1);             
+       fevalv[0] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
+         new FEVALUES (mapping, el, quadrature, flags));
+      }
+    else
+      {
+       fevalv.resize(el.n_base_elements());
+       for (unsigned int i=0;i<fevalv.size();++i)
+         {
+           fevalv[i] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
+             new FEVALUES (mapping, el.base_element(i), quadrature, flags));
+         }
+      }
+    n_components = el.n_components();
+  }
+  
+
   template <int dim, int spacedim>
   inline const FEValuesBase<dim, spacedim>&
   IntegrationInfo<dim,spacedim>::fe_values() const
index dadacaa15db199a289ecb8a9e8666e71a06a677f..b22fbf024c8792bda375abf5d7b77d215779fcc9 100644 (file)
@@ -119,35 +119,6 @@ namespace MeshWorker
   }
   
 
-  template<int dim, int sdim>
-  template <class FEVALUES>
-  void
-  IntegrationInfo<dim,sdim>::initialize(
-    const FiniteElement<dim,sdim>& el,
-    const Mapping<dim,sdim>& mapping,
-    const Quadrature<FEVALUES::integral_dimension>& quadrature,
-    const UpdateFlags flags,
-    const BlockInfo* block_info)
-  {
-    if (block_info == 0 || block_info->local().size() == 0)
-      {
-       fevalv.resize(1);             
-       fevalv[0] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
-         new FEVALUES (mapping, el, quadrature, flags));
-      }
-    else
-      {
-       fevalv.resize(el.n_base_elements());
-       for (unsigned int i=0;i<fevalv.size();++i)
-         {
-           fevalv[i] = boost::shared_ptr<FEValuesBase<dim,sdim> > (
-             new FEVALUES (mapping, el.base_element(i), quadrature, flags));
-         }
-      }
-    n_components = el.n_components();
-  }
-  
-
   template<int dim, int sdim>
   void
   IntegrationInfo<dim,sdim>::initialize_data(
@@ -212,8 +183,9 @@ namespace MeshWorker
 
 
   template<int dim, int sdim>
+  template <typename number>
   void
-  IntegrationInfo<dim,sdim>::fill_local_data(const DoFInfo<dim, sdim>& info, bool split_fevalues)
+  IntegrationInfo<dim,sdim>::fill_local_data(const DoFInfo<dim, sdim, number>& info, bool split_fevalues)
   {
      if (split_fevalues)
        {
index c53b15956ab20150e00f8ad8d81020ce53df8672..83c742746fc8b4d7be92445d70d10bcf442a4a88 100644 (file)
@@ -26,25 +26,30 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace MeshWorker
 {
+  template class IntegrationInfo<deal_II_dimension, deal_II_dimension>;
+  
   template class LocalResults<float>;
   template class DoFInfo<deal_II_dimension,deal_II_dimension,float>;
+  template void IntegrationInfo<deal_II_dimension>::fill_local_data(
+    const DoFInfo<deal_II_dimension, deal_II_dimension, float>&, bool);
   
   template class LocalResults<double>;
   template class DoFInfo<deal_II_dimension,deal_II_dimension,double>;
-  template class IntegrationInfo<deal_II_dimension>;
+  template void IntegrationInfo<deal_II_dimension>::fill_local_data(
+    const DoFInfo<deal_II_dimension, deal_II_dimension, double>&, bool);
   
-  template void IntegrationInfo<deal_II_dimension>
-  ::initialize<FEValues<deal_II_dimension> >(
-    const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-    const Quadrature<FEValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-  template void IntegrationInfo<deal_II_dimension>
-  ::initialize<FEFaceValues<deal_II_dimension> >(
-    const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-    const Quadrature<FEFaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
-  template void IntegrationInfo<deal_II_dimension>
-  ::initialize<FESubfaceValues<deal_II_dimension> >(
-    const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
-    const Quadrature<FESubfaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
+//   template void IntegrationInfo<deal_II_dimension>
+//   ::initialize<FEValues<deal_II_dimension> >(
+//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
+//     const Quadrature<FEValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
+//   template void IntegrationInfo<deal_II_dimension>
+//   ::initialize<FEFaceValues<deal_II_dimension> >(
+//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
+//     const Quadrature<FEFaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
+//   template void IntegrationInfo<deal_II_dimension>
+//   ::initialize<FESubfaceValues<deal_II_dimension> >(
+//     const FiniteElement<deal_II_dimension>&, const Mapping<deal_II_dimension>&,
+//     const Quadrature<FESubfaceValues<deal_II_dimension>::integral_dimension>&, const UpdateFlags, const BlockInfo*);
 }
 
 #endif
index 6a0a6db4194f5e71799a23b2a8a68a534c4ae861..93d26832fae22430cc0de32167f82420a9724002 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -14,6 +14,9 @@
 #define __deal2__matrix_block_h
 
 #include <base/config.h>
+#include <base/named_data.h>
+
+#include <boost/shared_ptr.hpp>
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -90,6 +93,41 @@ struct MatrixBlock
 };
 
 
+/**
+ * A vector of MatrixBlock, which is implemented using shared
+ * pointers, in order to allow for copying and rearranging. Each
+ * matrix block can be identified by name.
+ *
+ * @author Baerbel Janssen, Guido Kanschat, 2010
+ */
+template <class MATRIX>
+class MatrixBlockVector : public NamedData<boost::shared_ptr<MatrixBlock<MATRIX> > >
+{
+  public:
+                                    /**
+                                     * Add a new matrix block at the
+                                     * position <tt(row,column)</tt>
+                                     * in the block system.
+                                     */
+    void add(unsigned int row, unsigned int column, const std::string& name);
+                                    /**
+                                     * Access a constant reference to
+                                     * the block at position <i>i</i>.
+                                     */
+    const MatrixBlock<MATRIX>& block(unsigned int i) const;
+                                    /**
+                                     * Access the matrix at position
+                                     * <i>i</i> for read and write
+                                     * access.
+                                     */
+    MATRIX& matrix(unsigned int i);
+    
+    
+};
+
+
+//----------------------------------------------------------------------//
+
 template <class MATRIX>
 MatrixBlock<MATRIX>::MatrixBlock()
                :
@@ -113,6 +151,33 @@ MatrixBlock<MATRIX>::MatrixBlock(unsigned int i, unsigned int j)
                row(i), column(j)
 {}
 
+//----------------------------------------------------------------------//
+
+template <class MATRIX>
+inline void
+MatrixBlockVector<MATRIX>::add(unsigned int row, unsigned int column, const std::string& name)
+{
+  boost::shared_ptr<MatrixBlock<MATRIX> > p(new MatrixBlock<MATRIX>(row, column));
+  NamedData<boost::shared_ptr<MatrixBlock<MATRIX> > >::add(p, name);
+}
+
+
+template <class MATRIX>
+inline const MatrixBlock<MATRIX>&
+MatrixBlockVector<MATRIX>::block(unsigned int i) const
+{
+  return *this->read(i);
+}
+
+
+template <class MATRIX>
+inline MATRIX&
+MatrixBlockVector<MATRIX>::matrix(unsigned int i)
+{
+  return (*this)(i)->matrix;
+}
+
+
 
 DEAL_II_NAMESPACE_CLOSE
 

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.