]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
implement new block structure in DoFInfo
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 12 Mar 2012 13:43:15 +0000 (13:43 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 12 Mar 2012 13:43:15 +0000 (13:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@25263 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/meshworker/dof_info.h
deal.II/include/deal.II/meshworker/dof_info.templates.h
deal.II/include/deal.II/meshworker/simple.h

index 93055edddaa83955a3ea58402bbda25c6325865c..40addebb6c6ab8061e34f15d5b72edc7f5879ac9 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors
+//    Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -93,9 +93,19 @@ namespace MeshWorker
                                        */
       unsigned int sub_number;
 
-                                      /// The DoF indices of the current cell
+                                      /*
+                                       * The DoF indices of the
+                                       * current cell
+                                       */
       std::vector<unsigned int> indices;
 
+                                      /**
+                                       * The DoF indices on the
+                                       * current cell, organized by
+                                       * local blocks
+                                       */
+      std::vector<std::vector<unsigned int> > indices_by_block;
+      
                                       /**
                                        * Constructor setting the
                                        * #block_info pointer.
index faa56529cea6ac3be81e4011e19682e7c050ad57..d43af1b6493a769e4e37ccf10c7fce565f3820f7 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2009, 2010, 2011 by the deal.II authors
+//    Copyright (C) 2009, 2010, 2011, 2012 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -22,7 +22,11 @@ namespace MeshWorker
   template <int dim, int spacedim, typename number>
   DoFInfo<dim,spacedim,number>::DoFInfo(const BlockInfo& info)
                  : block_info(&info, typeid(*this).name())
-  {}
+  {
+    indices_by_block.resize(info.local().size());
+    for (unsigned int i=0;i<indices_by_block.size();++i)
+      indices_by_block[i].resize(info.local().block_size(i));
+  }
 
 
   template <int dim, int spacedim, typename number>
@@ -44,6 +48,14 @@ namespace MeshWorker
          {
            indices_org.resize(c->get_fe().dofs_per_cell);
            c->get_dof_indices(indices_org);
+           for (unsigned int i=0;i<indices.size();++i)
+             {
+               const std::pair<unsigned int, unsigned int>
+                 bi = block_info->local().global_to_local(i);
+               indices_by_block[bi.first][bi.second] = indices_org[i];
+             }
+                                            // Remove this after
+                                            // changing block codes
            for (unsigned int i=0;i<indices.size();++i)
              indices[this->block_info->renumber(i)] = indices_org[i];
          }
index 56f90d101053b80d9a4fcfae2071a7bc186527c1..6e4078053498b89080c5c72cc355983c6b69a3d1 100644 (file)
@@ -115,6 +115,16 @@ namespace MeshWorker
  * this class can be used in a MeshWorker::loop() to assemble the cell
  * and face matrices into the global matrix.
  *
+ * The assembler can handle two different types of local data. First,
+ * by default, the obvious choice of taking a single local matrix with
+ * dimensions equal to the number of degrees of freedom of the
+ * cell. Alternatively, a local block structure can be initialized
+ * with initialize_local_blocks(). After this, the local data will be
+ * arranged as an array of n by n FullMatrix blocks, which are
+ * ordered lexicographically in DoFInfo. Note that
+ * initialize_local_blocks() has to be called before initialize_info()
+ * to take the desired effect.
+ *
  * @todo On locally refined meshes, a ConstraintMatrix should be used
  * to automatically eliminate hanging nodes.
  *
@@ -140,9 +150,45 @@ namespace MeshWorker
                                          */
        void initialize(MATRIX& m);
                                         /**
-                                         * Initialize the constraints.
+                                         * Initialize the
+                                         * constraints for laer use
+                                         * by the assemble functions. 
                                          */
         void initialize(const ConstraintMatrix& constraints);
+       
+                                        /**
+                                         * Store information on the
+                                         * local block structure. If
+                                         * the assembler is
+                                         * inititialized with this
+                                         * function,
+                                         * initialize_info() will
+                                         * generate one local matrix
+                                         * for each block row and
+                                         * column, whch will be
+                                         * numbered
+                                         * lexicographically, row by
+                                         * row.
+                                         *
+                                         * In spite of using local
+                                         * block structure, all
+                                         * blocks will be enteres
+                                         * into the same global
+                                         * matrix, disregarding any
+                                         * global block structure.
+                                         *
+                                         * @note The argument of this
+                                         * function will be copied
+                                         * into the member object
+                                         * #local_indices. Thus,
+                                         * every subsequent change in
+                                         * the block structure must
+                                         * be initialzied or will not
+                                         * be used by the assembler.
+                                         */
+       
+       void initialize_local_blocks(const BlockIndices& local_indices);
+       
                                         /**
                                          * Initialize the local data
                                          * in the
@@ -191,10 +237,21 @@ namespace MeshWorker
                                          */
        SmartPointer<MATRIX,MatrixSimple<MATRIX> > matrix;
                                         /**
-                                         * A pointer to the object containing constraints.
+                                         * A pointer to the object
+                                         * containing constraints.
                                          */
         SmartPointer<const ConstraintMatrix,MatrixSimple<MATRIX> > constraints;
 
+                                        /**
+                                         * The object containing the
+                                         * local block structure. Set
+                                         * by
+                                         * initialize_local_blocks()
+                                         * and used by assembling
+                                         * functions.
+                                         */
+       BlockIndices local_indices;
+       
                                         /**
                                          * The smallest positive
                                          * number that will be
@@ -579,7 +636,15 @@ namespace MeshWorker
       constraints = &c;
     }
 
+    
+    template <class MATRIX>
+    inline void
+    MatrixSimple<MATRIX>::initialize_local_blocks(const BlockIndices& b)
+    {
+      local_indices = b;
+    }
 
+    
     template <class MATRIX >
     template <class DOFINFO>
     inline void
@@ -598,7 +663,7 @@ namespace MeshWorker
     {
       AssertDimension(M.m(), i1.size());
       AssertDimension(M.n(), i2.size());
-
+      
       if(constraints == 0)
        {
          for (unsigned int j=0; j<i1.size(); ++j)
@@ -617,7 +682,11 @@ namespace MeshWorker
     inline void
     MatrixSimple<MATRIX>::assemble(const DOFINFO& info)
     {
-      assemble(info.matrix(0,false).matrix, info.indices, info.indices);
+      if (local_indices.size() == 0)
+       assemble(info.matrix(0,false).matrix, info.indices, info.indices);
+      else
+       {
+       }
     }
 
 
@@ -627,12 +696,18 @@ namespace MeshWorker
     MatrixSimple<MATRIX>::assemble(const DOFINFO& info1,
                                   const DOFINFO& info2)
     {
-      assemble(info1.matrix(0,false).matrix, info1.indices, info1.indices);
-      assemble(info1.matrix(0,true).matrix, info1.indices, info2.indices);
-      assemble(info2.matrix(0,false).matrix, info2.indices, info2.indices);
-      assemble(info2.matrix(0,true).matrix, info2.indices, info1.indices);
+      if (local_indices.size() == 0)
+       {
+         assemble(info1.matrix(0,false).matrix, info1.indices, info1.indices);
+         assemble(info1.matrix(0,true).matrix, info1.indices, info2.indices);
+         assemble(info2.matrix(0,false).matrix, info2.indices, info2.indices);
+         assemble(info2.matrix(0,true).matrix, info2.indices, info1.indices);
+       }
+      else
+       {
+       }
     }
-
+  
 
 //----------------------------------------------------------------------//
 

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.