]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix orientation of 1D periodic case 1649/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 22 Sep 2015 16:04:22 +0000 (18:04 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 22 Sep 2015 16:04:22 +0000 (18:04 +0200)
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/periodicity_1d.cc [new file with mode: 0644]
tests/grid/periodicity_1d.output [new file with mode: 0644]

index 96394ba5844c2cb4c932e912f6ba96ccb0298923..dce6a7f6c7abd0c0de17b1c0c2d63047fd62f606 100644 (file)
@@ -2662,7 +2662,7 @@ next_cell:
     static inline std::bitset<3> lookup (const MATCH_T &)
     {
       // The 1D case is trivial
-      return 4; // [true ,false,false]
+      return 1; // [true ,false,false]
     }
   };
 
index 27ab3513ddf45df6dd89a5d8a892c074232ad852..e6a85a182782d5a24bd2961faea7839192889abf 100644 (file)
@@ -96,7 +96,7 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
   template
     std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> >::type>
     compute_active_cell_halo_layer (const parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> &,
-                                    const std_cxx11::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, 
+                                    const std_cxx11::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension,
                                                                                                                parallel::distributed::Triangulation<deal_II_dimension, deal_II_space_dimension> >::type&)> &);
 
   template
@@ -311,26 +311,22 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
                                                 const Tensor<1,deal_II_space_dimension> &,
                                                 const FullMatrix<double> &);
 
-    #if deal_II_dimension >= 2
-
-      template
-      void collect_periodic_faces<X> (const X &,
-                                      const types::boundary_id,
-                                      const types::boundary_id,
-                                      const int,
-                                      std::vector<PeriodicFacePair<X::cell_iterator> > &,
-                                      const Tensor<1,X::space_dimension> &,
-                                      const FullMatrix<double> &);
-
-      template
-      void collect_periodic_faces<X> (const X &,
-                                      const types::boundary_id,
-                                      const int,
-                                      std::vector<PeriodicFacePair<X::cell_iterator> > &,
-                                      const Tensor<1,X::space_dimension> &,
-                                      const FullMatrix<double> &);
+    template
+    void collect_periodic_faces<X> (const X &,
+                                    const types::boundary_id,
+                                    const types::boundary_id,
+                                    const int,
+                                    std::vector<PeriodicFacePair<X::cell_iterator> > &,
+                                    const Tensor<1,X::space_dimension> &,
+                                    const FullMatrix<double> &);
 
-    #endif
+    template
+    void collect_periodic_faces<X> (const X &,
+                                    const types::boundary_id,
+                                    const int,
+                                    std::vector<PeriodicFacePair<X::cell_iterator> > &,
+                                    const Tensor<1,X::space_dimension> &,
+                                    const FullMatrix<double> &);
 
   \}
 #endif
@@ -379,4 +375,3 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
 \}
 #endif
 }
-
diff --git a/tests/grid/periodicity_1d.cc b/tests/grid/periodicity_1d.cc
new file mode 100644 (file)
index 0000000..3806b84
--- /dev/null
@@ -0,0 +1,70 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2007 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Check periodic boundary faces
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_tools.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/base/logstream.h>
+
+#include <fstream>
+#include <string>
+
+
+
+template <int dim>
+void check()
+{
+  Triangulation<dim> tr;
+  Point<dim> p0, p1;
+  for (unsigned int d=0; d<dim; ++d)
+    p1[d] = 1+d;
+  std::vector<unsigned int> refinements(dim, 3);
+  GridGenerator::subdivided_hyper_rectangle(tr, refinements, p0, p1, true);
+  tr.refine_global(4-dim);
+
+  std::vector<GridTools::PeriodicFacePair<typename Triangulation<dim>::cell_iterator> > periodic_faces;
+  for (unsigned int d=0; d<dim; ++d)
+    GridTools::collect_periodic_faces(tr, 2*d, 2*d+1, d, periodic_faces);
+
+  deallog << "Test run in " << dim << " dimensions" << std::endl;
+  for (unsigned int i=0; i<periodic_faces.size(); ++i)
+    {
+      deallog << periodic_faces[i].cell[0]->index() << " "
+              << periodic_faces[i].cell[1]->index() << " "
+              << periodic_faces[i].face_idx[0] << " "
+              << periodic_faces[i].face_idx[1] << " "
+              << periodic_faces[i].orientation
+              << std::endl;
+    }
+}
+
+
+int main()
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check<1>();
+  check<2>();
+  check<3>();
+}
diff --git a/tests/grid/periodicity_1d.output b/tests/grid/periodicity_1d.output
new file mode 100644 (file)
index 0000000..7f620fc
--- /dev/null
@@ -0,0 +1,38 @@
+
+DEAL::Test run in 1 dimensions
+DEAL::0 2 0 1 001
+DEAL::Test run in 2 dimensions
+DEAL::0 2 0 1 001
+DEAL::3 5 0 1 001
+DEAL::6 8 0 1 001
+DEAL::0 6 2 3 001
+DEAL::1 7 2 3 001
+DEAL::2 8 2 3 001
+DEAL::Test run in 3 dimensions
+DEAL::0 2 0 1 001
+DEAL::3 5 0 1 001
+DEAL::6 8 0 1 001
+DEAL::9 11 0 1 001
+DEAL::12 14 0 1 001
+DEAL::15 17 0 1 001
+DEAL::18 20 0 1 001
+DEAL::21 23 0 1 001
+DEAL::24 26 0 1 001
+DEAL::0 6 2 3 001
+DEAL::1 7 2 3 001
+DEAL::2 8 2 3 001
+DEAL::9 15 2 3 001
+DEAL::10 16 2 3 001
+DEAL::11 17 2 3 001
+DEAL::18 24 2 3 001
+DEAL::19 25 2 3 001
+DEAL::20 26 2 3 001
+DEAL::0 18 4 5 001
+DEAL::1 19 4 5 001
+DEAL::2 20 4 5 001
+DEAL::3 21 4 5 001
+DEAL::4 22 4 5 001
+DEAL::5 23 4 5 001
+DEAL::6 24 4 5 001
+DEAL::7 25 4 5 001
+DEAL::8 26 4 5 001

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.