]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add broken test
authorTimo Heister <timo.heister@gmail.com>
Fri, 24 May 2013 12:56:20 +0000 (12:56 +0000)
committerTimo Heister <timo.heister@gmail.com>
Fri, 24 May 2013 12:56:20 +0000 (12:56 +0000)
git-svn-id: https://svn.dealii.org/trunk@29566 0785d39b-7218-0410-832d-ea1e28bc413d

tests/mpi/mg_04.cc [new file with mode: 0644]
tests/mpi/mg_04/ncpu_1/cmp/generic [new file with mode: 0644]
tests/mpi/mg_04/ncpu_10/cmp/generic [new file with mode: 0644]
tests/mpi/mg_04/ncpu_2/cmp/generic [new file with mode: 0644]
tests/mpi/mg_04/ncpu_4/cmp/generic [new file with mode: 0644]

diff --git a/tests/mpi/mg_04.cc b/tests/mpi/mg_04.cc
new file mode 100644 (file)
index 0000000..0c7e531
--- /dev/null
@@ -0,0 +1,138 @@
+//---------------------------------------------------------------------------
+//    $Id: p4est_data_out_01.cc 25200 2012-03-02 21:11:21Z heister $
+//    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.
+//
+//---------------------------------------------------------------------------
+
+
+// distributing FE_Q multigrid DoFs fails
+/*
+*/
+
+#include "../tests.h"
+#include "coarse_grid_common.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/numerics/data_out.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+#include <fstream>
+
+template<int dim>
+void output(parallel::distributed::Triangulation<dim> & tr)
+{
+  const std::string filename = ("mg_04/mesh." +
+                                Utilities::int_to_string
+                                (tr.locally_owned_subdomain(), 4) +
+                                ".vtu");
+
+
+}
+
+template<int dim>
+void test()
+{
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+  if (myid == 0)
+    deallog << "hyper_cube" << std::endl;
+
+  parallel::distributed::Triangulation<dim> tr(MPI_COMM_WORLD,
+                                              Triangulation<dim>::none,
+                                              parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+  GridGenerator::hyper_ball (tr);
+  //GridGenerator::hyper_cube(tr);
+  tr.refine_global(1);
+  DoFHandler<dim> dofh(tr);
+
+  output(tr);
+
+  static const FE_Q<dim> fe(1);
+  dofh.distribute_dofs (fe);
+  dofh.distribute_mg_dofs (fe);
+
+  {
+   std::ofstream file((std::string("mg_04/ncpu_") + Utilities::int_to_string(Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD)) + "/dat." + Utilities::int_to_string(myid)).c_str());
+    file << "**** proc " << myid << std::endl;
+          for (unsigned int lvl=0;lvl<tr.n_levels();++lvl)
+            {
+              file << "level " << lvl << ": ";
+              typename DoFHandler<dim>::cell_iterator
+              cell = dofh.begin(lvl),
+              endc = dofh.end(lvl);
+
+              for (;cell!=endc;++cell)
+                {
+              std::vector<unsigned int> dofs(fe.n_dofs_per_cell());
+              cell->get_mg_dof_indices(dofs);
+
+              for (unsigned int i=0;i<dofs.size();++i)
+                if (dofs[i]==numbers::invalid_dof_index)
+                  file << "- ";
+                else
+                  file << dofs[i] << " ";
+              file << " | ";
+                }
+              file << std::endl;
+            }
+  }
+
+  MPI_Barrier(MPI_COMM_WORLD);
+
+  if (myid==0)
+    {
+      for (unsigned int i=0;i<Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);++i)
+        {
+          cat_file((std::string("mg_04/ncpu_") + Utilities::int_to_string(Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD)) + "/dat." + Utilities::int_to_string(i)).c_str());
+        }
+    }      
+
+  if (myid==0)
+    deallog << "OK" << std::endl;
+}
+
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile(output_file_for_mpi("mg_04").c_str());
+      deallog.attach(logfile);
+      deallog.depth_console(0);
+      deallog.threshold_double(1.e-10);
+
+      deallog.push("2d");
+      test<2>();
+      deallog.pop();
+    }
+  else
+    test<2>();
+
+}
diff --git a/tests/mpi/mg_04/ncpu_1/cmp/generic b/tests/mpi/mg_04/ncpu_1/cmp/generic
new file mode 100644 (file)
index 0000000..9140b59
--- /dev/null
@@ -0,0 +1,8 @@
+
+DEAL:0:2d::hyper_cube
+**** proc 0
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: 0  | 1  | 2  | 3  | 4  | 5  | 6  | 7  | 8  | 9  | 10  | 11  | 12  | 13  | 14  | 15  | 
+
+DEAL:0:2d::OK
diff --git a/tests/mpi/mg_04/ncpu_10/cmp/generic b/tests/mpi/mg_04/ncpu_10/cmp/generic
new file mode 100644 (file)
index 0000000..a3c35cc
--- /dev/null
@@ -0,0 +1,41 @@
+
+DEAL:0:2d::hyper_cube
+**** proc 0
+level 0: -  | 
+
+**** proc 1
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: 0  | 1  | 2  | 3  | 4  | -  | 6  | -  | 8  | 9  | -  | -  | 12  | -  | -  | -  | 
+
+**** proc 2
+level 0: -  | 
+
+**** proc 3
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | 1  | -  | 3  | 4  | 5  | 6  | 7  | -  | 9  | -  | -  | 12  | 13  | -  | -  | 
+
+**** proc 4
+level 0: -  | 
+
+**** proc 5
+level 0: -  | 
+
+**** proc 6
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | 2  | 3  | -  | -  | 6  | -  | 8  | 9  | 10  | 11  | 12  | -  | 14  | -  | 
+
+**** proc 7
+level 0: -  | 
+
+**** proc 8
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | -  | 3  | -  | -  | 6  | 7  | -  | 9  | -  | 11  | 12  | 13  | 14  | 15  | 
+
+**** proc 9
+level 0: -  | 
+
+DEAL:0:2d::OK
diff --git a/tests/mpi/mg_04/ncpu_2/cmp/generic b/tests/mpi/mg_04/ncpu_2/cmp/generic
new file mode 100644 (file)
index 0000000..d548817
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:0:2d::hyper_cube
+**** proc 0
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: 0  | 1  | 2  | 3  | 4  | -  | 6  | -  | 8  | 9  | -  | -  | 12  | -  | -  | -  | 
+
+**** proc 1
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | 1  | -  | 3  | 4  | 5  | 6  | 7  | -  | 9  | -  | -  | 12  | 13  | -  | -  | 
+
+**** proc 2
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | 2  | 3  | -  | -  | 6  | -  | 8  | 9  | 10  | 11  | 12  | -  | 14  | -  | 
+
+**** proc 3
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | -  | 3  | -  | -  | 6  | 7  | -  | 9  | -  | 11  | 12  | 13  | 14  | 15  | 
+
+DEAL:0:2d::OK
diff --git a/tests/mpi/mg_04/ncpu_4/cmp/generic b/tests/mpi/mg_04/ncpu_4/cmp/generic
new file mode 100644 (file)
index 0000000..d548817
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL:0:2d::hyper_cube
+**** proc 0
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: 0  | 1  | 2  | 3  | 4  | -  | 6  | -  | 8  | 9  | -  | -  | 12  | -  | -  | -  | 
+
+**** proc 1
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | 1  | -  | 3  | 4  | 5  | 6  | 7  | -  | 9  | -  | -  | 12  | 13  | -  | -  | 
+
+**** proc 2
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | 2  | 3  | -  | -  | 6  | -  | 8  | 9  | 10  | 11  | 12  | -  | 14  | -  | 
+
+**** proc 3
+level 0: 0  | 
+level 1: 0  | 1  | 2  | 3  | 
+level 2: -  | -  | -  | 3  | -  | -  | 6  | 7  | -  | 9  | -  | 11  | 12  | 13  | 14  | 15  | 
+
+DEAL:0:2d::OK

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.