From: Martin Kronbichler Date: Wed, 27 Jul 2016 11:16:17 +0000 (+0200) Subject: Fix mg dof indices on objects not in standard orientation. X-Git-Tag: v8.5.0-rc1~837^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ecff93bc6e9d1f83981541d9460f4f96e0da3140;p=dealii.git Fix mg dof indices on objects not in standard orientation. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 923bcbf99f..af59acd918 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -171,7 +171,7 @@ inconvenience this causes. (Jonathan Robey, 2016/07/21) -
  • New: Added GridGenerator::quarter_hyper_ball() to generate the +
  • New: Added GridGenerator::quarter_hyper_ball() to generate the intersection of a hyper ball with the positive orthant relative to its center.
    @@ -388,7 +388,14 @@ inconvenience this causes.

    Specific improvements

      -
    1. Improved: Allow for including dofs for individual components on +
    2. Fixed: Level indices for geometric multigrid queried through + DoFAccessor::get_mg_dof_indices() would return wrong indices on lines + and faces in non-standard orientation in 3D. This is now fixed. +
      + (Martin Kronbichler, 2016/07/27) +
    3. + +
    4. Improved: Allow for including dofs for individual components on boundary in DoFTools::make_vertex_patches().
      (Ryan Grove, Daniel Arndt, 2016/07/21) diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index 21f279c8a0..be5d11813c 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -1460,11 +1460,18 @@ namespace internal for (unsigned int line = 0; line < GeometryInfo<3>::lines_per_cell; ++line) for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof) - accessor.line(line)->set_mg_dof_index(level, dof, *next++); + accessor.line (line)->set_mg_dof_index + (level, fe.adjust_line_dof_index_for_line_orientation(dof,accessor.line_orientation(line)), *next++); for (unsigned int quad = 0; quad < GeometryInfo<3>::quads_per_cell; ++quad) for (unsigned int dof = 0; dof < fe.dofs_per_quad; ++dof) - accessor.quad(quad)->set_mg_dof_index(level, dof, *next++); + accessor.quad (quad)->set_mg_dof_index + (level, fe.adjust_quad_dof_index_for_face_orientation + (dof, + accessor.face_orientation(quad), + accessor.face_flip(quad), + accessor.face_rotation(quad)), + *next++); for (unsigned int dof = 0; dof < fe.dofs_per_hex; ++dof) accessor.set_mg_dof_index(level, dof, *next++); @@ -1883,11 +1890,17 @@ namespace internal for (unsigned int line = 0; line < GeometryInfo<3>::lines_per_cell; ++line) for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof) - *next++ = accessor.line (line)->mg_dof_index (level, dof); + *next++ = accessor.line (line)->mg_dof_index + (level, accessor.get_fe(fe_index).adjust_line_dof_index_for_line_orientation(dof,accessor.line_orientation(line))); for (unsigned int quad = 0; quad < GeometryInfo<3>::quads_per_cell; ++quad) for (unsigned int dof = 0; dof < fe.dofs_per_quad; ++dof) - *next++ = accessor.quad (quad)->mg_dof_index (level, dof); + *next++ = accessor.quad (quad)->mg_dof_index + (level, accessor.get_fe(fe_index).adjust_quad_dof_index_for_face_orientation + (dof, + accessor.face_orientation(quad), + accessor.face_flip(quad), + accessor.face_rotation(quad))); for (unsigned int dof = 0; dof < fe.dofs_per_hex; ++dof) *next++ = accessor.mg_dof_index (level, dof); diff --git a/tests/multigrid/dof_04.cc b/tests/multigrid/dof_04.cc new file mode 100644 index 0000000000..5feefaa0b4 --- /dev/null +++ b/tests/multigrid/dof_04.cc @@ -0,0 +1,76 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// check DoFAccessor::get_mg_dof_indices on meshes that are not in standard +// orientation by comparing to DoFAccessor::get_dof_indices + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + + +template +void check() +{ + // need cubic polynomials that have two dofs on lines + FE_Q fe(3); + + Triangulation tr; + if (dim > 1) + GridGenerator::hyper_shell(tr, Point(), 0.5, 1, 12); + else + GridGenerator::hyper_cube(tr); + + DoFHandler dof(tr); + dof.distribute_dofs(fe); + dof.distribute_mg_dofs(fe); + + std::vector dof_indices(fe.dofs_per_cell); + std::vector mg_dof_indices(fe.dofs_per_cell); + for (typename DoFHandler::active_cell_iterator cell=dof.begin_active(); + cell != dof.end(); ++cell) + { + cell->get_dof_indices(dof_indices); + cell->get_mg_dof_indices(mg_dof_indices); + bool has_error = false; + // dof indices should have the same order on both the mg dofs and the + // usual dofs becuse there is only one level + for (unsigned int i=0; icenter() << ": "; + for (unsigned int i=0; i (); + check<2> (); + check<3> (); +} diff --git a/tests/multigrid/dof_04.output b/tests/multigrid/dof_04.output new file mode 100644 index 0000000000..8b58a4b83f --- /dev/null +++ b/tests/multigrid/dof_04.output @@ -0,0 +1,4 @@ + +DEAL::1D OK +DEAL::2D OK +DEAL::3D OK diff --git a/tests/multigrid/dof_05.cc b/tests/multigrid/dof_05.cc new file mode 100644 index 0000000000..86e2a21778 --- /dev/null +++ b/tests/multigrid/dof_05.cc @@ -0,0 +1,52 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// check that distributed_mg_dofs works correctly on non-standard oriented +// meshes when used in parallel + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + + +template +void check() +{ + // need cubic polynomials that have two dofs on lines + FE_Q fe(3); + + parallel::distributed::Triangulation tr(MPI_COMM_WORLD, Triangulation::none, + parallel::distributed::Triangulation::construct_multigrid_hierarchy); + GridGenerator::hyper_shell(tr, Point(), 0.5, 1, 12); + tr.refine_global(1); + + DoFHandler dof(tr); + dof.distribute_dofs(fe); + dof.distribute_mg_dofs(fe); + deallog << dim << "D OK" << std::endl; +} + +int main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi(argc, argv); + mpi_initlog(); + check<2> (); + check<3> (); +} diff --git a/tests/multigrid/dof_05.with_mpi=true.with_p4est=true.mpirun=3.output b/tests/multigrid/dof_05.with_mpi=true.with_p4est=true.mpirun=3.output new file mode 100644 index 0000000000..ad8782de68 --- /dev/null +++ b/tests/multigrid/dof_05.with_mpi=true.with_p4est=true.mpirun=3.output @@ -0,0 +1,3 @@ + +DEAL::2D OK +DEAL::3D OK