]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
test number_cache for shared_tria
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 29 Apr 2014 22:36:03 +0000 (22:36 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 29 Apr 2014 22:36:03 +0000 (22:36 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_sharedtria@32859 0785d39b-7218-0410-832d-ea1e28bc413d

tests/sharedtria/dof_01.cc [new file with mode: 0644]
tests/sharedtria/dof_01.mpirun=3.output [new file with mode: 0644]
tests/sharedtria/dof_01.output [new file with mode: 0644]

diff --git a/tests/sharedtria/dof_01.cc b/tests/sharedtria/dof_01.cc
new file mode 100644 (file)
index 0000000..e844631
--- /dev/null
@@ -0,0 +1,151 @@
+// ---------------------------------------------------------------------
+// $Id: dof_handler_number_cache.cc 31761 2013-11-22 14:42:37Z heister $
+//
+// Copyright (C) 2008 - 2013 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check number cache for shared_tria
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/distributed/shared_tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/intergrid_map.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+
+#include <fstream>
+#include <cstdlib>
+#include <numeric>
+
+
+template<int dim>
+void test()
+{
+  parallel::shared::Triangulation<dim>
+  triangulation (MPI_COMM_WORLD);
+
+  FESystem<dim> fe (FE_Q<dim>(3),2,
+                    FE_DGQ<dim>(1),1);
+
+  DoFHandler<dim> dof_handler (triangulation);
+
+  GridGenerator::hyper_cube(triangulation);
+  triangulation.refine_global (2);
+
+  const unsigned int n_refinements[] = { 0, 4, 3, 2 };
+  for (unsigned int i=0; i<n_refinements[dim]; ++i)
+    {
+      // refine one-fifth of cells randomly
+      std::vector<bool> flags (triangulation.n_active_cells(), false);
+      for (unsigned int k=0; k<flags.size()/5 + 1; ++k)
+        flags[Testing::rand() % flags.size()] = true;
+      // make sure there's at least one that
+      // will be refined
+      flags[0] = true;
+
+      // refine triangulation
+      unsigned int index=0;
+      for (typename Triangulation<dim>::active_cell_iterator
+           cell = triangulation.begin_active();
+           cell != triangulation.end(); ++cell)
+          {
+            if (flags[index])
+              cell->set_refine_flag();
+            ++index;
+          }
+
+      Assert (index <= triangulation.n_active_cells(), ExcInternalError());
+
+      // flag all other cells for coarsening
+      // (this should ensure that at least
+      // some of them will actually be
+      // coarsened)
+      index=0;
+      for (typename Triangulation<dim>::active_cell_iterator
+           cell = triangulation.begin_active();
+           cell != triangulation.end(); ++cell)
+          {
+            if (!flags[index])
+              cell->set_coarsen_flag();
+            ++index;
+          }
+
+      triangulation.execute_coarsening_and_refinement ();
+      dof_handler.distribute_dofs (fe);
+
+      deallog
+       << "n_dofs: " << dof_handler.n_dofs() << std::endl
+       << "n_locally_owned_dofs: " << dof_handler.n_locally_owned_dofs() << std::endl;
+
+      deallog << "n_locally_owned_dofs_per_processor: ";
+      std::vector<unsigned int> v = dof_handler.n_locally_owned_dofs_per_processor();
+      unsigned int sum = 0;
+      for (unsigned int i=0;i<v.size();++i)
+       {
+         deallog << v[i] << " ";
+         sum += v[i];
+       }
+      deallog << " sum: " << sum << std::endl;
+
+      Assert(dof_handler.n_locally_owned_dofs() == dof_handler.n_locally_owned_dofs_per_processor()[triangulation.locally_owned_subdomain()], ExcInternalError());
+      
+
+      Assert( dof_handler.n_locally_owned_dofs() == dof_handler.locally_owned_dofs().n_elements(), ExcInternalError());
+
+      const unsigned int N = dof_handler.n_dofs();
+
+      Assert (dof_handler.n_locally_owned_dofs() <= N,
+              ExcInternalError());
+      Assert (std::accumulate (dof_handler.n_locally_owned_dofs_per_processor().begin(),
+                               dof_handler.n_locally_owned_dofs_per_processor().end(),
+                               0U) == N,
+              ExcInternalError());
+
+      IndexSet all (N);
+      for (unsigned int i=0;
+          i<dof_handler.locally_owned_dofs_per_processor().size(); ++i)
+       {
+         IndexSet intersect = all & dof_handler.locally_owned_dofs_per_processor()[i];
+         Assert(intersect.n_elements()==0, ExcInternalError());
+
+         all.add_indices(dof_handler.locally_owned_dofs_per_processor()[i]);
+       }
+
+      Assert(all == complete_index_set(N), ExcInternalError());
+    }
+}
+
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  MPILogInitAll all;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/sharedtria/dof_01.mpirun=3.output b/tests/sharedtria/dof_01.mpirun=3.output
new file mode 100644 (file)
index 0000000..892356f
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL:0:2d::n_dofs: 818
+DEAL:0:2d::n_locally_owned_dofs: 818
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 818  sum: 818
+DEAL:0:2d::n_dofs: 1874
+DEAL:0:2d::n_locally_owned_dofs: 1874
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 1874  sum: 1874
+DEAL:0:2d::n_dofs: 2800
+DEAL:0:2d::n_locally_owned_dofs: 2800
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 2800  sum: 2800
+DEAL:0:3d::n_dofs: 13282
+DEAL:0:3d::n_locally_owned_dofs: 13282
+DEAL:0:3d::n_locally_owned_dofs_per_processor: 13282  sum: 13282
+DEAL:0:3d::n_dofs: 42960
+DEAL:0:3d::n_locally_owned_dofs: 42960
+DEAL:0:3d::n_locally_owned_dofs_per_processor: 42960  sum: 42960
diff --git a/tests/sharedtria/dof_01.output b/tests/sharedtria/dof_01.output
new file mode 100644 (file)
index 0000000..d074284
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL:0:2d::n_dofs: 818
+DEAL:0:2d::n_locally_owned_dofs: 818
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 818  sum: 818
+DEAL:0:2d::n_dofs: 1754
+DEAL:0:2d::n_locally_owned_dofs: 1754
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 1754  sum: 1754
+DEAL:0:2d::n_dofs: 3056
+DEAL:0:2d::n_locally_owned_dofs: 3056
+DEAL:0:2d::n_locally_owned_dofs_per_processor: 3056  sum: 3056
+DEAL:0:3d::n_dofs: 13282
+DEAL:0:3d::n_locally_owned_dofs: 13282
+DEAL:0:3d::n_locally_owned_dofs_per_processor: 13282  sum: 13282
+DEAL:0:3d::n_dofs: 41826
+DEAL:0:3d::n_locally_owned_dofs: 41826
+DEAL:0:3d::n_locally_owned_dofs_per_processor: 41826  sum: 41826

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.