/**
- * Test p-multigrid for a uniformly refined mesh both for simplex and
- * hypercube mesh.
+ * Test p-multigrid for a uniformly and locally refined meshes both for
+ * simplex and hypercube mesh.
*/
#include "multigrid_util.h"
void
test(const unsigned int n_refinements,
const unsigned int fe_degree_fine,
- const bool do_simplex_mesh)
+ const bool do_simplex_mesh,
+ const unsigned int mesh_type)
{
using VectorType = LinearAlgebra::distributed::Vector<Number>;
GridGenerator::subdivided_hyper_cube_with_simplices(tria, 2);
else
GridGenerator::subdivided_hyper_cube(tria, 2);
- tria.refine_global(n_refinements);
+
+ if (mesh_type == 0)
+ {
+ tria.refine_global(n_refinements);
+ }
+ else if (mesh_type == 1)
+ {
+ for (unsigned int i = 1; i < n_refinements; i++)
+ {
+ for (auto cell : tria.active_cell_iterators())
+ if (cell->is_locally_owned())
+ {
+ bool flag = true;
+ for (int d = 0; d < dim; d++)
+ if (cell->center()[d] > 0.5)
+ flag = false;
+ if (flag)
+ cell->set_refine_flag();
+ }
+ tria.execute_coarsening_and_refinement();
+ }
+ }
+ else
+ AssertThrow(false, ExcNotImplemented());
const auto level_degrees =
MGTransferGlobalCoarseningTools::create_polynomial_coarsening_sequence(
constraint.reinit(locally_relevant_dofs);
VectorTools::interpolate_boundary_values(
*mapping, dof_handler, 0, Functions::ZeroFunction<dim>(), constraint);
+ DoFTools::make_hanging_node_constraints(dof_handler, constraint);
constraint.close();
// set up operator
operators,
transfer);
+ constraints[max_level].distribute(dst);
+
deallog << dim << " " << fe_degree_fine << " " << n_refinements << " "
<< (do_simplex_mesh ? "tri " : "quad") << " "
<< solver_control.last_step() << std::endl;
deallog.precision(8);
- for (unsigned int n_refinements = 2; n_refinements <= 4; ++n_refinements)
- for (unsigned int degree = 2; degree <= 4; ++degree)
- test<2>(n_refinements, degree, false /*quadrilateral*/);
+ for (unsigned int mesh_type = 0; mesh_type < 2; ++mesh_type)
+ {
+ for (unsigned int n_refinements = 2; n_refinements <= 4; ++n_refinements)
+ for (unsigned int degree = 2; degree <= 4; ++degree)
+ test<2>(n_refinements, degree, false /*quadrilateral*/, mesh_type);
- for (unsigned int n_refinements = 2; n_refinements <= 4; ++n_refinements)
- for (unsigned int degree = 2; degree <= 2; ++degree)
- test<2>(n_refinements, degree, true /*triangle*/);
+ for (unsigned int n_refinements = 2; n_refinements <= 4; ++n_refinements)
+ for (unsigned int degree = 2; degree <= 2; ++degree)
+ test<2>(n_refinements, degree, true /*triangle*/, mesh_type);
+ }
}