From 039896eeaaeb16532b99a2e78e2a323ace7f3ff1 Mon Sep 17 00:00:00 2001 From: guido Date: Fri, 2 Jun 2000 19:15:14 +0000 Subject: [PATCH] Reduced sparsity pattern for DG-Multigrid git-svn-id: https://svn.dealii.org/trunk@2990 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/multigrid/mg_dof_tools.h | 24 ++++ .../deal.II/source/multigrid/mg_dof_tools.cc | 107 ++++++++++++++++++ 2 files changed, 131 insertions(+) diff --git a/deal.II/deal.II/include/multigrid/mg_dof_tools.h b/deal.II/deal.II/include/multigrid/mg_dof_tools.h index d4da86b441..210f2a6e2f 100644 --- a/deal.II/deal.II/include/multigrid/mg_dof_tools.h +++ b/deal.II/deal.II/include/multigrid/mg_dof_tools.h @@ -64,6 +64,30 @@ class MGDoFTools SparsityPattern &sparsity, const unsigned int level); + /** + * 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 MGDoFHandler &dof, + SparsityPattern &sparsity, + const unsigned int level, + const FullMatrix& int_mask, + const FullMatrix& flux_mask); }; diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index 20672ab899..f151345986 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -15,6 +15,7 @@ #include #include #include +#include #include #include #include @@ -101,6 +102,105 @@ MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler &dof, } + +template +void +MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler &dof, + SparsityPattern &sparsity, + const unsigned int level, + const FullMatrix& int_mask, + const FullMatrix& flux_mask) +{ + const unsigned int n_dofs = dof.n_dofs(level); + 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); + MGDoFHandler::cell_iterator cell = dof.begin(level), + endc = dof.end(level); + + + 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_mg_dof_indices (dofs_on_this_cell); + // make sparsity pattern for this cell + for (unsigned int i=0; i::faces_per_cell; + ++face) + { + MGDoFHandler::face_iterator cell_face = cell->face(face); + if (cell_face->user_flag_set ()) + continue; + + if (! cell->at_boundary (face) ) + { + MGDoFHandler::cell_iterator neighbor = cell->neighbor(face); + + if (neighbor->level() < cell->level()) + continue; + + unsigned int neighbor_face = cell->neighbor_of_neighbor(face); + + neighbor->get_mg_dof_indices (dofs_on_other_cell); + for (unsigned int i=0; iface(neighbor_face)->set_user_flag (); + } + } + } +} + + + + // explicit instantiations template void MGDoFTools::make_sparsity_pattern (const MGDoFHandler &, SparsityPattern &, @@ -109,3 +209,10 @@ template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler 1 +template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler &, + SparsityPattern &, + const unsigned int, + const FullMatrix&, + const FullMatrix&); +#endif -- 2.39.5