]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix a snafu from last night: the function I moved into an internal namespace is still...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 8 Jan 2009 16:01:33 +0000 (16:01 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 8 Jan 2009 16:01:33 +0000 (16:01 +0000)
git-svn-id: https://svn.dealii.org/trunk@18136 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/multigrid/mg_dof_tools.cc

index b76e1b40c26410f8788af5f7ff834cd0b86da5bc..1a928aad0d1dfe440a4bbed9fdc99db8a2f9652b 100644 (file)
@@ -1768,6 +1768,20 @@ class DoFTools
     convert_couplings_to_blocks (const DoFHandler<dim,spacedim>& dof_handler,
                                 const Table<2, Coupling>& table_by_component,
                                 std::vector<Table<2,Coupling> >& tables_by_block);
+
+                                    /**
+                                     * Given a finite element and a table how
+                                     * the vector components of it couple
+                                     * with each other, compute and return a
+                                     * table that describes how the
+                                     * individual shape functions couple with
+                                     * each other.
+                                     */
+    template <int dim, int spacedim>
+    static
+    Table<2,Coupling>
+    dof_couplings_from_component_couplings (const FiniteElement<dim,spacedim> &fe,
+                                           const Table<2,Coupling> &component_couplings);
     
                                     /**
                                      * Exception
@@ -1813,20 +1827,6 @@ class DoFTools
                                      * Exception
                                      */
     DeclException0 (ExcInvalidBoundaryIndicator);
-
-  private:
-
-                                    /**
-                                     * Compute coupling of dofs from
-                                     * coupling of components.
-                                     */
-    template <int dim, int spacedim>
-    static void compute_dof_couplings (
-      Table<2,Coupling>& dof_couplings,
-      const Table<2,Coupling>& component_couplings,
-      const FiniteElement<dim,spacedim>& fe);
-
-    friend class MGTools;
 };
 
 
index 58e319be5eaf71e069dff58928ac096300a0c051..1c38e132480e814e96963fdd1ce00b67219fb50c 100644 (file)
@@ -681,65 +681,57 @@ DoFTools::make_flux_sparsity_pattern (
 #endif
 
 
-namespace internal
+template <int dim, int spacedim>
+Table<2,DoFTools::Coupling>
+DoFTools::dof_couplings_from_component_couplings
+(const FiniteElement<dim,spacedim> &fe,
+ const Table<2,Coupling> &component_couplings)
 {
-  namespace
-  {
-    template <int dim, int spacedim>
-    void
-    compute_dof_couplings(Table<2,DoFTools::Coupling>&        dof_couplings,
-                         const Table<2,DoFTools::Coupling>&  component_couplings,
-                         const FiniteElement<dim,spacedim>& 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()));
+  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;
+  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));
+  Table<2,DoFTools::Coupling> dof_couplings (n_dofs, n_dofs);
 
-      for (unsigned int i=0; i<n_dofs; ++i)
+  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 ii
-           = (fe.is_primitive(i) ?
-              fe.system_to_component_index(i).first
+         const unsigned int jj
+           = (fe.is_primitive(j) ?
+              fe.system_to_component_index(j).first
               :
-              (std::find (fe.get_nonzero_components(i).begin(),
-                          fe.get_nonzero_components(i).end(),
+              (std::find (fe.get_nonzero_components(j).begin(),
+                          fe.get_nonzero_components(j).end(),
                           true)
                -
-               fe.get_nonzero_components(i).begin())
+               fe.get_nonzero_components(j).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());          
+         Assert (jj < fe.n_components(),
+                 ExcInternalError());          
 
-             dof_couplings(i,j) = component_couplings(ii,jj);
-           }
+         dof_couplings(i,j) = component_couplings(ii,jj);
        }
     }
-  }
+  return dof_couplings;
 }
 
 
@@ -778,11 +770,9 @@ DoFTools::make_flux_sparsity_pattern (
   typename DH::active_cell_iterator cell = dof.begin_active(),
                                    endc = dof.end();
   
-  Table<2,Coupling> int_dof_mask(total_dofs, total_dofs);
-  Table<2,Coupling> flux_dof_mask(total_dofs, total_dofs);
-  
-  internal::compute_dof_couplings(int_dof_mask, int_mask, fe);
-  internal::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  const Table<2,Coupling>
+    int_dof_mask  = dof_couplings_from_component_couplings(fe, int_mask),
+    flux_dof_mask = dof_couplings_from_component_couplings(fe, flux_mask);
   
   for (unsigned int i=0; i<total_dofs; ++i)
     for (unsigned int f=0; f<GeometryInfo<DH::dimension>::faces_per_cell;++f)
index 8d04955e8c35bf673ace8a4ebf1e828dbc0510b9..742e9d4d6de6a1c919daa45ca4d0b39e78fd3804 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -685,11 +685,9 @@ MGTools::make_flux_sparsity_pattern (
   typename MGDoFHandler<dim,spacedim>::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);
-
-  DoFTools::compute_dof_couplings(int_dof_mask, int_mask, fe);
-  DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  const Table<2,DoFTools::Coupling>
+    int_dof_mask  = DoFTools::dof_couplings_from_component_couplings(fe, int_mask),
+    flux_dof_mask = DoFTools::dof_couplings_from_component_couplings(fe, flux_mask);
   
   for (unsigned int i=0; i<total_dofs; ++i)
     for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
@@ -864,11 +862,10 @@ MGTools::make_flux_sparsity_pattern_edge (
   Table<2,bool> support_on_face(dofs_per_cell, GeometryInfo<dim>::faces_per_cell);
   
   typename MGDoFHandler<dim,spacedim>::cell_iterator cell = dof.begin(level),
-                                           endc = dof.end(level);
-
-   Table<2,DoFTools::Coupling> flux_dof_mask(dofs_per_cell, dofs_per_cell);
+                                                    endc = dof.end(level);
 
-   DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+  const Table<2,DoFTools::Coupling> flux_dof_mask
+    = DoFTools::dof_couplings_from_component_couplings(fe, flux_mask);
   
   for (unsigned int i=0; i<dofs_per_cell; ++i)
     for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)

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.