From 6ba5e2cf869b10f589da36b3e9c389a7eaafae68 Mon Sep 17 00:00:00 2001 From: kayser-herold Date: Sat, 29 Jul 2006 23:18:13 +0000 Subject: [PATCH] Added a testcase for the hanging_node_constraints in the hp-case and modified test step-3c to work properly in the hp-setting. git-svn-id: https://svn.dealii.org/trunk@13511 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/hp/Makefile | 3 +- tests/hp/hp_hanging_nodes_01.cc | 178 ++++++++++++++++++++++++++++++++ tests/hp/step-3c.cc | 9 +- 3 files changed, 186 insertions(+), 4 deletions(-) create mode 100644 tests/hp/hp_hanging_nodes_01.cc diff --git a/tests/hp/Makefile b/tests/hp/Makefile index f96dfe52f4..908f7e6709 100644 --- a/tests/hp/Makefile +++ b/tests/hp/Makefile @@ -32,7 +32,8 @@ tests_x = hp_dof_handler \ step-* \ n_active_fe_indices \ integrate_difference \ - hp_*_dof_identities_* + hp_*_dof_identities_* \ + hp_hanging_nodes_* # from above list of regular expressions, generate the real set of # tests diff --git a/tests/hp/hp_hanging_nodes_01.cc b/tests/hp/hp_hanging_nodes_01.cc new file mode 100644 index 0000000000..369860485b --- /dev/null +++ b/tests/hp/hp_hanging_nodes_01.cc @@ -0,0 +1,178 @@ +/* $Id$ */ +/* Version: $Name$ */ +/* */ +/* Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2006 by the deal.II authors */ +/* */ +/* This file is subject to QPL and may not be distributed */ +/* without copyright and license information. Please refer */ +/* to the file deal.II/doc/license.html for the text and */ +/* further information on this license. */ + + +#include +#include +#include +#include +#include + +#include +#include +#include + +#include + +#include +#include + +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include +#include +#include + + +// This a test for the hp capable version of the make_hanging_node_constraints +// method. It uses a triangulation with one refined element beside an +// unrefined element to create the constraints for this configuration. + +template +int generate_grid (Triangulation &tria) +{ + Point p1, + p2; + std::vector sub_div; + + // Define a rectangular shape + for (unsigned int d=0; d < dim; ++d) + { + p1(d) = 0; + p2(d) = (d == 0) ? 2.0 : 1.0; + sub_div.push_back ( (d == 0) ? 2 : 1); + } + GridGenerator::subdivided_hyper_rectangle (tria, sub_div, p1, p2, true); + + // Refine the first cell. + tria.begin_active ()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + + return (0); +} + + + +template +void test_constraints (hp::FECollection &fe_coll) +{ + Triangulation tria; + + // Setup a rectangular domain + // where one cell is h-refined, + // while the other cell is + // unrefined. Furthermore every cell + // gets a different active_fe_index. + // This should serve as a testcase + // for the hanging node constraints. + generate_grid (tria); + + // Now assign increasing + // active_fe_indices to + // the different cells. + hp::DoFHandler dof_handler (tria); + typename hp::DoFHandler::active_cell_iterator cell = dof_handler.begin_active (), + endc = dof_handler.end (); + unsigned int fe_indx = 0; + for (; cell != endc; ++cell) + { + cell->set_active_fe_index (fe_indx); + ++fe_indx; + } + + // Distribute DoFs; + dof_handler.distribute_dofs (fe_coll); + deallog << "DoFs: " << dof_handler.n_dofs () << std::endl; + + // Create the constraints. + ConstraintMatrix constraint_matrix; + + DoFTools::make_hanging_node_constraints (dof_handler, + constraint_matrix); + + // Output the constraints + constraint_matrix.print (deallog.get_file_stream ()); +} + + +template +void test_constraints_old (FiniteElement &fe) +{ + Triangulation tria; + + // Setup a rectangular domain + // where one cell is h-refined, + // while the other cell is + // unrefined. Furthermore every cell + // gets a different active_fe_index. + // This should serve as a testcase + // for the hanging node constraints. + generate_grid (tria); + + // Now assign increasing + // active_fe_indices to + // the different cells. + DoFHandler dof_handler (tria); + typename DoFHandler::active_cell_iterator cell = dof_handler.begin_active (), + endc = dof_handler.end (); + + // Distribute DoFs; + dof_handler.distribute_dofs (fe); + deallog << "DoFs: " << dof_handler.n_dofs () << std::endl; + + // Create the constraints. + ConstraintMatrix constraint_matrix; + + DoFTools::make_hanging_node_constraints (dof_handler, + constraint_matrix); + + // Output the constraints + constraint_matrix.print (deallog.get_file_stream ()); +} + +int main () +{ + std::ofstream logfile("hp_hanging_nodes_01/output"); + logfile.precision(2); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + FE_Q<2> fe_1 (1); + FE_Q<2> fe_2 (2); + FE_Q<2> fe_3 (3); + + hp::FECollection<2> fe_coll2; + fe_coll2.push_back (fe_3); + fe_coll2.push_back (fe_2); + fe_coll2.push_back (fe_2); + + fe_coll2.push_back (fe_2); + fe_coll2.push_back (fe_3); + + test_constraints<2> (fe_coll2); + + test_constraints_old<2> (fe_1); + + return 0; +} diff --git a/tests/hp/step-3c.cc b/tests/hp/step-3c.cc index 96cc6bb583..b2737bffc2 100644 --- a/tests/hp/step-3c.cc +++ b/tests/hp/step-3c.cc @@ -144,12 +144,15 @@ void LaplaceProblem::make_grid_and_dofs () void LaplaceProblem::assemble_system () { - hp::QCollection<2> quadrature_formula(QGauss<2>(6)); + hp::QCollection<2> quadrature_formula; + + for (unsigned int p = 0; p < fe.size (); ++p) + quadrature_formula.push_back (QGauss<2> (p + 2)); + hp::FEValues<2> x_fe_values (fe, quadrature_formula, update_values | update_gradients | update_JxW_values); const unsigned int max_dofs_per_cell = fe.max_dofs_per_cell (); - const unsigned int n_q_points = quadrature_formula[0].n_quadrature_points; FullMatrix cell_matrix (max_dofs_per_cell, max_dofs_per_cell); Vector cell_rhs (max_dofs_per_cell); @@ -169,7 +172,7 @@ void LaplaceProblem::assemble_system () cell_rhs = 0; const unsigned int dofs_per_cell = cell->get_fe ().dofs_per_cell; - + const unsigned int n_q_points = quadrature_formula[cell->active_fe_index ()].n_quadrature_points; for (unsigned int i=0; i