]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
additional functions for DG-elements
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 20 Jun 2000 22:39:59 +0000 (22:39 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 20 Jun 2000 22:39:59 +0000 (22:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@3038 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 210f2a6e2f6939200ed762551d6d2843f5d97d53..224483b997c498093bea386ec25f812ce5c2ffb5 100644 (file)
@@ -27,7 +27,7 @@
  * All member functions are static, so there is no need to create an
  * object of class #MGDoFTools#.
  *
- * @author Wolfgang Bangerth and others, 1999
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2000
  */
 class MGDoFTools 
 {
@@ -64,6 +64,19 @@ class MGDoFTools
                                SparsityPattern         &sparsity,
                                const unsigned int       level);
 
+                                    /**
+                                     * Create sparsity pattern for
+                                     * the fluxes at refinement
+                                     * edges. 
+                                     * @see{make_flux_sparsity_pattern}
+                                     * @see{DoFTools}
+                                     */
+    template <int dim>
+    static void
+    make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof_handler,
+                                    SparsityPattern         &sparsity,
+                                    const unsigned int       level);
+
                                     /**
                                      * This function does the same as
                                      * the other with the same name,
index d24dc56122636fab1e2572eabe35282d3ce2d78c..6cbb5db9caeacbf3f92e9361135b06ebe58eea9b 100644 (file)
@@ -406,13 +406,14 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                      for (unsigned int i=0; i<total_dofs; ++i)
                        {
                          for (unsigned int j=0; j<total_dofs; ++j)
-                           if (flux_dof_mask[i][j])
-                             {
+                           {
+                             if (flux_dof_mask[i][j])
                                sparsity.add (dofs_on_this_cell[i],
                                              dofs_on_other_cell[j]);
+                             if (flux_dof_mask[j][i])
                                sparsity.add (dofs_on_other_cell[i],
                                              dofs_on_this_cell[j]);
-                             }
+                           }
                        }
                      sub_neighbor->face(neighbor_face)->set_user_flag ();
                    }
@@ -421,10 +422,11 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
                  for (unsigned int i=0; i<total_dofs; ++i)
                    {
                      for (unsigned int j=0; j<total_dofs; ++j)
-                       if (flux_dof_mask[i][j])
-                         {
+                       {
+                         if (flux_dof_mask[i][j])
                            sparsity.add (dofs_on_this_cell[i],
                                          dofs_on_other_cell[j]);
+                         if (flux_dof_mask[j][i])
                            sparsity.add (dofs_on_other_cell[i],
                                          dofs_on_this_cell[j]);
                        }
index f151345986873c70d9148d4d99403d2a9079f5f4..6da824e9b2477961dc68833516205fbde7fddec5 100644 (file)
@@ -48,6 +48,8 @@ void MGDoFTools::make_sparsity_pattern (const MGDoFHandler<dim> &dof,
     }
 }
 
+
+
 template<int dim>
 void
 MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
@@ -103,6 +105,59 @@ MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
 
 
 
+template<int dim>
+void
+MGDoFTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof,
+                                            SparsityPattern       &sparsity,
+                                            const unsigned int level)
+{
+  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);
+  
+  Assert (sparsity.n_rows() == coarse_dofs,
+         ExcDimensionMismatch (sparsity.n_rows(), coarse_dofs));
+  Assert (sparsity.n_cols() == fine_dofs,
+         ExcDimensionMismatch (sparsity.n_cols(), fine_dofs));
+
+  const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
+  vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
+  vector<unsigned int> dofs_on_other_cell(dofs_per_cell);
+  MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
+                                  endc = dof.end(level);
+  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)
+       {
+         if ( (! cell->at_boundary(face)) &&
+              (static_cast<unsigned int>(cell->neighbor_level(face)) != level) )
+           {
+             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)
+                   {
+                     sparsity.add (dofs_on_this_cell[i],
+                                   dofs_on_other_cell[j]);
+                     sparsity.add (dofs_on_this_cell[j],
+                                   dofs_on_other_cell[i]);
+                   }
+               }
+           }
+       } 
+    }
+}
+
+
+
 template<int dim>
 void
 MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
@@ -209,6 +264,10 @@ template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<deal_II
                                                      SparsityPattern &,
                                                      const unsigned int);
 
+template void MGDoFTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<deal_II_dimension> &,
+                                                          SparsityPattern &,
+                                                          const unsigned int);
+
 #if deal_II_dimension > 1
 template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<deal_II_dimension> &,
                                                      SparsityPattern &,

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.