From: Guido Kanschat Date: Fri, 26 May 2000 02:37:10 +0000 (+0000) Subject: optimized flux sparsity X-Git-Tag: v8.0.0~20456 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9e24b41391a30b88146e5b531edb9991f7eb9915;p=dealii.git optimized flux sparsity git-svn-id: https://svn.dealii.org/trunk@2950 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index f3c7c88a98..79baa24559 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -284,6 +284,30 @@ class DoFTools static void make_flux_sparsity_pattern (const DoFHandler &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 + static void + make_flux_sparsity_pattern (const DoFHandler &dof, + SparsityPattern &sparsity, + const FullMatrix& int_mask, + const FullMatrix& flux_mask); /** * Make up the constraints which diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 530936b80f..d24dc56122 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -223,10 +223,6 @@ void DoFTools::make_boundary_sparsity_pattern (const DoFHandler& 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 void DoFTools::make_flux_sparsity_pattern (const DoFHandler &dof, @@ -319,6 +315,129 @@ DoFTools::make_flux_sparsity_pattern (const DoFHandler &dof, +template +void +DoFTools::make_flux_sparsity_pattern (const DoFHandler &dof, + SparsityPattern &sparsity, + const FullMatrix& int_mask, + const FullMatrix& 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 dofs_on_this_cell(total_dofs); + vector dofs_on_other_cell(total_dofs); + DoFHandler::active_cell_iterator cell = dof.begin_active(), + endc = dof.end(); + + + vector > int_dof_mask(total_dofs, + vector(total_dofs, false)); + vector > flux_dof_mask(total_dofs, + vector(total_dofs, false)); + for (unsigned int i=0; i& > (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::faces_per_cell; + ++face) + { + DoFHandler::face_iterator cell_face = cell->face(face); + if (cell_face->user_flag_set ()) + continue; + + if (! cell->at_boundary (face) ) + { + DoFHandler::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::subfaces_per_face; + ++sub_nr) + { + DoFHandler::cell_iterator sub_neighbor + = neighbor-> + child(GeometryInfo::child_cell_on_face(neighbor_face, sub_nr)); + + sub_neighbor->get_dof_indices (dofs_on_other_cell); + for (unsigned int i=0; iface(neighbor_face)->set_user_flag (); + } + } else { + neighbor->get_dof_indices (dofs_on_other_cell); + for (unsigned int i=0; iface(neighbor_face)->set_user_flag (); + } + } + } + } +} + + + #if deal_II_dimension == 1 template <> @@ -1362,6 +1481,11 @@ template void DoFTools::make_flux_sparsity_pattern (const DoFHandler& dof, SparsityPattern &sparsity); template void +DoFTools::make_flux_sparsity_pattern (const DoFHandler& dof, + SparsityPattern &, + const FullMatrix&, + const FullMatrix&); +template void DoFTools::make_boundary_sparsity_pattern (const DoFHandler& dof, const vector &, SparsityPattern &);