]> https://gitweb.dealii.org/ - dealii.git/commitdiff
FullMatrix arguments replaced by more obvious enum
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 9 Mar 2005 03:05:35 +0000 (03:05 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 9 Mar 2005 03:05:35 +0000 (03:05 +0000)
git-svn-id: https://svn.dealii.org/trunk@10052 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/include/multigrid/mg_dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/multigrid/mg_dof_tools.cc
tests/bits/dof_tools_18a.cc
tests/bits/dof_tools_18b.cc
tests/bits/dof_tools_18c.cc
tests/bits/dof_tools_18d.cc

index fc1f71d3371aa8224f88b511b4dfdc86a47fb17e..49065cbc18e918fa3458f3cdc57ffff8b8bc8a7d 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998 - 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 #include <set>
 #include <map>
 
+template<int dim, class T> class Table;
 class SparsityPattern;
 template <typename number> class Vector;
-template <typename number> class FullMatrix;
 template <int dim> class Function;
 template <int dim> class Point;
+template <int dim> class FiniteElement;
 template <int dim> class DoFHandler;
 template <int dim> class MGDoFHandler;
 class ConstraintMatrix;
@@ -141,11 +142,30 @@ template <int dim> class Mapping;
  * don't do it here) and not relevant in this context.
  *
  *
- * @author Wolfgang Bangerth and others, 1998, 1999, 2000, 2001
+ * @author Wolfgang Bangerth, Guido Kanschat and others, 1998 - 2005
  */
 class DoFTools
 {
   public:
+                                    /**
+                                     * The flags used in tables by certain
+                                     * <tt>make_*_pattern</tt>
+                                     * functions to determine whether
+                                     * two components of the solution
+                                     * couple in the interior of mesh
+                                     * cells or at the boundary.
+                                     */
+    enum Coupling
+    {
+                                          /// The two components do not couple
+         none,
+                                          /// The two components do couple
+         always,
+/// The two components couple only if their shape functions can be nonzero on this face
+         nonzero
+    };
+    
+    
                                     /**
                                      * Locate non-zero entries of the
                                      * system matrix.
@@ -471,8 +491,8 @@ class DoFTools
     static void
     make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                SparsityPattern       &sparsity,
-                               const FullMatrix<double>& int_mask,
-                               const FullMatrix<double>& flux_mask);
+                               const Table<2,Coupling>& int_mask,
+                               const Table<2,Coupling>& flux_mask);
     
                                     /**
                                      * Make up the constraints which
@@ -1386,6 +1406,18 @@ class DoFTools
                                 std::vector<std::map<unsigned int, float> > &weights,
                                 const typename DoFHandler<dim>::active_cell_iterator &begin,
                                 const typename DoFHandler<dim>::active_cell_iterator &end);
+
+                                    /**
+                                     * Compute coupling of dofs from
+                                     * coupling of components.
+                                     */
+    template <int dim>
+    static void compute_dof_couplings (
+      Table<2,Coupling>& dof_couplings,
+      const Table<2,Coupling>& component_couplings,
+      const FiniteElement<dim>& fe);
+
+    friend class MGTools;
 };
 
 
index 746b3ac80d0530e985a0a2b71e7eb6da26075b05..fdd172fe83d48bf7f6c714591ad0474e3ea1148d 100644 (file)
@@ -1,4 +1,4 @@
-//----------------------------  mg_dof_tools.h  ---------------------------
+//---------------------------------------------------------------------------
 //    $Id$
 //    Version: $Name$
 //
@@ -9,11 +9,12 @@
 //    to the file deal.II/doc/license.html for the  text  and
 //    further information on this license.
 //
-//----------------------------  mg_dof_tools.h  ---------------------------
+//---------------------------------------------------------------------------
 #ifndef __deal2__mg_dof_tools_h
 #define __deal2__mg_dof_tools_h
 
 #include <base/config.h>
+#include <dofs/dof_tools.h>
 #include <multigrid/mg_dof_handler.h>
 
 #include <vector>
@@ -36,9 +37,9 @@ template <class number> class FullMatrix;
  * for more information.
  *
  * All member functions are static, so there is no need to create an
- * object of class @p MGTools.
+ * object of class MGTools.
  *
- * @author Wolfgang Bangerth, Guido Kanschat, 1999-2003
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999 - 2005
  */
 class MGTools
 {
@@ -113,8 +114,27 @@ class MGTools
     make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
                                SparsityPattern       &sparsity,
                                const unsigned int level,
-                               const FullMatrix<double> &int_mask,
-                               const FullMatrix<double> &flux_mask);
+                               const Table<2,DoFTools::Coupling> &int_mask,
+                               const Table<2,DoFTools::Coupling> &flux_mask);
+
+                                    /**
+                                     * Create sparsity pattern for
+                                     * the fluxes at refinement
+                                     * edges. The matrix maps a
+                                     * function of the fine level
+                                     * space @p level to the coarser
+                                     * space. This is the version
+                                     * restricting the pattern to the
+                                     * elements actually needed.
+                                     *
+                                     * make_flux_sparsity_pattern()
+                                     */
+    template <int dim, class SparsityPattern>
+    static void
+    make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof_handler,
+                                    SparsityPattern         &sparsity,
+                                    const unsigned int       level,
+                                    const Table<2,DoFTools::Coupling> &flux_mask);
 
                                     /**
                                      * Count the dofs component-wise
index e626c0131f3208123b9f91ab9e3e1d916b006756..4a3a297fdb9c68867d9faa09ac519024f7fa2da7 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -15,6 +15,7 @@
 #include <base/multithread_info.h>
 #include <base/thread_management.h>
 #include <base/quadrature_lib.h>
+#include <base/table.h>
 #include <grid/tria.h>
 #include <grid/tria_iterator.h>
 #include <grid/intergrid_map.h>
@@ -472,10 +473,10 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<1> &dof,
 
 template <int dim, class SparsityPattern>
 void
-DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
-                                     SparsityPattern       &sparsity,
-                                     const FullMatrix<double> &int_mask,
-                                     const FullMatrix<double> &flux_mask)
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim>dof,
+                                     SparsityPattern&       sparsity,
+                                     const Table<2,Coupling>& int_mask,
+                                     const Table<2,Coupling>& flux_mask)
 {
   const unsigned int n_dofs = dof.n_dofs();
   const FiniteElement<dim>& fe = dof.get_fe();
@@ -485,14 +486,14 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
          ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
   Assert (sparsity.n_cols() == n_dofs,
          ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
-  Assert (int_mask.m() == n_comp,
-         ExcDimensionMismatch (int_mask.m(), n_comp));
-  Assert (int_mask.n() == n_comp,
-         ExcDimensionMismatch (int_mask.n(), n_comp));
-  Assert (flux_mask.m() == n_comp,
-         ExcDimensionMismatch (flux_mask.m(), n_comp));
-  Assert (flux_mask.n() == n_comp,
-         ExcDimensionMismatch (flux_mask.n(), n_comp));
+  Assert (int_mask.n_rows() == n_comp,
+         ExcDimensionMismatch (int_mask.n_rows(), n_comp));
+  Assert (int_mask.n_cols() == n_comp,
+         ExcDimensionMismatch (int_mask.n_cols(), n_comp));
+  Assert (flux_mask.n_rows() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+  Assert (flux_mask.n_cols() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
   
   const unsigned int total_dofs = fe.dofs_per_cell;
   std::vector<unsigned int> dofs_on_this_cell(total_dofs);
@@ -501,56 +502,16 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
   
   typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
                                                 endc = dof.end();
-
-
-//TODO[GK]: would someone explain what the letters 'f', 'a', and 'n' mean?
-//Even better would, of course, be to change this to an enum
-  Table<2,char> int_dof_mask(total_dofs, total_dofs);
-  Table<2,char> flux_dof_mask(total_dofs, total_dofs);
-
+  
+  Table<2,Coupling> int_dof_mask(total_dofs, total_dofs);
+  Table<2,Coupling> flux_dof_mask(total_dofs, total_dofs);
+  
+  compute_dof_couplings(int_dof_mask, int_mask, fe);
+  compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  
   for (unsigned int i=0; i<total_dofs; ++i)
-    {
-      const unsigned int ii
-        = (dof.get_fe().is_primitive(i) ?
-           dof.get_fe().system_to_component_index(i).first
-           :
-           (std::find (dof.get_fe().get_nonzero_components(i).begin(),
-                       dof.get_fe().get_nonzero_components(i).end(),
-                       true)
-            -
-           dof.get_fe().get_nonzero_components(i).begin())
-           );
-      Assert (ii < dof.get_fe().n_components(),
-              ExcInternalError());
-
-      for (unsigned int j=0; j<total_dofs; ++j)
-       {
-          const unsigned int jj
-            = (dof.get_fe().is_primitive(j) ?
-               dof.get_fe().system_to_component_index(j).first
-               :
-               (std::find (dof.get_fe().get_nonzero_components(j).begin(),
-                           dof.get_fe().get_nonzero_components(j).end(),
-                           true)
-                -
-                dof.get_fe().get_nonzero_components(j).begin())
-               );
-          Assert (jj < dof.get_fe().n_components(),
-                  ExcInternalError());          
-
-//TODO[GK]: the documentation only says that non-zeroness counts, but here
-//individual values are tested. this should be documented. Even better, rather
-//than passing in a matrix of doubles, we should pass in a Table<2,some-enum>
-         if (int_mask(ii,jj) != 0)
-           int_dof_mask(i,j) = 'f';
-         if (flux_mask(ii,jj) == 1.)
-           flux_dof_mask(i,j) = 'a';
-         if (flux_mask(ii,jj) == 2.)
-           flux_dof_mask(i,j) = 'f';
-       }
-      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
-       support_on_face(i,f) = fe.has_support_on_face(i,f);
-    }
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+      support_on_face(i,f) = fe.has_support_on_face(i,f);
   
                                   // Clear user flags because we will
                                   // need them. But first we save
@@ -570,7 +531,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                       // make sparsity pattern for this cell
       for (unsigned int i=0; i<total_dofs; ++i)
        for (unsigned int j=0; j<total_dofs; ++j)
-         if (int_dof_mask(i,j))
+         if (int_dof_mask(i,j) != none)
            sparsity.add (dofs_on_this_cell[i],
                          dofs_on_this_cell[j]);
 
@@ -593,13 +554,13 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                    {
                      const bool j_non_zero_i = support_on_face (j, face);
                      
-                     if (flux_dof_mask(i,j) == 'f')
+                     if (flux_dof_mask(i,j) == always)
                         sparsity.add (dofs_on_this_cell[i],
                                       dofs_on_this_cell[j]);
-                     if (flux_dof_mask(i,j) == 'a')
-                        if (i_non_zero_i && j_non_zero_i)
-                          sparsity.add (dofs_on_this_cell[i],
-                                        dofs_on_this_cell[j]);
+                     if (flux_dof_mask(i,j) == nonzero
+                         && i_non_zero_i && j_non_zero_i)
+                       sparsity.add (dofs_on_this_cell[i],
+                                     dofs_on_this_cell[j]);
                    }
                }
            }
@@ -611,10 +572,10 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                               // by coarser cells
              if (neighbor->level() < cell->level())
                continue;
-
+             
              const unsigned int
                 neighbor_face = cell->neighbor_of_neighbor(face);
-
+             
              if (neighbor->has_children())
                {
                  for (unsigned int sub_nr = 0;
@@ -634,11 +595,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                            {
                              const bool j_non_zero_i = support_on_face (j, face);
                              const bool j_non_zero_e  =support_on_face (j, neighbor_face);
-
-                             if (flux_dof_mask(i,j)==0)
-                               flux_dof_mask(i,j) = 'n';
-                         
-                             if (flux_dof_mask(i,j) == 'f')
+                             if (flux_dof_mask(i,j) == always)
                                {
                                  sparsity.add (dofs_on_this_cell[i],
                                                dofs_on_other_cell[j]);
@@ -649,7 +606,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                  sparsity.add (dofs_on_other_cell[i],
                                                dofs_on_other_cell[j]);
                                }
-                             if (flux_dof_mask(i,j) == 'a')
+                             if (flux_dof_mask(i,j) == nonzero)
                                {
                                  if (i_non_zero_i && j_non_zero_e)
                                    sparsity.add (dofs_on_this_cell[i],
@@ -665,7 +622,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                                  dofs_on_other_cell[j]);
                                }
                              
-                             if (flux_dof_mask(j,i) == 'f')
+                             if (flux_dof_mask(j,i) == always)
                                {
                                  sparsity.add (dofs_on_this_cell[j],
                                                dofs_on_other_cell[i]);
@@ -676,7 +633,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                  sparsity.add (dofs_on_other_cell[j],
                                                dofs_on_other_cell[i]);
                                }
-                             if (flux_dof_mask(j,i) == 'a')
+                             if (flux_dof_mask(j,i) == nonzero)
                                {
                                  if (j_non_zero_i && i_non_zero_e)
                                    sparsity.add (dofs_on_this_cell[j],
@@ -707,10 +664,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                        {
                          const bool j_non_zero_i = support_on_face (j, face);
                          const bool j_non_zero_e = support_on_face (j, neighbor_face);
-                         if (flux_dof_mask(i,j)==0)
-                           flux_dof_mask(i,j) = 'n';
-                         
-                         if (flux_dof_mask(i,j) == 'f')
+                         if (flux_dof_mask(i,j) == always)
                            {
                              sparsity.add (dofs_on_this_cell[i],
                                            dofs_on_other_cell[j]);
@@ -721,7 +675,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                              sparsity.add (dofs_on_other_cell[i],
                                            dofs_on_other_cell[j]);
                            }
-                         if (flux_dof_mask(i,j) == 'a')
+                         if (flux_dof_mask(i,j) == nonzero)
                            {
                              if (i_non_zero_i && j_non_zero_e)
                                sparsity.add (dofs_on_this_cell[i],
@@ -737,7 +691,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                                              dofs_on_other_cell[j]); 
                            }
 
-                         if (flux_dof_mask(j,i) == 'f')
+                         if (flux_dof_mask(j,i) == always)
                            {
                              sparsity.add (dofs_on_this_cell[j],
                                            dofs_on_other_cell[i]);
@@ -748,7 +702,7 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                              sparsity.add (dofs_on_other_cell[j],
                                            dofs_on_other_cell[i]);
                            }
-                         if (flux_dof_mask(j,i) == 'a')
+                         if (flux_dof_mask(j,i) == nonzero)
                            {
                              if (j_non_zero_i && i_non_zero_e)
                                sparsity.add (dofs_on_this_cell[j],
@@ -2817,6 +2771,62 @@ DoFTools::map_dofs_to_support_points (const Mapping<dim>       &mapping,
 }
 
 
+template <int dim>
+void
+DoFTools::compute_dof_couplings(
+  Table<2,Coupling>&        dof_couplings,
+  const Table<2,Coupling>&  component_couplings,
+  const FiniteElement<dim>& fe)
+{
+  Assert(component_couplings.n_rows() == fe.n_components(),
+        ExcDimensionMismatch(component_couplings.n_rows(),
+                             fe.n_components()));
+  Assert(component_couplings.n_cols() == fe.n_components(),
+        ExcDimensionMismatch(component_couplings.n_cols(),
+                             fe.n_components()));
+  
+  const unsigned int n_dofs = fe.dofs_per_cell;
+
+  Assert(dof_couplings.n_rows() == n_dofs,
+        ExcDimensionMismatch(dof_couplings.n_rows(), n_dofs));
+  Assert(dof_couplings.n_cols() == n_dofs,
+        ExcDimensionMismatch(dof_couplings.n_cols(), n_dofs));
+
+  for (unsigned int i=0; i<n_dofs; ++i)
+    {
+      const unsigned int ii
+        = (fe.is_primitive(i) ?
+           fe.system_to_component_index(i).first
+           :
+           (std::find (fe.get_nonzero_components(i).begin(),
+                       fe.get_nonzero_components(i).end(),
+                       true)
+            -
+           fe.get_nonzero_components(i).begin())
+           );
+      Assert (ii < fe.n_components(),
+              ExcInternalError());
+
+      for (unsigned int j=0; j<n_dofs; ++j)
+       {
+          const unsigned int jj
+            = (fe.is_primitive(j) ?
+               fe.system_to_component_index(j).first
+               :
+               (std::find (fe.get_nonzero_components(j).begin(),
+                           fe.get_nonzero_components(j).end(),
+                           true)
+                -
+                fe.get_nonzero_components(j).begin())
+               );
+          Assert (jj < fe.n_components(),
+                  ExcInternalError());          
+
+         dof_couplings(i,j) = component_couplings(ii,jj);
+       }
+    }
+}
+
 
 // explicit instantiations
 template void
@@ -3001,26 +3011,26 @@ template void
 DoFTools::make_flux_sparsity_pattern<deal_II_dimension,SparsityPattern>
 (const DoFHandler<deal_II_dimension>& dof,
  SparsityPattern    &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
 template void
 DoFTools::make_flux_sparsity_pattern<deal_II_dimension,CompressedSparsityPattern>
 (const DoFHandler<deal_II_dimension>& dof,
  CompressedSparsityPattern    &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
 template void
 DoFTools::make_flux_sparsity_pattern<deal_II_dimension,BlockSparsityPattern>
 (const DoFHandler<deal_II_dimension>& dof,
  BlockSparsityPattern    &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
 template void
 DoFTools::make_flux_sparsity_pattern<deal_II_dimension,CompressedBlockSparsityPattern>
 (const DoFHandler<deal_II_dimension>& dof,
  CompressedBlockSparsityPattern    &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
 #endif
 
 
@@ -3138,3 +3148,10 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension>
 (const Mapping<deal_II_dimension>       &mapping,
  const DoFHandler<deal_II_dimension>    &dof_handler,
  std::vector<Point<deal_II_dimension> > &support_points);
+
+template
+void
+DoFTools::compute_dof_couplings<deal_II_dimension>(
+  Table<2,Coupling>&       dof_couplings,
+  const Table<2,Coupling>& component_couplings,
+  const FiniteElement<deal_II_dimension>&      fe);
index 268613979e4a737b8b5ca1c1ffb997162bf4bcab..6c0d3ee1f9864e519eb3f1b02a74f12f92698732 100644 (file)
@@ -1,15 +1,15 @@
-//----------------------------  mg_dof_tools.cc  ---------------------------
+//---------------------------------------------------------------------------
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998 - 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 //    to the file deal.II/doc/license.html for the  text  and
 //    further information on this license.
 //
-//----------------------------  mg_dof_tools.cc  ---------------------------
+//---------------------------------------------------------------------------
 
 
 #include <base/multithread_info.h>
@@ -181,54 +181,52 @@ MGTools::make_flux_sparsity_pattern_edge (
 
 
 
-//TODO[GK]: Replace FullMatrix by Table<2,char>
+
+
+
 template<int dim, class SparsityPattern>
 void
 MGTools::make_flux_sparsity_pattern (
   const MGDoFHandler<dim> &dof,
   SparsityPattern       &sparsity,
   const unsigned int level,
-  const FullMatrix<double> &int_mask,
-  const FullMatrix<double> &flux_mask)
+  const Table<2,DoFTools::Coupling> &int_mask,
+  const Table<2,DoFTools::Coupling> &flux_mask)
 {
+  const FiniteElement<dim>& fe = dof.get_fe();
   const unsigned int n_dofs = dof.n_dofs(level);
-  const unsigned int n_comp = dof.get_fe().n_components();
+  const unsigned int n_comp = fe.n_components();
   
   Assert (sparsity.n_rows() == n_dofs,
          ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
   Assert (sparsity.n_cols() == n_dofs,
          ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
-  Assert (int_mask.m() == n_comp,
-         ExcDimensionMismatch (int_mask.m(), n_comp));
-  Assert (int_mask.n() == n_comp,
-         ExcDimensionMismatch (int_mask.n(), n_comp));
-  Assert (flux_mask.m() == n_comp,
-         ExcDimensionMismatch (flux_mask.m(), n_comp));
-  Assert (flux_mask.n() == n_comp,
-         ExcDimensionMismatch (flux_mask.n(), n_comp));
+  Assert (int_mask.n_rows() == n_comp,
+         ExcDimensionMismatch (int_mask.n_rows(), n_comp));
+  Assert (int_mask.n_cols() == n_comp,
+         ExcDimensionMismatch (int_mask.n_cols(), n_comp));
+  Assert (flux_mask.n_rows() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+  Assert (flux_mask.n_cols() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
   
-  const unsigned int total_dofs = dof.get_fe().dofs_per_cell;
+  const unsigned int total_dofs = fe.dofs_per_cell;
   std::vector<unsigned int> dofs_on_this_cell(total_dofs);
   std::vector<unsigned int> dofs_on_other_cell(total_dofs);
+  Table<2,bool> support_on_face(total_dofs, GeometryInfo<dim>::faces_per_cell);
+  
   typename MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
                                            endc = dof.end(level);
 
+  Table<2,DoFTools::Coupling> int_dof_mask(total_dofs, total_dofs);
+  Table<2,DoFTools::Coupling> flux_dof_mask(total_dofs, total_dofs);
 
-  std::vector<std::vector<bool> > int_dof_mask(total_dofs,
-                                std::vector<bool>(total_dofs, false));
-  std::vector<std::vector<bool> > flux_dof_mask(total_dofs,
-                                std::vector<bool>(total_dofs, false));
+  DoFTools::compute_dof_couplings(int_dof_mask, int_mask, fe);
+  DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  
   for (unsigned int i=0; i<total_dofs; ++i)
-    for (unsigned int j=0; j<total_dofs; ++j)
-      {
-       unsigned int ii = dof.get_fe().system_to_component_index(i).first;
-       unsigned int jj = dof.get_fe().system_to_component_index(j).first;
-       
-       if (int_mask(ii,jj) != 0)
-         int_dof_mask[i][j] = true;
-       if (flux_mask(ii,jj) != 0)
-         flux_dof_mask[i][j] = true;
-      }
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+      support_on_face(i,f) = fe.has_support_on_face(i,f);
   
                                   // Clear user flags because we will
                                   // need them. But first we save
@@ -248,7 +246,7 @@ MGTools::make_flux_sparsity_pattern (
                                       // make sparsity pattern for this cell
       for (unsigned int i=0; i<total_dofs; ++i)
        for (unsigned int j=0; j<total_dofs; ++j)
-         if (int_dof_mask[i][j])
+         if (int_dof_mask[i][j] != DoFTools::none)
            sparsity.add (dofs_on_this_cell[i],
                          dofs_on_this_cell[j]);
 
@@ -261,39 +259,191 @@ MGTools::make_flux_sparsity_pattern (
          if (cell_face->user_flag_set ())
            continue;
 
-         if (! cell->at_boundary (face) )
+         if (cell->at_boundary (face) )
+           {
+             for (unsigned int i=0; i<total_dofs; ++i)
+               {
+                 const bool i_non_zero_i = support_on_face (i, face);
+                 for (unsigned int j=0; j<total_dofs; ++j)
+                   {
+                     const bool j_non_zero_i = support_on_face (j, face);
+                     
+                     if (flux_dof_mask(i,j) == DoFTools::always)
+                        sparsity.add (dofs_on_this_cell[i],
+                                      dofs_on_this_cell[j]);
+                     if (flux_dof_mask(i,j) == DoFTools::nonzero
+                         && i_non_zero_i && j_non_zero_i)
+                       sparsity.add (dofs_on_this_cell[i],
+                                     dofs_on_this_cell[j]);
+                   }
+               }
+           }
+         else
            {
              typename MGDoFHandler<dim>::cell_iterator
                neighbor = cell->neighbor(face);
-
+             
              if (neighbor->level() < cell->level())
                continue;
-
+             
              unsigned int neighbor_face = cell->neighbor_of_neighbor(face);
-
+             
              neighbor->get_mg_dof_indices (dofs_on_other_cell);
              for (unsigned int i=0; i<total_dofs; ++i)
                {
+                 const bool i_non_zero_i = support_on_face (i, face);
+                 const bool i_non_zero_e = support_on_face (i, neighbor_face);
                  for (unsigned int j=0; j<total_dofs; ++j)
-                   if (flux_dof_mask[i][j])
-                     {
-                       sparsity.add (dofs_on_this_cell[i],
-                                     dofs_on_other_cell[j]);
-                       sparsity.add (dofs_on_other_cell[i],
-                                     dofs_on_this_cell[j]);
-                     }
+                   {
+                     const bool j_non_zero_i = support_on_face (j, face);
+                     const bool j_non_zero_e = support_on_face (j, neighbor_face);
+                     if (flux_dof_mask(i,j) == DoFTools::always)
+                       {
+                         sparsity.add (dofs_on_this_cell[i],
+                                       dofs_on_other_cell[j]);
+                         sparsity.add (dofs_on_other_cell[i],
+                                       dofs_on_this_cell[j]);
+                         sparsity.add (dofs_on_this_cell[i],
+                                       dofs_on_this_cell[j]);
+                         sparsity.add (dofs_on_other_cell[i],
+                                       dofs_on_other_cell[j]);
+                       }
+                     if (flux_dof_mask(i,j) == DoFTools::nonzero)
+                       {
+                         if (i_non_zero_i && j_non_zero_e)
+                           sparsity.add (dofs_on_this_cell[i],
+                                         dofs_on_other_cell[j]);
+                         if (i_non_zero_e && j_non_zero_i)
+                           sparsity.add (dofs_on_other_cell[i],
+                                         dofs_on_this_cell[j]);
+                         if (i_non_zero_i && j_non_zero_i)
+                           sparsity.add (dofs_on_this_cell[i],
+                                         dofs_on_this_cell[j]);
+                         if (i_non_zero_e && j_non_zero_e)
+                           sparsity.add (dofs_on_other_cell[i],
+                                         dofs_on_other_cell[j]); 
+                       }
+                     
+                     if (flux_dof_mask(j,i) == DoFTools::always)
+                       {
+                         sparsity.add (dofs_on_this_cell[j],
+                                       dofs_on_other_cell[i]);
+                         sparsity.add (dofs_on_other_cell[j],
+                                       dofs_on_this_cell[i]);
+                         sparsity.add (dofs_on_this_cell[j],
+                                       dofs_on_this_cell[i]);
+                         sparsity.add (dofs_on_other_cell[j],
+                                       dofs_on_other_cell[i]);
+                       }
+                     if (flux_dof_mask(j,i) == DoFTools::nonzero)
+                       {
+                         if (j_non_zero_i && i_non_zero_e)
+                           sparsity.add (dofs_on_this_cell[j],
+                                         dofs_on_other_cell[i]);
+                         if (j_non_zero_e && i_non_zero_i)
+                           sparsity.add (dofs_on_other_cell[j],
+                                         dofs_on_this_cell[i]);
+                         if (j_non_zero_i && i_non_zero_i)
+                           sparsity.add (dofs_on_this_cell[j],
+                                         dofs_on_this_cell[i]);
+                         if (j_non_zero_e && i_non_zero_e)
+                           sparsity.add (dofs_on_other_cell[j],
+                                         dofs_on_other_cell[i]); 
+                       }
+                   }
                }
              neighbor->face(neighbor_face)->set_user_flag (); 
            }
        }
     }
-
+  
                                   // finally restore the user flags
   const_cast<Triangulation<dim> &>(dof.get_tria()).load_user_flags(user_flags);
 }
 
 
 
+template<int dim, class SparsityPattern>
+void
+MGTools::make_flux_sparsity_pattern_edge (
+  const MGDoFHandler<dim> &dof,
+  SparsityPattern       &sparsity,
+  const unsigned int level,
+  const Table<2,DoFTools::Coupling> &flux_mask)
+{
+  const FiniteElement<dim>& fe = dof.get_fe();
+  const unsigned int n_comp = fe.n_components();
+  
+  Assert ((level>=1) && (level<dof.get_tria().n_levels()),
+         ExcIndexRange(level, 1, dof.get_tria().n_levels()));
+
+  const unsigned int fine_dofs = dof.n_dofs(level);
+  const unsigned int coarse_dofs = dof.n_dofs(level-1);
+  
+  // Matrix maps from fine level to coarse level
+
+  Assert (sparsity.n_rows() == coarse_dofs,
+         ExcDimensionMismatch (sparsity.n_rows(), coarse_dofs));
+  Assert (sparsity.n_cols() == fine_dofs,
+         ExcDimensionMismatch (sparsity.n_cols(), fine_dofs));
+  Assert (flux_mask.n_rows() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+  Assert (flux_mask.n_cols() == n_comp,
+         ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
+
+  const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
+  std::vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
+  std::vector<unsigned int> dofs_on_other_cell(dofs_per_cell);
+  Table<2,bool> support_on_face(total_dofs, GeometryInfo<dim>::faces_per_cell);
+  
+  typename MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
+                                           endc = dof.end(level);
+
+   Table<2,DoFTools::Coupling> flux_dof_mask(total_dofs, total_dofs);
+
+   DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  
+  for (unsigned int i=0; i<total_dofs; ++i)
+    for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+      support_on_face(i,f) = fe.has_support_on_face(i,f);
+  
+  for (; cell!=endc; ++cell)
+    {
+      cell->get_mg_dof_indices (dofs_on_this_cell);
+                                      // Loop over all interior neighbors
+      for (unsigned int face = 0;
+          face < GeometryInfo<dim>::faces_per_cell;
+          ++face)
+       {
+         // Neighbor is coarser
+
+         if ( (! cell->at_boundary(face)) &&
+              (static_cast<unsigned int>(cell->neighbor_level(face)) != level) )
+           {
+             typename MGDoFHandler<dim>::cell_iterator
+               neighbor = cell->neighbor(face);
+             neighbor->get_mg_dof_indices (dofs_on_other_cell);
+
+             for (unsigned int i=0; i<dofs_per_cell; ++i)
+               {
+                 for (unsigned int j=0; j<dofs_per_cell; ++j)
+                   {
+                     if (flux_dof_mask(i,j) != DoFTools::none)
+                       {
+                         sparsity.add (dofs_on_other_cell[i],
+                                       dofs_on_this_cell[j]);
+                         sparsity.add (dofs_on_other_cell[j],
+                                       dofs_on_this_cell[i]);
+                       }
+                   }
+               }
+           }
+       } 
+    }
+}
+
+
+
 template <int dim>
 void
 MGTools::count_dofs_per_component (
@@ -567,29 +717,29 @@ template void
 MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
                                       SparsityPattern &,
                                       const unsigned int,
-                                      const FullMatrix<double>&,
-                                      const FullMatrix<double>&);
+                                      const Table<2,DoFTools::Coupling>&,
+                                      const Table<2,DoFTools::Coupling>&);
 
 template void
 MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
                                       CompressedSparsityPattern &,
                                       const unsigned int,
-                                      const FullMatrix<double>&,
-                                      const FullMatrix<double>&);
+                                      const Table<2,DoFTools::Coupling>&,
+                                      const Table<2,DoFTools::Coupling>&);
 
 template void
 MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
                                       BlockSparsityPattern &,
                                       const unsigned int,
-                                      const FullMatrix<double>&,
-                                      const FullMatrix<double>&);
+                                      const Table<2,DoFTools::Coupling>&,
+                                      const Table<2,DoFTools::Coupling>&);
 
 template void
 MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
                                       CompressedBlockSparsityPattern &,
                                       const unsigned int,
-                                      const FullMatrix<double>&,
-                                      const FullMatrix<double>&);
+                                      const Table<2,DoFTools::Coupling>&,
+                                      const Table<2,DoFTools::Coupling>&);
 
 #endif
 
index 37a20e72594105f8d4a23d0b1923914d8e7d3988..d3cd963dc7b1ba4a5bd5c97238d123f8138299e3 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 2003, 2004 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -29,13 +29,13 @@ std::string output_file_name = "dof_tools_18a.output";
                                    // something)
 void
 make_masks (const unsigned int n,
-            FullMatrix<double> &m1,
-            FullMatrix<double> &m2)
+            Table<2,DoFTools::Coupling> &m1,
+            Table<2,DoFTools::Coupling> &m2)
 {
   m1.reinit (n,n);
   m2.reinit (n,n);
   for (unsigned int i=0; i<n; ++i)
-    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
 }
 
   
@@ -48,7 +48,8 @@ template <int dim>
 void
 check_this (const DoFHandler<dim> &dof_handler)
 {
-  FullMatrix<double> mask_int, mask_ext;
+  Table<2,DoFTools::Coupling> mask_int;
+  Table<2,DoFTools::Coupling> mask_ext;
   make_masks (dof_handler.get_fe().n_components(),
               mask_int, mask_ext);
   
index 5d8709d283eac6c81d7c68507fff5829a2713753..6c30eb26fa4090ad2ce3dea4282d65b0fb184bac 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 2003, 2004 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -28,13 +28,13 @@ std::string output_file_name = "dof_tools_18b.output";
                                    // something)
 void
 make_masks (const unsigned int n,
-            FullMatrix<double> &m1,
-            FullMatrix<double> &m2)
+            Table<2, DoFTools::Coupling> &m1,
+            Table<2, DoFTools::Coupling> &m2)
 {
   m1.reinit (n,n);
   m2.reinit (n,n);
   for (unsigned int i=0; i<n; ++i)
-    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
 }
 
   
@@ -48,7 +48,8 @@ template <int dim>
 void
 check_this (const DoFHandler<dim> &dof_handler)
 {
-  FullMatrix<double> mask_int, mask_ext;
+  Table<2, DoFTools::Coupling> mask_int;
+  Table<2, DoFTools::Coupling> mask_ext;
   make_masks (dof_handler.get_fe().n_components(),
               mask_int, mask_ext);
   
index e01a87a614c0b1214de20ddd4f3a953f94f53cc5..3686a17ea9979aea21c1e14a06ca578c4a7275aa 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 2003, 2004 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -28,13 +28,13 @@ std::string output_file_name = "dof_tools_18c.output";
                                    // something)
 void
 make_masks (const unsigned int n,
-            FullMatrix<double> &m1,
-            FullMatrix<double> &m2)
+            Table<2, DoFTools::Coupling> &m1,
+            Table<2, DoFTools::Coupling> &m2)
 {
   m1.reinit (n,n);
   m2.reinit (n,n);
   for (unsigned int i=0; i<n; ++i)
-    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
 }
 
   
@@ -56,7 +56,8 @@ check_this (const DoFHandler<dim> &dof_handler)
   if (dof_handler.get_fe().is_primitive() != true)
     return;
   
-  FullMatrix<double> mask_int, mask_ext;
+  Table<2, DoFTools::Coupling> mask_int;
+  Table<2, DoFTools::Coupling> mask_ext;
   make_masks (dof_handler.get_fe().n_components(),
               mask_int, mask_ext);
   
index 9837505c36543c75c6496a825132aa9ec345db70..2592a6687309528a9b0490f48fa0923b2219d79c 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$ 
 //
-//    Copyright (C) 2003, 2004 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -28,13 +28,13 @@ std::string output_file_name = "dof_tools_18d.output";
                                    // something)
 void
 make_masks (const unsigned int n,
-            FullMatrix<double> &m1,
-            FullMatrix<double> &m2)
+            Table<2, DoFTools::Coupling> &m1,
+            Table<2, DoFTools::Coupling> &m2)
 {
   m1.reinit (n,n);
   m2.reinit (n,n);
   for (unsigned int i=0; i<n; ++i)
-    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+    m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
 }
 
   
@@ -57,7 +57,8 @@ check_this (const DoFHandler<dim> &dof_handler)
   if (dof_handler.get_fe().is_primitive() != true)
     return;
   
-  FullMatrix<double> mask_int, mask_ext;
+  Table<2, DoFTools::Coupling> mask_int;
+  Table<2, DoFTools::Coupling> mask_ext;
   make_masks (dof_handler.get_fe().n_components(),
               mask_int, mask_ext);
   

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.