]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add test case mpi/tria_signals_06: self test of tria_signals test driver. 1743/head
authorLei Qiao <qiaol618@gmail.com>
Thu, 15 Oct 2015 17:10:08 +0000 (12:10 -0500)
committerLei Qiao <qiaol618@gmail.com>
Thu, 22 Oct 2015 14:56:38 +0000 (09:56 -0500)
tests/mpi/tria_signals_06.cc [new file with mode: 0644]
tests/mpi/tria_signals_06.mpirun=4.output [new file with mode: 0644]

diff --git a/tests/mpi/tria_signals_06.cc b/tests/mpi/tria_signals_06.cc
new file mode 100644 (file)
index 0000000..cda82fc
--- /dev/null
@@ -0,0 +1,202 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/base/std_cxx11/bind.h>
+
+#include <fstream>
+#include <ostream>
+
+// This is a test for test driver itself. Signals on refinement is disabled,
+// so non-positive n_active_cell_gap is expected. Thus we can make sure the
+// SignalListener class is not just report zero results blindly, and both
+// signals on refinement and coarsening are working.
+//
+// However, if this test failed, it is not necessarily the signals are wrong.
+// The reason can also be that the adaptation mechanism is changed.
+//
+// This test is modified from tria_signals_04.
+
+template<int dim, int spacedim>
+class SignalListener
+{
+public:
+  SignalListener(Triangulation<dim, spacedim> &tria_in)
+    :
+    n_active_cells(tria_in.n_active_cells()),
+    tria(tria_in)
+  {
+    // tria_in.signals.post_refinement_on_cell.connect
+    // (std_cxx11::bind (&SignalListener<dim, spacedim>::count_on_refine,
+    //                   this,
+    //                   std_cxx11::placeholders::_1));
+
+    tria_in.signals.pre_coarsening_on_cell.connect
+    (std_cxx11::bind (&SignalListener<dim, spacedim>::count_on_coarsen,
+                      this,
+                      std_cxx11::placeholders::_1));
+  }
+
+  int n_active_cell_gap()
+  {
+    return (n_active_cells -
+            static_cast<int> (tria.n_active_cells()));
+  }
+
+private:
+  void count_on_refine(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
+  {
+    n_active_cells += cell->n_children();
+    --n_active_cells;
+
+    return;
+  }
+
+  void count_on_coarsen(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
+  {
+    ++n_active_cells;
+    n_active_cells -= cell->n_children();
+
+    return;
+  }
+
+  int n_active_cells;
+  const Triangulation<dim, spacedim> &tria;
+};
+
+
+template<int dim, int spacedim>
+void test()
+{
+  typedef parallel::distributed::Triangulation<dim, spacedim> TriaType;
+
+  {
+    const std::string prefix = Utilities::int_to_string (dim, 1) +
+                               "d-" +
+                               Utilities::int_to_string (spacedim, 1)
+                               + "d";
+    deallog.push(prefix.c_str());
+  }
+
+  // Option dealii::Triangulation<dim, spacedim>::maximum_smoothing can't
+  // run in parallel at the time that this test is created.
+  TriaType tria(MPI_COMM_WORLD,
+                typename dealii::Triangulation<dim, spacedim>::MeshSmoothing
+                (dealii::Triangulation<dim, spacedim>::smoothing_on_refinement |
+                 dealii::Triangulation<dim, spacedim>::smoothing_on_coarsening));
+
+
+  GridGenerator::hyper_cube(tria);
+  SignalListener<dim, spacedim> count_cell_via_signal(tria);
+
+  tria.refine_global(1);
+
+
+  // The following loop is borrowed from p4est_3d_refine_01 with some modifications.
+  for  (int n_loop = 0;
+        // Terminate loop on global information to prevent premature termination
+        // on only part of processors. (n_loop < 20) is just a passive safety to
+        // avoid infinite loop.
+        (tria.n_global_active_cells() < 20000) && (n_loop < 20);
+        ++n_loop)
+    {
+      std::vector<bool> flags (tria.n_active_cells(), false);
+
+      // Refine one fifth of all cells each time (but at least one).
+      // Note that only the own marked cells will be refined.
+      // But refine flags on own cells could be effected by flags on ghost cells
+      // through mesh smoothing.
+      for (unsigned int i=0; i<tria.n_active_cells() / 3 + 1; ++i)
+        {
+          const unsigned int x = Testing::rand() % flags.size();
+          flags[x] = true;
+        }
+
+      unsigned int index=0;
+      unsigned int locals=0;
+
+      for (typename Triangulation<dim, spacedim>::active_cell_iterator
+           cell = tria.begin_active();
+           cell != tria.end(); ++cell, ++index)
+        if (flags[index])
+          {
+            if (cell->is_locally_owned())
+              ++locals;
+            cell->set_refine_flag();
+          }
+
+      if (locals > 5)
+        {
+          // Coarsen some cells randomly only if we have enough local cells
+          // marked to be refined
+          std::fill(flags.begin(), flags.end(), false);
+          for (unsigned int i=0; i<tria.n_active_cells() / 3; ++i)
+            {
+              const unsigned int x = Testing::rand() % flags.size();
+              flags[x] = true;
+            }
+
+          index=0;
+          for (typename Triangulation<dim, spacedim>::active_cell_iterator
+               cell = tria.begin_active();
+               cell != tria.end(); ++cell, ++index)
+            if (flags[index] && ! cell->refine_flag_set())
+              {
+                cell->set_coarsen_flag();
+              }
+        }
+
+      tria.execute_coarsening_and_refinement ();
+
+      deallog << "n_loop: " << n_loop
+              << ", n_cell_gap: "
+              << count_cell_via_signal.n_active_cell_gap() << std::endl;
+    }
+
+  deallog.pop();
+  return;
+}
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, /* int max_num_threads */ 1);
+  MPILogInitAll log;
+
+  // parallel::distributed::Triangulation<1, spacedim> is not valid.
+  {
+    const int dim = 2;
+    const int spacedim = 2;
+    test<dim,spacedim> ();
+  }
+
+  {
+    const int dim = 2;
+    const int spacedim = 3;
+    test<dim,spacedim> ();
+  }
+
+  {
+    const int dim = 3;
+    const int spacedim = 3;
+    test<dim,spacedim> ();
+  }
+
+  return (0);
+}
+
diff --git a/tests/mpi/tria_signals_06.mpirun=4.output b/tests/mpi/tria_signals_06.mpirun=4.output
new file mode 100644 (file)
index 0000000..daeb541
--- /dev/null
@@ -0,0 +1,135 @@
+
+DEAL:0:2d-2d::n_loop: 0, n_cell_gap: -9
+DEAL:0:2d-2d::n_loop: 1, n_cell_gap: -24
+DEAL:0:2d-2d::n_loop: 2, n_cell_gap: -51
+DEAL:0:2d-2d::n_loop: 3, n_cell_gap: -84
+DEAL:0:2d-2d::n_loop: 4, n_cell_gap: -150
+DEAL:0:2d-2d::n_loop: 5, n_cell_gap: -303
+DEAL:0:2d-2d::n_loop: 6, n_cell_gap: -564
+DEAL:0:2d-2d::n_loop: 7, n_cell_gap: -1059
+DEAL:0:2d-2d::n_loop: 8, n_cell_gap: -1953
+DEAL:0:2d-2d::n_loop: 9, n_cell_gap: -3609
+DEAL:0:2d-2d::n_loop: 10, n_cell_gap: -7116
+DEAL:0:2d-3d::n_loop: 0, n_cell_gap: -3
+DEAL:0:2d-3d::n_loop: 1, n_cell_gap: -3
+DEAL:0:2d-3d::n_loop: 2, n_cell_gap: -9
+DEAL:0:2d-3d::n_loop: 3, n_cell_gap: -30
+DEAL:0:2d-3d::n_loop: 4, n_cell_gap: -48
+DEAL:0:2d-3d::n_loop: 5, n_cell_gap: -72
+DEAL:0:2d-3d::n_loop: 6, n_cell_gap: -132
+DEAL:0:2d-3d::n_loop: 7, n_cell_gap: -219
+DEAL:0:2d-3d::n_loop: 8, n_cell_gap: -369
+DEAL:0:2d-3d::n_loop: 9, n_cell_gap: -711
+DEAL:0:2d-3d::n_loop: 10, n_cell_gap: -1413
+DEAL:0:2d-3d::n_loop: 11, n_cell_gap: -2628
+DEAL:0:2d-3d::n_loop: 12, n_cell_gap: -4875
+DEAL:0:2d-3d::n_loop: 13, n_cell_gap: -9405
+DEAL:0:3d-3d::n_loop: 0, n_cell_gap: -21
+DEAL:0:3d-3d::n_loop: 1, n_cell_gap: -63
+DEAL:0:3d-3d::n_loop: 2, n_cell_gap: -189
+DEAL:0:3d-3d::n_loop: 3, n_cell_gap: -455
+DEAL:0:3d-3d::n_loop: 4, n_cell_gap: -1603
+DEAL:0:3d-3d::n_loop: 5, n_cell_gap: -5439
+DEAL:0:3d-3d::n_loop: 6, n_cell_gap: -19705
+
+DEAL:1:2d-2d::n_loop: 0, n_cell_gap: -9
+DEAL:1:2d-2d::n_loop: 1, n_cell_gap: -27
+DEAL:1:2d-2d::n_loop: 2, n_cell_gap: -72
+DEAL:1:2d-2d::n_loop: 3, n_cell_gap: -120
+DEAL:1:2d-2d::n_loop: 4, n_cell_gap: -219
+DEAL:1:2d-2d::n_loop: 5, n_cell_gap: -414
+DEAL:1:2d-2d::n_loop: 6, n_cell_gap: -687
+DEAL:1:2d-2d::n_loop: 7, n_cell_gap: -1227
+DEAL:1:2d-2d::n_loop: 8, n_cell_gap: -2244
+DEAL:1:2d-2d::n_loop: 9, n_cell_gap: -4041
+DEAL:1:2d-2d::n_loop: 10, n_cell_gap: -7761
+DEAL:1:2d-3d::n_loop: 0, n_cell_gap: 0
+DEAL:1:2d-3d::n_loop: 1, n_cell_gap: 0
+DEAL:1:2d-3d::n_loop: 2, n_cell_gap: 0
+DEAL:1:2d-3d::n_loop: 3, n_cell_gap: -30
+DEAL:1:2d-3d::n_loop: 4, n_cell_gap: -48
+DEAL:1:2d-3d::n_loop: 5, n_cell_gap: -81
+DEAL:1:2d-3d::n_loop: 6, n_cell_gap: -138
+DEAL:1:2d-3d::n_loop: 7, n_cell_gap: -273
+DEAL:1:2d-3d::n_loop: 8, n_cell_gap: -420
+DEAL:1:2d-3d::n_loop: 9, n_cell_gap: -858
+DEAL:1:2d-3d::n_loop: 10, n_cell_gap: -1569
+DEAL:1:2d-3d::n_loop: 11, n_cell_gap: -2913
+DEAL:1:2d-3d::n_loop: 12, n_cell_gap: -5373
+DEAL:1:2d-3d::n_loop: 13, n_cell_gap: -10128
+DEAL:1:3d-3d::n_loop: 0, n_cell_gap: -21
+DEAL:1:3d-3d::n_loop: 1, n_cell_gap: -77
+DEAL:1:3d-3d::n_loop: 2, n_cell_gap: -238
+DEAL:1:3d-3d::n_loop: 3, n_cell_gap: -602
+DEAL:1:3d-3d::n_loop: 4, n_cell_gap: -1981
+DEAL:1:3d-3d::n_loop: 5, n_cell_gap: -6608
+DEAL:1:3d-3d::n_loop: 6, n_cell_gap: -22239
+
+
+DEAL:2:2d-2d::n_loop: 0, n_cell_gap: 0
+DEAL:2:2d-2d::n_loop: 1, n_cell_gap: -33
+DEAL:2:2d-2d::n_loop: 2, n_cell_gap: -69
+DEAL:2:2d-2d::n_loop: 3, n_cell_gap: -129
+DEAL:2:2d-2d::n_loop: 4, n_cell_gap: -207
+DEAL:2:2d-2d::n_loop: 5, n_cell_gap: -411
+DEAL:2:2d-2d::n_loop: 6, n_cell_gap: -669
+DEAL:2:2d-2d::n_loop: 7, n_cell_gap: -1188
+DEAL:2:2d-2d::n_loop: 8, n_cell_gap: -2169
+DEAL:2:2d-2d::n_loop: 9, n_cell_gap: -3954
+DEAL:2:2d-2d::n_loop: 10, n_cell_gap: -7644
+DEAL:2:2d-3d::n_loop: 0, n_cell_gap: 0
+DEAL:2:2d-3d::n_loop: 1, n_cell_gap: 0
+DEAL:2:2d-3d::n_loop: 2, n_cell_gap: -9
+DEAL:2:2d-3d::n_loop: 3, n_cell_gap: -24
+DEAL:2:2d-3d::n_loop: 4, n_cell_gap: -54
+DEAL:2:2d-3d::n_loop: 5, n_cell_gap: -99
+DEAL:2:2d-3d::n_loop: 6, n_cell_gap: -162
+DEAL:2:2d-3d::n_loop: 7, n_cell_gap: -276
+DEAL:2:2d-3d::n_loop: 8, n_cell_gap: -426
+DEAL:2:2d-3d::n_loop: 9, n_cell_gap: -786
+DEAL:2:2d-3d::n_loop: 10, n_cell_gap: -1479
+DEAL:2:2d-3d::n_loop: 11, n_cell_gap: -2769
+DEAL:2:2d-3d::n_loop: 12, n_cell_gap: -5133
+DEAL:2:2d-3d::n_loop: 13, n_cell_gap: -9693
+DEAL:2:3d-3d::n_loop: 0, n_cell_gap: -21
+DEAL:2:3d-3d::n_loop: 1, n_cell_gap: -77
+DEAL:2:3d-3d::n_loop: 2, n_cell_gap: -210
+DEAL:2:3d-3d::n_loop: 3, n_cell_gap: -616
+DEAL:2:3d-3d::n_loop: 4, n_cell_gap: -1967
+DEAL:2:3d-3d::n_loop: 5, n_cell_gap: -6657
+DEAL:2:3d-3d::n_loop: 6, n_cell_gap: -22134
+
+
+DEAL:3:2d-2d::n_loop: 0, n_cell_gap: -9
+DEAL:3:2d-2d::n_loop: 1, n_cell_gap: -33
+DEAL:3:2d-2d::n_loop: 2, n_cell_gap: -57
+DEAL:3:2d-2d::n_loop: 3, n_cell_gap: -108
+DEAL:3:2d-2d::n_loop: 4, n_cell_gap: -171
+DEAL:3:2d-2d::n_loop: 5, n_cell_gap: -312
+DEAL:3:2d-2d::n_loop: 6, n_cell_gap: -555
+DEAL:3:2d-2d::n_loop: 7, n_cell_gap: -1065
+DEAL:3:2d-2d::n_loop: 8, n_cell_gap: -1995
+DEAL:3:2d-2d::n_loop: 9, n_cell_gap: -3666
+DEAL:3:2d-2d::n_loop: 10, n_cell_gap: -7161
+DEAL:3:2d-3d::n_loop: 0, n_cell_gap: 0
+DEAL:3:2d-3d::n_loop: 1, n_cell_gap: 0
+DEAL:3:2d-3d::n_loop: 2, n_cell_gap: -9
+DEAL:3:2d-3d::n_loop: 3, n_cell_gap: -21
+DEAL:3:2d-3d::n_loop: 4, n_cell_gap: -39
+DEAL:3:2d-3d::n_loop: 5, n_cell_gap: -69
+DEAL:3:2d-3d::n_loop: 6, n_cell_gap: -120
+DEAL:3:2d-3d::n_loop: 7, n_cell_gap: -237
+DEAL:3:2d-3d::n_loop: 8, n_cell_gap: -363
+DEAL:3:2d-3d::n_loop: 9, n_cell_gap: -720
+DEAL:3:2d-3d::n_loop: 10, n_cell_gap: -1377
+DEAL:3:2d-3d::n_loop: 11, n_cell_gap: -2619
+DEAL:3:2d-3d::n_loop: 12, n_cell_gap: -4908
+DEAL:3:2d-3d::n_loop: 13, n_cell_gap: -9354
+DEAL:3:3d-3d::n_loop: 0, n_cell_gap: -21
+DEAL:3:3d-3d::n_loop: 1, n_cell_gap: -77
+DEAL:3:3d-3d::n_loop: 2, n_cell_gap: -238
+DEAL:3:3d-3d::n_loop: 3, n_cell_gap: -616
+DEAL:3:3d-3d::n_loop: 4, n_cell_gap: -1764
+DEAL:3:3d-3d::n_loop: 5, n_cell_gap: -6321
+DEAL:3:3d-3d::n_loop: 6, n_cell_gap: -20755
+

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.