]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
make_flux_sparsity_pattern in DofTools
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 3 Dec 1999 02:49:15 +0000 (02:49 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 3 Dec 1999 02:49:15 +0000 (02:49 +0000)
git-svn-id: https://svn.dealii.org/trunk@1970 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_handler.h
deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_handler.cc
deal.II/deal.II/source/dofs/dof_tools.cc

index 4eea08f89cc4c805880084e10757b3d5e9876d54..6bf15cb273f50756b052346affe04da75d1693fe 100644 (file)
@@ -449,7 +449,7 @@ class DoFHandler  :  public Subscriptor,
                                      * have to compress the matrix yourself,
                                      * using #SparseMatrixStruct::compress()#.
                                      */
-    void make_sparsity_pattern (SparseMatrixStruct &) const; 
+    void make_sparsity_pattern (SparseMatrixStruct &) const;
 
                                     /**
                                      * This function does mistly the same as
@@ -494,7 +494,7 @@ class DoFHandler  :  public Subscriptor,
                                      */
     void make_sparsity_pattern (const vector<vector<bool> > &mask,
                                SparseMatrixStruct          &) const;
-    
+
                                     /**
                                      * Write the sparsity structure of the
                                      * matrix composed of the basis functions
index f0248902415e3fb15b9e394db9b20ddd42edd530..5a0e3afb564a83437062ef7397462730ac572920 100644 (file)
 class DoFTools
 {
   public:
+                                    /**
+                                     * Generate sparsity pattern for fluxes.
+                                     * This is a replacement of the
+                                     * function
+                                     * #make_sparsity_pattern# for
+                                     * discontinuous methods. Since
+                                     * the fluxes include couplings
+                                     * between neighboring elements,
+                                     * the normal couplings and these
+                                     * extra matrix entries are
+                                     * considered.
+                                     */
+    template<int dim>
+    static void make_flux_sparsity_pattern (const DoFHandler<dim>&,
+                                    SparseMatrixStruct &);
+    
                                     /**
                                      * Extract the indices of the degrees
                                      * of freedom belonging to certain
index 7b480a0e4bccc8f11333ac0fb7aa4fc5538f6eae..b7b9589a09d9f43d0a46874c88dba7a696971596 100644 (file)
@@ -1974,7 +1974,8 @@ void DoFHandler<3>::make_hanging_node_constraints (ConstraintMatrix &constraints
 
 
 template <int dim>
-void DoFHandler<dim>::make_sparsity_pattern (SparseMatrixStruct &sparsity) const {
+void DoFHandler<dim>::make_sparsity_pattern (SparseMatrixStruct &sparsity) const
+{
   Assert (selected_fe != 0, ExcNoFESelected());
   Assert (sparsity.n_rows() == n_dofs(),
          ExcDifferentDimensions (sparsity.n_rows(), n_dofs()));
@@ -1993,10 +1994,8 @@ void DoFHandler<dim>::make_sparsity_pattern (SparseMatrixStruct &sparsity) const
        for (unsigned int j=0; j<dofs_per_cell; ++j)
          sparsity.add (dofs_on_this_cell[i],
                        dofs_on_this_cell[j]);
-    };
-};
-
-
+    }
+}
 
 template <int dim>
 void
index 7847191faacf0c0b27a22678a8ac3542f3a00e81..8a30654c59d4a271fa0269ff185bf24dca6f5c2c 100644 (file)
@@ -9,8 +9,76 @@
 #include <fe/fe.h>
 #include <fe/fe_system.h>
 #include <basic/dof_tools.h>
+#include <lac/sparsematrix.h>
 
+#if deal_II_dimension == 1
 
+//TODO: Implement make_flux_sparsity_pattern in 1d
+
+template <>
+void
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<1>&,
+                                     SparseMatrixStruct &)
+{
+  Assert(false, ExcNotImplemented());
+}
+
+#else
+
+//TODO: Check this function for potential of optimization. (G)
+
+template<int dim>
+void
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim>& dof,
+                                     SparseMatrixStruct &sparsity)
+{
+  const unsigned int n_dofs = dof.n_dofs();
+  
+  Assert (sparsity.n_rows() == n_dofs,
+         ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
+  Assert (sparsity.n_cols() == n_dofs,
+         ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
+
+  const unsigned int total_dofs = dof.get_fe().dofs_per_cell;
+  vector<int> dofs_on_this_cell(total_dofs);
+  vector<int> dofs_on_other_cell(total_dofs);
+  DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
+                                       endc = dof.end();
+  for (; cell!=endc; ++cell) 
+    {
+      cell->get_dof_indices (dofs_on_this_cell);
+                                      // make sparsity pattern for this cell
+      for (unsigned int i=0; i<total_dofs; ++i)
+       for (unsigned int j=0; j<total_dofs; ++j)
+         sparsity.add (dofs_on_this_cell[i],
+                       dofs_on_this_cell[j]);
+
+                                      // Loop over all interior neighbors
+      for (unsigned int face = 0;
+          face < GeometryInfo<dim>::faces_per_cell;
+          ++face)
+       {
+         if (cell->face(face)->boundary_indicator() == 255)
+           {
+             DoFHandler<dim>::active_cell_iterator neighbor = cell->neighbor(face);
+             neighbor->get_dof_indices (dofs_on_other_cell);
+                                              // only add one direction
+                                              // The other is taken care of
+                                              // by neighbor.
+             for (unsigned int i=0; i<total_dofs; ++i)
+               {
+                 for (unsigned int j=0; j<total_dofs; ++j)
+                   {
+                     sparsity.add (dofs_on_this_cell[i],
+                                   dofs_on_other_cell[j]);
+                   }
+               }
+           }
+       } 
+    }
+}
+
+#endif
 
 template<int dim>
 void
@@ -79,6 +147,10 @@ DoFTools::extract_level_dofs(const unsigned int       level,
 
 
 // explicit instantiations
+#if deal_II_dimension > 1
+template void DoFTools::make_flux_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
+                                                   SparseMatrixStruct &sparsity);
+#endif
 template void DoFTools::extract_dofs(const DoFHandler<deal_II_dimension>& dof,
                                     const vector<bool>& local_select,
                                     vector<bool>& selected_dofs);

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.