From: bangerth Date: Mon, 7 Apr 2014 20:42:02 +0000 (+0000) Subject: Add a test that currently fails. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1d7b03d889ab8818c7643a4fd569e44697029707;p=dealii-svn.git Add a test that currently fails. git-svn-id: https://svn.dealii.org/trunk@32729 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/mpi/torus.cc b/tests/mpi/torus.cc new file mode 100644 index 0000000000..e61edf5dea --- /dev/null +++ b/tests/mpi/torus.cc @@ -0,0 +1,118 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2014 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. +// +// --------------------------------------------------------------------- + +/* + * Author: Guido Kanschat, Texas A&M University, 2009 + */ + +#include "../tests.h" +#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 +#include +#include + +#include +#include +#include +#include + +#include +#include + + +typedef parallel::distributed::Triangulation<2,3>::cell_iterator cell_iterator; +DeclException1(ExcMissingCell, + cell_iterator, + << "Trying to find cell " << arg1 << " but it doesn't appear to be in the list"); + +int main(int argc, char *argv[]) +{ + dealii::Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv); + MPILogInitAll log; + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + { + static std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + } + + parallel::distributed::Triangulation<2,3> triangulation(MPI_COMM_WORLD, + typename Triangulation<2,3>::MeshSmoothing + (Triangulation<2,3>::smoothing_on_refinement | + Triangulation<2,3 >::smoothing_on_coarsening)); + GridGenerator::torus(triangulation, 1, 0.2); + + // create a set of all cells, and insert all cells into it + std::set::cell_iterator> cells; + for (parallel::distributed::Triangulation<2,3>::cell_iterator cell= triangulation.begin(0); + cell!=triangulation.end(0); + ++cell) + { + cells.insert (cell); + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "Adding cell " << cell << std::endl; + } + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "List contains " << cells.size() << " items" << std::endl; + + // verify that every cell is in there + for(parallel::distributed::Triangulation<2,3>::cell_iterator cell= triangulation.begin(0); + cell!=triangulation.end(0); + ++cell) + Assert (cells.find(cell)!=cells.end(), + ExcMissingCell(cell)); + + // refine triangulation and verify that every coarse mesh cell is in there + triangulation.refine_global(2); + + if (Utilities::MPI::this_mpi_process (MPI_COMM_WORLD) == 0) + deallog << "List contains " << cells.size() << " items" << std::endl; + for(parallel::distributed::Triangulation<2,3>::cell_iterator cell= triangulation.begin(0); + cell!=triangulation.end(0); + ++cell) + Assert (cells.find(cell)!=cells.end(), + ExcMissingCell(cell)); +} diff --git a/tests/mpi/torus.mpirun=2.output b/tests/mpi/torus.mpirun=2.output new file mode 100644 index 0000000000..decf2f7cfd --- /dev/null +++ b/tests/mpi/torus.mpirun=2.output @@ -0,0 +1,125 @@ + +DEAL:0::Element: FE_DGQ<2>(2) +DEAL:0::Step 0 +DEAL:0::Triangulation 12 cells, 2 levels +DEAL:0::DoFHandler 108 dofs, level dofs 27 108 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 20.6576 +DEAL:0:cg::Convergence step 15 value 8.11943e-13 +DEAL:0::Estimate 1.10045 +DEAL:0::Step 1 +DEAL:0::Triangulation 18 cells, 3 levels +DEAL:0::DoFHandler 162 dofs, level dofs 27 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.1045 +DEAL:0:cg::Convergence step 16 value 8.12709e-13 +DEAL:0::Estimate 1.03046 +DEAL:0::Step 2 +DEAL:0::Triangulation 27 cells, 4 levels +DEAL:0::DoFHandler 243 dofs, level dofs 27 108 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.3244 +DEAL:0:cg::Convergence step 17 value 6.05005e-13 +DEAL:0::Estimate 0.836764 +DEAL:0::Step 3 +DEAL:0::Triangulation 39 cells, 5 levels +DEAL:0::DoFHandler 351 dofs, level dofs 27 108 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.4584 +DEAL:0:cg::Convergence step 16 value 3.55455e-13 +DEAL:0::Estimate 0.644552 +DEAL:0::Step 4 +DEAL:0::Triangulation 54 cells, 6 levels +DEAL:0::DoFHandler 486 dofs, level dofs 27 108 180 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.5437 +DEAL:0:cg::Convergence step 16 value 8.06366e-13 +DEAL:0::Estimate 0.474987 +DEAL:0::Step 5 +DEAL:0::Triangulation 84 cells, 7 levels +DEAL:0::DoFHandler 756 dofs, level dofs 27 108 324 216 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 25.5435 +DEAL:0:cg::Convergence step 17 value 8.48262e-13 +DEAL:0::Estimate 0.356044 + +DEAL:1::Element: FE_DGQ<2>(2) +DEAL:1::Step 0 +DEAL:1::Triangulation 12 cells, 2 levels +DEAL:1::DoFHandler 108 dofs, level dofs 27 108 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 20.6576 +DEAL:1:cg::Convergence step 15 value 8.11943e-13 +DEAL:1::Estimate 1.10045 +DEAL:1::Step 1 +DEAL:1::Triangulation 18 cells, 3 levels +DEAL:1::DoFHandler 162 dofs, level dofs 27 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.1045 +DEAL:1:cg::Convergence step 16 value 8.12709e-13 +DEAL:1::Estimate 1.03046 +DEAL:1::Step 2 +DEAL:1::Triangulation 27 cells, 4 levels +DEAL:1::DoFHandler 243 dofs, level dofs 27 108 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.3244 +DEAL:1:cg::Convergence step 17 value 6.05005e-13 +DEAL:1::Estimate 0.836764 +DEAL:1::Step 3 +DEAL:1::Triangulation 39 cells, 5 levels +DEAL:1::DoFHandler 351 dofs, level dofs 27 108 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.4584 +DEAL:1:cg::Convergence step 16 value 3.55455e-13 +DEAL:1::Estimate 0.644552 +DEAL:1::Step 4 +DEAL:1::Triangulation 54 cells, 6 levels +DEAL:1::DoFHandler 486 dofs, level dofs 27 108 180 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.5437 +DEAL:1:cg::Convergence step 16 value 8.06366e-13 +DEAL:1::Estimate 0.474987 +DEAL:1::Step 5 +DEAL:1::Triangulation 84 cells, 7 levels +DEAL:1::DoFHandler 756 dofs, level dofs 27 108 324 216 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 25.5435 +DEAL:1:cg::Convergence step 17 value 8.48262e-13 +DEAL:1::Estimate 0.356044 +