From 8a7ecf83f881e86bd587423efaed022ac1c1ab7d Mon Sep 17 00:00:00 2001 From: heister Date: Wed, 24 Aug 2011 20:39:36 +0000 Subject: [PATCH] add failing no_flux test (need to fix DoFRenumbering::hierarchical() first git-svn-id: https://svn.dealii.org/trunk@24188 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/mpi/no_flux_constraints_02.cc | 226 ++++++++++++++++++ .../no_flux_constraints_02/ncpu_1/cmp/generic | 0 .../no_flux_constraints_02/ncpu_9/cmp/generic | 14 ++ 3 files changed, 240 insertions(+) create mode 100644 tests/mpi/no_flux_constraints_02.cc create mode 100644 tests/mpi/no_flux_constraints_02/ncpu_1/cmp/generic create mode 100644 tests/mpi/no_flux_constraints_02/ncpu_9/cmp/generic diff --git a/tests/mpi/no_flux_constraints_02.cc b/tests/mpi/no_flux_constraints_02.cc new file mode 100644 index 0000000000..219d903588 --- /dev/null +++ b/tests/mpi/no_flux_constraints_02.cc @@ -0,0 +1,226 @@ +//--------------------------------------------------------------------------- +// $Id: no_flux_constraints.cc 24175 2011-08-24 12:14:21Z kronbichler $ +// Version: $Name$ +// +// Copyright (C) 2009, 2010 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. +// +//--------------------------------------------------------------------------- + + +// check ConstraintMatrix for a distributed mesh on hyper shell with both +// hanging nodes and no-normal-flux constraints on the outer boundary of the +// shell. +// fails because DoFRenumbering::hierarchical does not work correctly with +// p4est right now. Also see test renumber_z_order_02 + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +template +void test() +{ + Assert (dim == 3, ExcNotImplemented()); + unsigned int myid = Utilities::System::get_this_mpi_process (MPI_COMM_WORLD); + unsigned int numprocs = Utilities::System::get_n_mpi_processes (MPI_COMM_WORLD); + + parallel::distributed::Triangulation triangulation(MPI_COMM_WORLD); + + // create hypershell mesh and refine some + // cells. use some large numbers which make + // round-off influence more pronounced + const double R0 = 6371000.-2890000.; /* m */ + const double R1 = 6371000.- 35000.; /* m */ + GridGenerator::hyper_shell (triangulation, + Point(), + R0, + R1, + 96, true); + triangulation.refine_global (2); + + for (typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + cell != triangulation.end(); ++cell) + if (!cell->is_ghost() && !cell->is_artificial()) + if (cell->center()[2] > 0.75 * R1) + { + cell->set_refine_flag(); + for (unsigned int c=0; c<8; ++c) + cell->parent()->child(c)->set_refine_flag(); + } + + { + unsigned int n_flagged_cells = 0; + for (typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + cell != triangulation.end(); ++cell) + if (!cell->is_ghost() && !cell->is_artificial()) + if (cell->refine_flag_set()) + ++n_flagged_cells; + + unsigned int global_f_c = 0; + MPI_Allreduce (&n_flagged_cells, &global_f_c, 1, MPI_UNSIGNED, + MPI_SUM, MPI_COMM_WORLD); + if (myid == 0) + deallog << "# flagged cells = " << global_f_c << std::endl; + } + + triangulation.prepare_coarsening_and_refinement(); + triangulation.execute_coarsening_and_refinement (); + + if (myid == 0) + deallog << "#cells = " << triangulation.n_global_active_cells() + << std::endl; + + // create FE_System and fill in no-normal flux + // conditions on boundary 1 (outer) + static const FESystem fe(FE_Q (1), dim); + DoFHandler dofh(triangulation); + dofh.distribute_dofs (fe); + DoFRenumbering::hierarchical(dofh); + + if (myid == 0) + deallog << "#dofs = " << dofh.locally_owned_dofs().size() + << std::endl; + + IndexSet relevant_set; + DoFTools::extract_locally_relevant_dofs (dofh, relevant_set); + + ConstraintMatrix constraints; + constraints.reinit(relevant_set); + DoFTools::make_hanging_node_constraints (dofh, constraints); + std::set no_normal_flux_boundaries; + no_normal_flux_boundaries.insert (1); + const unsigned int degree = 1; +/* VectorTools::compute_no_normal_flux_constraints (dofh, 0, + no_normal_flux_boundaries, + constraints, + MappingQ(degree));*/ + constraints.close(); + + if (myid==0) + system("rm -rf no_flux_constraints_02/cm_?.dot"); + + MPI_Barrier(MPI_COMM_WORLD); + + { //write the constraintmatrix to a file on each cpu + char fname[] = "no_flux_constraints_02/cm_0.dot"; + fname[26]+=myid; + std::ofstream file(fname); + constraints.print(file); + } + MPI_Barrier(MPI_COMM_WORLD); + sleep(1); + if (myid==0) + { + //sort and merge the constraint matrices on proc 0, generate a checksum + //and output that into the deallog + system("cat no_flux_constraints_02/cm_?.dot|sort -n|uniq >no_flux_constraints_02/cm"); + system("md5sum no_flux_constraints_02/cm >no_flux_constraints_02/cm.check"); + { + std::ifstream file("no_flux_constraints_02/cm.check"); + std::string str; + while (!file.eof()) + { + std::getline(file, str); + deallog << str << std::endl; + } + } + } + // print the number of constraints. since + // processors might write info in different + // orders, copy all numbers to root processor + std::vector n_constraints_glob (numprocs); + unsigned int n_constraints = constraints.n_constraints(); + MPI_Gather (&n_constraints, 1, MPI_UNSIGNED, + &n_constraints_glob[0], 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD); + if (myid == 0) + for (unsigned int i=0; i local_dof_indices (dofs_per_cell); + Vector local_vector (dofs_per_cell); + for (unsigned int i=0; i::active_cell_iterator + cell = dofh.begin_active(), + endc = dofh.end(); + for (; cell!=endc; ++cell) + if (cell->subdomain_id() == triangulation.locally_owned_subdomain()) + { + cell->get_dof_indices (local_dof_indices); + constraints.distribute_local_to_global (local_vector, + local_dof_indices, + vector); + } + vector.compress (Add); + } + + // now check that no entries were generated + // for constrained entries on the locally + // owned range. + const std::pair range = vector.local_range(); + for (unsigned int i=range.first; i(); + deallog.pop(); + } + else + test<3>(); + } + + return 0; +} diff --git a/tests/mpi/no_flux_constraints_02/ncpu_1/cmp/generic b/tests/mpi/no_flux_constraints_02/ncpu_1/cmp/generic new file mode 100644 index 0000000000..e69de29bb2 diff --git a/tests/mpi/no_flux_constraints_02/ncpu_9/cmp/generic b/tests/mpi/no_flux_constraints_02/ncpu_9/cmp/generic new file mode 100644 index 0000000000..ba60b2ea3d --- /dev/null +++ b/tests/mpi/no_flux_constraints_02/ncpu_9/cmp/generic @@ -0,0 +1,14 @@ + +DEAL:0:3d::# flagged cells = 224 +DEAL:0:3d::#cells = 7712 +DEAL:0:3d::#dofs = 27018 +DEAL:0:3d::#constraints on 0: 168 +DEAL:0:3d::#constraints on 1: 1017 +DEAL:0:3d::#constraints on 2: 617 +DEAL:0:3d::#constraints on 3: 1096 +DEAL:0:3d::#constraints on 4: 626 +DEAL:0:3d::#constraints on 5: 1032 +DEAL:0:3d::#constraints on 6: 695 +DEAL:0:3d::#constraints on 7: 907 +DEAL:0:3d::#constraints on 8: 308 +DEAL:0:3d::OK -- 2.39.5