]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
optimized flux sparsity
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 May 2000 02:37:10 +0000 (02:37 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 May 2000 02:37:10 +0000 (02:37 +0000)
git-svn-id: https://svn.dealii.org/trunk@2950 0785d39b-7218-0410-832d-ea1e28bc413d

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

index f3c7c88a987afc57862ded6cdc45860fb7bd0109..79baa245594038ac6126358cce42f5cf9150b0f0 100644 (file)
@@ -284,6 +284,30 @@ class DoFTools
     static void
     make_flux_sparsity_pattern (const DoFHandler<dim> &dof_handler,
                                SparsityPattern       &sparsity_pattern);
+
+                                    /**
+                                     * This function does the same as
+                                     * the other with the same name,
+                                     * but it gets two additional
+                                     * coefficient matrices. A matrix
+                                     * entry will only be generated
+                                     * for two basis functions, if
+                                     * there is a non-zero entry
+                                     * linking their associated
+                                     * components in the coefficient
+                                     * matrix.
+                                     *
+                                     * There is one matrix for
+                                     * couplings in a cell and one
+                                     * for the couplings occuring in
+                                     * fluxes.
+                                     */
+    template <int dim>
+    static void
+    make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
+                               SparsityPattern       &sparsity,
+                               const FullMatrix<double>& int_mask,
+                               const FullMatrix<double>& flux_mask);
     
                                     /**
                                      * Make up the constraints which
index 530936b80f3a2243e165e73f9d3d00e219d41dbe..d24dc56122636fab1e2572eabe35282d3ce2d78c 100644 (file)
@@ -223,10 +223,6 @@ void DoFTools::make_boundary_sparsity_pattern (const DoFHandler<dim>& dof,
 };
 
 
-//TODO: Check this function for potential of optimization. (G)
-// sometimes, a shape function might be zero on an edge, but we
-// need more information from the finite element used. This should be
-// left for future improvements.
 template<int dim>
 void
 DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
@@ -319,6 +315,129 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
 
 
 
+template<int dim>
+void
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
+                                     SparsityPattern       &sparsity,
+                                     const FullMatrix<double>& int_mask,
+                                     const FullMatrix<double>& flux_mask)
+{
+  const unsigned int n_dofs = dof.n_dofs();
+  const unsigned int n_comp = dof.get_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));
+  
+  const unsigned int total_dofs = dof.get_fe().dofs_per_cell;
+  vector<unsigned int> dofs_on_this_cell(total_dofs);
+  vector<unsigned int> dofs_on_other_cell(total_dofs);
+  DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
+                                       endc = dof.end();
+
+
+  vector<vector<bool> > int_dof_mask(total_dofs,
+                                vector<bool>(total_dofs, false));
+  vector<vector<bool> > flux_dof_mask(total_dofs,
+                                vector<bool>(total_dofs, false));
+  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;
+      }
+  
+  (const_cast<Triangulation<dim>& > (dof.get_tria())).clear_user_flags();
+  
+  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)
+         if (int_dof_mask[i][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)
+       {
+         DoFHandler<dim>::face_iterator cell_face = cell->face(face);
+         if (cell_face->user_flag_set ())
+           continue;
+
+         if (! cell->at_boundary (face) )
+           {
+             DoFHandler<dim>::cell_iterator neighbor = cell->neighbor(face);
+                                              // Refinement edges are taken care of
+                                              // by coarser cells
+             if (neighbor->level() < cell->level())
+               continue;
+
+             unsigned int neighbor_face = cell->neighbor_of_neighbor(face);
+
+             if (neighbor->has_children())
+               {
+                 for (unsigned int sub_nr = 0;
+                      sub_nr != GeometryInfo<dim>::subfaces_per_face;
+                      ++sub_nr)
+                   {
+                     DoFHandler<dim>::cell_iterator sub_neighbor
+                       = neighbor->
+                       child(GeometryInfo<dim>::child_cell_on_face(neighbor_face, sub_nr));
+
+                     sub_neighbor->get_dof_indices (dofs_on_other_cell);
+                     for (unsigned int i=0; i<total_dofs; ++i)
+                       {
+                         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]);
+                             }
+                       }
+                     sub_neighbor->face(neighbor_face)->set_user_flag ();
+                   }
+               } else {
+                 neighbor->get_dof_indices (dofs_on_other_cell);
+                 for (unsigned int i=0; i<total_dofs; ++i)
+                   {
+                     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]);
+                       }
+                   }
+                 neighbor->face(neighbor_face)->set_user_flag (); 
+               }
+           } 
+       }
+    }
+}
+
+
+
 #if deal_II_dimension == 1
 
 template <>
@@ -1362,6 +1481,11 @@ template void
 DoFTools::make_flux_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
                                      SparsityPattern    &sparsity);
 template void
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
+                                     SparsityPattern    &,
+                                     const FullMatrix<double>&,
+                                     const FullMatrix<double>&);
+template void
 DoFTools::make_boundary_sparsity_pattern (const DoFHandler<deal_II_dimension>& dof,
                                          const vector<unsigned int>  &,
                                          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.