]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix disabling bug in hp::DoFHandler::n_boundary_dofs().
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 22 Jul 2013 02:46:30 +0000 (02:46 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 22 Jul 2013 02:46:30 +0000 (02:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@30085 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/source/hp/dof_handler.cc
tests/hp/n_boundary_dofs.cc [new file with mode: 0644]
tests/hp/n_boundary_dofs/cmp/generic [new file with mode: 0644]

index 1ff6e2ef50c0ab653be706462977738788f478f9..a47b3af016b1516b2c63787b7aba5faa680a9db8 100644 (file)
@@ -164,6 +164,12 @@ this function.
 
 <ol>
 
+<li>Fixed: hp::DoFHandler::n_boundary_dofs() had a bug that always led
+to a failed assertion. This is now fixed.
+<br>
+(Wolfgang Bangerth, 2013/07/21)
+</li>
+
 <li>Fixed: VectorTools::project has an option to first project onto the
 boundary. However, the implementation of this option ignored the mapping
 that is provided to the function. This is now fixed.
index c87c3f5dcb79ec70af2e116ceaadfec33aabeac9..0b4b865540e5cc5c2fccc76ea9250f35c066e5a3 100644 (file)
@@ -1955,7 +1955,8 @@ namespace hp
             const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
             dofs_on_face.resize (dofs_per_face);
 
-            cell->face(f)->get_dof_indices (dofs_on_face);
+            cell->face(f)->get_dof_indices (dofs_on_face,
+                                           cell->active_fe_index());
             for (unsigned int i=0; i<dofs_per_face; ++i)
               boundary_dofs.insert(dofs_on_face[i]);
           };
diff --git a/tests/hp/n_boundary_dofs.cc b/tests/hp/n_boundary_dofs.cc
new file mode 100644 (file)
index 0000000..a31c7f8
--- /dev/null
@@ -0,0 +1,86 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// 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 the consistency of the number cache of DoFHandler for a sequential
+// object. Like deal.II/n_boundary_dofs but for an hp object
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/grid/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/hp/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 <deal.II/hp/fe_collection.h>
+
+#include <fstream>
+#include <cstdlib>
+
+
+template<int dim>
+void test()
+{
+  Triangulation<dim> triangulation;
+  GridGenerator::hyper_cube(triangulation);
+  triangulation.refine_global (1);
+
+  hp::FECollection<dim> fe;
+  for (unsigned int i=1; i<=GeometryInfo<dim>::max_children_per_cell; ++i)
+    fe.push_back (FE_Q<dim>(i));
+
+  hp::DoFHandler<dim> dof_handler (triangulation);
+
+  unsigned int index=0;
+  for (typename hp::DoFHandler<dim>::active_cell_iterator
+        cell = dof_handler.begin_active();
+       cell != dof_handler.end(); ++cell, ++index)
+    cell->set_active_fe_index (index);
+
+  dof_handler.distribute_dofs (fe);
+
+  const unsigned int N = dof_handler.n_boundary_dofs();
+  deallog << N << std::endl;
+}
+
+
+int main()
+{
+  std::ofstream logfile("n_boundary_dofs/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  deallog.push("1d");
+  test<1>();
+  deallog.pop();
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/hp/n_boundary_dofs/cmp/generic b/tests/hp/n_boundary_dofs/cmp/generic
new file mode 100644 (file)
index 0000000..4c6e557
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL:1d::2
+DEAL:2d::20
+DEAL:3d::698

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.