From: Daniel Arndt Date: Sun, 20 May 2018 11:33:44 +0000 (+0200) Subject: Simplify DoFTools::make_flux_sparsity_pattern and test more cases X-Git-Tag: v9.1.0-rc1~1124^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F6631%2Fhead;p=dealii.git Simplify DoFTools::make_flux_sparsity_pattern and test more cases --- diff --git a/source/dofs/dof_tools_sparsity.cc b/source/dofs/dof_tools_sparsity.cc index d4f6f9babc..c7e392b1eb 100644 --- a/source/dofs/dof_tools_sparsity.cc +++ b/source/dofs/dof_tools_sparsity.cc @@ -983,21 +983,36 @@ namespace DoFTools std::vector dofs_on_this_cell(DoFTools::max_dofs_per_cell(dof)); std::vector dofs_on_other_cell(DoFTools::max_dofs_per_cell(dof)); - std::vector > int_dof_mask (fe.size()); - - int_dof_mask = dof_couplings_from_component_couplings (fe, int_mask); - - // Convert the int_dof_mask to bool_int_dof_mask so we can pass it - // to constraints.add_entries_local_to_global() - std::vector > bool_int_dof_mask (fe.size()); + const unsigned int n_components = fe.n_components(); + AssertDimension (int_mask.size(0), n_components); + AssertDimension (int_mask.size(1), n_components); + AssertDimension (flux_mask.size(0), n_components); + AssertDimension (flux_mask.size(1), n_components); + + // note that we also need to set the respective entries if flux_mask + // says so. this is necessary since we need to consider all degrees of + // freedom on a cell for interior faces. + Table<2,Coupling> int_and_flux_mask (n_components, n_components); + for (unsigned int c1=0; c1 > int_and_flux_dof_mask + = dof_couplings_from_component_couplings (fe, int_and_flux_mask); + + // Convert int_and_flux_dof_mask to bool_int_and_flux_dof_mask so we + // can pass it to constraints.add_entries_local_to_global() + std::vector > bool_int_and_flux_dof_mask (fe.size()); for (unsigned int f=0; f(fe[f].dofs_per_cell,fe[f].dofs_per_cell)); - bool_int_dof_mask[f].fill (false); + bool_int_and_flux_dof_mask[f].reinit (TableIndices<2>(fe[f].dofs_per_cell, + fe[f].dofs_per_cell)); + bool_int_and_flux_dof_mask[f].fill (false); for (unsigned int i=0; iget_fe().dofs_per_cell); cell->get_dof_indices (dofs_on_this_cell); - // make sparsity pattern for this cell + // make sparsity pattern for this cell also taking into account + // the couplings due to face contributions on the same cell constraints.add_entries_local_to_global (dofs_on_this_cell, sparsity, keep_constrained_dofs, - bool_int_dof_mask[cell->active_fe_index()]); - // Loop over all interior neighbors + bool_int_and_flux_dof_mask[cell->active_fe_index()]); + + // Loop over interior faces for (unsigned int face = 0; face < GeometryInfo::faces_per_cell; ++face) @@ -1029,35 +1046,6 @@ namespace DoFTools const bool periodic_neighbor = cell->has_periodic_neighbor (face); - for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) - { - const unsigned int ii - = (cell->get_fe().is_primitive(i) ? - cell->get_fe().system_to_component_index(i).first - : - cell->get_fe().get_nonzero_components(i).first_selected_component() - ); - - Assert (ii < cell->get_fe().n_components(), ExcInternalError()); - - for (unsigned int j=0; jget_fe().dofs_per_cell; ++j) - { - const unsigned int jj - = (cell->get_fe().is_primitive(j) ? - cell->get_fe().system_to_component_index(j).first - : - cell->get_fe().get_nonzero_components(j).first_selected_component() - ); - - Assert (jj < cell->get_fe().n_components(), ExcInternalError()); - - if ((flux_mask(ii,jj) == always) - || - (flux_mask(ii,jj) == nonzero)) - sparsity.add (dofs_on_this_cell[i], - dofs_on_this_cell[j]); - } - } if ((!cell->at_boundary(face)) || periodic_neighbor) { typename dealii::hp::DoFHandler::level_cell_iterator diff --git a/tests/hp/flux_sparsity_02.cc b/tests/hp/flux_sparsity_02.cc index f9757fac9f..aab2d77f90 100644 --- a/tests/hp/flux_sparsity_02.cc +++ b/tests/hp/flux_sparsity_02.cc @@ -73,14 +73,75 @@ check () dof_handler.begin_active()->set_active_fe_index(1); dof_handler.distribute_dofs(fe_collection); - DynamicSparsityPattern dsp(dof_handler.n_dofs(), - dof_handler.n_dofs()); - - DoFTools::make_flux_sparsity_pattern (dof_handler, dsp); - dsp.compress(); - - // Print sparsity pattern, we expect that all dofs couple with each other. - dsp.print(deallog.get_file_stream()); + deallog << "Dimension " << dim << std::endl; + + { + deallog << "cell and face coupling" << std::endl; + DynamicSparsityPattern dsp(dof_handler.n_dofs(), + dof_handler.n_dofs()); + + Table<2, DoFTools::Coupling> cell_coupling(1,1); + Table<2, DoFTools::Coupling> face_coupling(1,1); + cell_coupling(0,0) = DoFTools::always; + face_coupling(0,0) = DoFTools::always; + + DoFTools::make_flux_sparsity_pattern (dof_handler, dsp, cell_coupling, face_coupling); + dsp.compress(); + + // Print sparsity pattern, we expect that all dofs couple with each other. + dsp.print(deallog.get_file_stream()); + } + + { + deallog << "cell coupling" << std::endl; + DynamicSparsityPattern dsp(dof_handler.n_dofs(), + dof_handler.n_dofs()); + + Table<2, DoFTools::Coupling> cell_coupling(1,1); + Table<2, DoFTools::Coupling> face_coupling(1,1); + cell_coupling(0,0) = DoFTools::always; + face_coupling(0,0) = DoFTools::none; + + DoFTools::make_flux_sparsity_pattern (dof_handler, dsp, cell_coupling, face_coupling); + dsp.compress(); + + // Print sparsity pattern, we expect no couplings across the face. + dsp.print(deallog.get_file_stream()); + } + + { + deallog << "face coupling" << std::endl; + DynamicSparsityPattern dsp(dof_handler.n_dofs(), + dof_handler.n_dofs()); + + Table<2, DoFTools::Coupling> cell_coupling(1,1); + Table<2, DoFTools::Coupling> face_coupling(1,1); + cell_coupling(0,0) = DoFTools::none; + face_coupling(0,0) = DoFTools::always; + + DoFTools::make_flux_sparsity_pattern (dof_handler, dsp, cell_coupling, face_coupling); + dsp.compress(); + + // Print sparsity pattern, we expect that all dofs couple with each other. + dsp.print(deallog.get_file_stream()); + } + + { + deallog << "no coupling" << std::endl; + DynamicSparsityPattern dsp(dof_handler.n_dofs(), + dof_handler.n_dofs()); + + Table<2, DoFTools::Coupling> cell_coupling(1,1); + Table<2, DoFTools::Coupling> face_coupling(1,1); + cell_coupling(0,0) = DoFTools::none; + face_coupling(0,0) = DoFTools::none; + + DoFTools::make_flux_sparsity_pattern (dof_handler, dsp, cell_coupling, face_coupling); + dsp.compress(); + + // Print sparsity pattern, we expect no couplings. + dsp.print(deallog.get_file_stream()); + } } diff --git a/tests/hp/flux_sparsity_02.output b/tests/hp/flux_sparsity_02.output index 26fea0c985..98b1f3efa4 100644 --- a/tests/hp/flux_sparsity_02.output +++ b/tests/hp/flux_sparsity_02.output @@ -1,9 +1,51 @@ +DEAL:2d::Dimension 2 +DEAL:2d::cell and face coupling [0,0,1,2,3,4] [1,0,1,2,3,4] [2,0,1,2,3,4] [3,0,1,2,3,4] [4,0,1,2,3,4] +DEAL:2d::cell coupling +[0,0,1,2,3] +[1,0,1,2,3] +[2,0,1,2,3] +[3,0,1,2,3] +[4,4] +DEAL:2d::face coupling +[0,0,1,2,3,4] +[1,0,1,2,3,4] +[2,0,1,2,3,4] +[3,0,1,2,3,4] +[4,0,1,2,3,4] +DEAL:2d::no coupling +[0] +[1] +[2] +[3] +[4] +DEAL:3d::Dimension 3 +DEAL:3d::cell and face coupling +[0,0,1,2,3,4,5,6,7,8] +[1,0,1,2,3,4,5,6,7,8] +[2,0,1,2,3,4,5,6,7,8] +[3,0,1,2,3,4,5,6,7,8] +[4,0,1,2,3,4,5,6,7,8] +[5,0,1,2,3,4,5,6,7,8] +[6,0,1,2,3,4,5,6,7,8] +[7,0,1,2,3,4,5,6,7,8] +[8,0,1,2,3,4,5,6,7,8] +DEAL:3d::cell coupling +[0,0,1,2,3,4,5,6,7] +[1,0,1,2,3,4,5,6,7] +[2,0,1,2,3,4,5,6,7] +[3,0,1,2,3,4,5,6,7] +[4,0,1,2,3,4,5,6,7] +[5,0,1,2,3,4,5,6,7] +[6,0,1,2,3,4,5,6,7] +[7,0,1,2,3,4,5,6,7] +[8,8] +DEAL:3d::face coupling [0,0,1,2,3,4,5,6,7,8] [1,0,1,2,3,4,5,6,7,8] [2,0,1,2,3,4,5,6,7,8] @@ -13,3 +55,13 @@ [6,0,1,2,3,4,5,6,7,8] [7,0,1,2,3,4,5,6,7,8] [8,0,1,2,3,4,5,6,7,8] +DEAL:3d::no coupling +[0] +[1] +[2] +[3] +[4] +[5] +[6] +[7] +[8]