]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests
authorBruno Turcksin <bruno.turcksin@gmail.com>
Mon, 6 Aug 2018 13:42:14 +0000 (09:42 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Fri, 10 Aug 2018 14:08:18 +0000 (10:08 -0400)
16 files changed:
tests/all-headers/CMakeLists.txt
tests/cuda/create_mesh.h [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_01.output
tests/cuda/matrix_free_matrix_vector_02.cu [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_02.output [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_03.cu [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_03.output [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_04.output
tests/cuda/matrix_free_matrix_vector_05.cu [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_05.output [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_06.cu [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_06.output [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_09.cu [new file with mode: 0644]
tests/cuda/matrix_free_matrix_vector_09.output [new file with mode: 0644]
tests/cuda/matrix_vector_common.h
tests/cuda/matrix_vector_mf.h

index 1d84c7f4c5dc530eab320960acad994660f28c54..1fd8895ee3c0d9267f746939700775b23d7db87b 100644 (file)
@@ -51,6 +51,7 @@ FILE(GLOB_RECURSE _headers RELATIVE ${_include_dir}/deal.II
 #
 LIST(REMOVE_ITEM _headers "lac/cuda_atomic.h"
   "matrix_free/cuda_fe_evaluation.h"
+  "matrix_free/cuda_hanging_nodes_internal.h"
   "matrix_free/cuda_matrix_free.templates.h"
   "matrix_free/cuda_tensor_product_kernels.h")
 
diff --git a/tests/cuda/create_mesh.h b/tests/cuda/create_mesh.h
new file mode 100644 (file)
index 0000000..d22c71d
--- /dev/null
@@ -0,0 +1,132 @@
+//-----------------------  create_mesh.h  -----------------------------
+//    Version: $Name$
+//
+//-----------------------  create_mesh.h  -----------------------------
+
+
+// this creates a mesh that contains cells of all different kinds detected in
+// the MatrixFree class: Use a mesh that consists of a square cell, then a
+// triangular cell, a parallelepiped cell and two trapezoidal cells, where one
+// is very close to being Cartesian (by 1e-8).
+
+#include <deal.II/base/point.h>
+
+#include <deal.II/grid/tria.h>
+
+#include <fstream>
+#include <iostream>
+
+
+void create_mesh(Triangulation<2> &tria, const double scale_grid = 1.)
+{
+  const unsigned int      dim = 2;
+  std::vector<Point<dim>> points(12);
+
+  // build the mesh layer by layer from points
+
+  // 1. cube cell
+  points[0] = Point<dim>(0, 0);
+  points[1] = Point<dim>(0, 1);
+  points[2] = Point<dim>(1, 0);
+  points[3] = Point<dim>(1, 1);
+
+  // 2. rectangular cell
+  points[4] = Point<dim>(3., 0);
+  points[5] = Point<dim>(3., 1);
+
+  // 3. parallelogram cell
+  points[6] = Point<dim>(5., 1.);
+  points[7] = Point<dim>(5., 2.);
+
+  // almost square cell (but trapezoidal by
+  // 1e-8)
+  points[8] = Point<dim>(6., 1.);
+  points[9] = Point<dim>(6., 2. + 1e-8);
+
+  // apparently trapezoidal cell
+  points[10] = Point<dim>(7., 1.4);
+  points[11] = Point<dim>(7.5, numbers::PI);
+
+  if (scale_grid != 1.)
+    for (unsigned int i = 0; i < points.size(); ++i)
+      points[i] *= scale_grid;
+
+
+  // connect the points to cells
+  std::vector<CellData<dim>> cells(5);
+  for (unsigned int i = 0; i < 5; ++i)
+    {
+      cells[i].vertices[0] = 0 + 2 * i;
+      cells[i].vertices[1] = 2 + 2 * i;
+      cells[i].vertices[2] = 1 + 2 * i;
+      cells[i].vertices[3] = 3 + 2 * i;
+      cells[i].material_id = 0;
+    }
+
+  tria.create_triangulation(points, cells, SubCellData());
+}
+
+
+
+void create_mesh(Triangulation<3> &tria, const double scale_grid = 1.)
+{
+  const unsigned int      dim = 3;
+  std::vector<Point<dim>> points(24);
+
+  // build the mesh layer by layer from points
+
+  // 1. cube cell
+  points[0] = Point<dim>(0, 0, 0);
+  points[1] = Point<dim>(0, 1., 0);
+  points[2] = Point<dim>(0, 0, 1);
+  points[3] = Point<dim>(0, 1., 1);
+  points[4] = Point<dim>(1., 0, 0);
+  points[5] = Point<dim>(1., 1., 0);
+  points[6] = Point<dim>(1., 0, 1);
+  points[7] = Point<dim>(1., 1., 1);
+
+  // 2. rectangular cell
+  points[8]  = Point<dim>(3., 0, 0);
+  points[9]  = Point<dim>(3., 1, 0);
+  points[10] = Point<dim>(3., 0, 1);
+  points[11] = Point<dim>(3., 1, 1);
+
+  // 3. parallelogram cell
+  points[12] = Point<dim>(5., 1., 1.);
+  points[13] = Point<dim>(5., 2., 1.);
+  points[14] = Point<dim>(5., 1., 2.);
+  points[15] = Point<dim>(5., 2., 2.);
+
+  // almost square cell (but trapezoidal by
+  // 1e-8 in y-direction)
+  points[16] = Point<dim>(6., 1., 1.);
+  points[17] = Point<dim>(6., 2. + 1e-8, 1.);
+  points[18] = Point<dim>(6., 1., 2.);
+  points[19] = Point<dim>(6., 2., 2.);
+
+  // apparently trapezoidal cell
+  points[20] = Point<dim>(7., 1.4, 1.2231);
+  points[21] = Point<dim>(7.5, numbers::PI, 1.334);
+  points[22] = Point<dim>(7., 1.5, 7.1);
+  points[23] = Point<dim>(7.5, 3.8, 2.99);
+
+  if (scale_grid != 1.)
+    for (unsigned int i = 0; i < points.size(); ++i)
+      points[i] *= scale_grid;
+
+  // connect the points to cells
+  std::vector<CellData<dim>> cells(5);
+  for (unsigned int i = 0; i < 5; ++i)
+    {
+      cells[i].vertices[0] = 0 + 4 * i;
+      cells[i].vertices[1] = 4 + 4 * i;
+      cells[i].vertices[2] = 1 + 4 * i;
+      cells[i].vertices[3] = 5 + 4 * i;
+      cells[i].vertices[4] = 2 + 4 * i;
+      cells[i].vertices[5] = 6 + 4 * i;
+      cells[i].vertices[6] = 3 + 4 * i;
+      cells[i].vertices[7] = 7 + 4 * i;
+      cells[i].material_id = 0;
+    }
+  tria.create_triangulation(points, cells, SubCellData());
+}
index 43ef13cd754847276775ad0edd4edc6128e9942f..ed73cb7cd98d5fbc36775ec0fef55eede460c983 100644 (file)
@@ -2,9 +2,6 @@
 DEAL:2d::Testing FE_Q<2>(1)
 DEAL:2d::Norm of difference: 0
 DEAL:2d::
-DEAL:2d::Testing FE_Q<2>(2)
-DEAL:2d::Norm of difference: 0
-DEAL:2d::
 DEAL:2d::Testing FE_Q<2>(3)
 DEAL:2d::Norm of difference: 0
 DEAL:2d::
diff --git a/tests/cuda/matrix_free_matrix_vector_02.cu b/tests/cuda/matrix_free_matrix_vector_02.cu
new file mode 100644 (file)
index 0000000..e42e6fb
--- /dev/null
@@ -0,0 +1,49 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.
+//
+// ---------------------------------------------------------------------
+
+// this function tests the correctness of the implementation of matrix free
+// matrix-vector products by comparing with the result of deal.II sparse
+// matrix. The mesh uses a hypercube mesh with no hanging nodes, but with zero
+// Dirichlet conditions.
+
+#include <deal.II/base/function.h>
+
+#include "../tests.h"
+
+std::ofstream logfile("output");
+
+#include "matrix_vector_common.h"
+
+
+template <int dim, int fe_degree>
+void
+test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(5 - dim);
+
+  FE_Q<dim>       fe(fe_degree);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe);
+  ConstraintMatrix constraints;
+  VectorTools::interpolate_boundary_values(dof,
+                                           0,
+                                           Functions::ZeroFunction<dim>(),
+                                           constraints);
+  constraints.close();
+
+  do_test<dim, fe_degree, double, fe_degree + 1>(dof, constraints);
+}
diff --git a/tests/cuda/matrix_free_matrix_vector_02.output b/tests/cuda/matrix_free_matrix_vector_02.output
new file mode 100644 (file)
index 0000000..ed73cb7
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL:2d::Testing FE_Q<2>(1)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:2d::Testing FE_Q<2>(3)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:3d::Testing FE_Q<3>(1)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(2)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(3)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
diff --git a/tests/cuda/matrix_free_matrix_vector_03.cu b/tests/cuda/matrix_free_matrix_vector_03.cu
new file mode 100644 (file)
index 0000000..cba3875
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// this function tests the correctness of the implementation of matrix free
+// matrix-vector products by comparing with the result of deal.II sparse
+// matrix. The mesh uses a hypercube mesh with hanging nodes (created by
+// randomly refining some cells) and zero Dirichlet conditions.
+
+#include <deal.II/base/function.h>
+
+#include "../tests.h"
+
+std::ofstream logfile("output");
+
+#include "matrix_vector_common.h"
+
+template <int dim, int fe_degree>
+void
+test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active(),
+                                                    endc = tria.end();
+  if (dim < 3 || fe_degree < 2)
+    tria.refine_global(2);
+  else
+    tria.refine_global(1);
+  tria.begin(tria.n_levels() - 1)->set_refine_flag();
+  tria.last()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  cell = tria.begin_active();
+  for (unsigned int i = 0; i < 10 - 3 * dim; ++i)
+    {
+      cell                 = tria.begin_active();
+      unsigned int counter = 0;
+      for (; cell != endc; ++cell, ++counter)
+        if (counter % (7 - i) == 0)
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement();
+    }
+
+  FE_Q<dim>       fe(fe_degree);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe);
+  AffineConstraints<double> constraints;
+  DoFTools::make_hanging_node_constraints(dof, constraints);
+
+  VectorTools::interpolate_boundary_values(dof,
+                                           0,
+                                           Functions::ZeroFunction<dim>(),
+                                           constraints);
+  constraints.close();
+
+  // Skip 2D tests with even fe_degree
+  if ((dim == 3) || ((fe_degree % 2) == 1))
+    do_test<dim, fe_degree, double, fe_degree + 1>(dof, constraints);
+}
diff --git a/tests/cuda/matrix_free_matrix_vector_03.output b/tests/cuda/matrix_free_matrix_vector_03.output
new file mode 100644 (file)
index 0000000..ab11bfe
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL:2d::Testing FE_Q<2>(1)
+DEAL:2d::Norm of difference: 0.
+DEAL:2d::
+DEAL:2d::Testing FE_Q<2>(3)
+DEAL:2d::Norm of difference: 0.
+DEAL:2d::
+DEAL:3d::Testing FE_Q<3>(1)
+DEAL:3d::Norm of difference: 0.
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(2)
+DEAL:3d::Norm of difference: 0.
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(3)
+DEAL:3d::Norm of difference: 0.
+DEAL:3d::
index 43ef13cd754847276775ad0edd4edc6128e9942f..ed73cb7cd98d5fbc36775ec0fef55eede460c983 100644 (file)
@@ -2,9 +2,6 @@
 DEAL:2d::Testing FE_Q<2>(1)
 DEAL:2d::Norm of difference: 0
 DEAL:2d::
-DEAL:2d::Testing FE_Q<2>(2)
-DEAL:2d::Norm of difference: 0
-DEAL:2d::
 DEAL:2d::Testing FE_Q<2>(3)
 DEAL:2d::Norm of difference: 0
 DEAL:2d::
diff --git a/tests/cuda/matrix_free_matrix_vector_05.cu b/tests/cuda/matrix_free_matrix_vector_05.cu
new file mode 100644 (file)
index 0000000..cdc08ea
--- /dev/null
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// this function tests the correctness of the implementation of matrix free
+// matrix-vector products by comparing with the result of deal.II sparse
+// matrix. The mesh uses a parallelogram mesh with hanging nodes (only cell
+// type: 1 = linear).
+
+#include "../tests.h"
+
+std::ofstream logfile("output");
+
+#include "matrix_vector_common.h"
+
+
+template <int dim, int fe_degree>
+void
+test()
+{
+  Triangulation<dim> tria;
+  Point<dim>         points[dim];
+  points[0][0] = 0.25;
+  points[0][1] = 0.123;
+  points[1][0] = 0.09983712334;
+  points[1][1] = 0.314159265358979;
+  if (dim == 3)
+    {
+      points[2][0] = 0.21;
+      points[2][2] = 0.4123;
+    }
+  GridGenerator::parallelepiped(tria, points);
+  typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active(),
+                                                    endc = tria.end();
+  for (; cell != endc; ++cell)
+    if (cell->center().norm() < 1e-8)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  cell = tria.begin_active();
+  for (; cell != endc; ++cell)
+    if (cell->center().norm() < 0.2)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  if (dim < 3)
+    tria.refine_global(2);
+  tria.begin(tria.n_levels() - 1)->set_refine_flag();
+  tria.last()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  for (int i = 0; i < 7 - 2 * fe_degree; ++i)
+    {
+      cell                 = tria.begin_active();
+      unsigned int counter = 0;
+      for (; cell != endc; ++cell, ++counter)
+        if (counter % (7 - i) == 0)
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement();
+    }
+
+  FE_Q<dim>       fe(fe_degree);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe);
+  AffineConstraints<double> constraints;
+  DoFTools::make_hanging_node_constraints(dof, constraints);
+  constraints.close();
+
+  // Skip 2D tests with even fe_degree
+  if ((dim == 3) || ((fe_degree % 2) == 1))
+    do_test<dim, fe_degree, double, fe_degree + 1>(dof, constraints);
+}
diff --git a/tests/cuda/matrix_free_matrix_vector_05.output b/tests/cuda/matrix_free_matrix_vector_05.output
new file mode 100644 (file)
index 0000000..ed73cb7
--- /dev/null
@@ -0,0 +1,16 @@
+
+DEAL:2d::Testing FE_Q<2>(1)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:2d::Testing FE_Q<2>(3)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:3d::Testing FE_Q<3>(1)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(2)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
+DEAL:3d::Testing FE_Q<3>(3)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
diff --git a/tests/cuda/matrix_free_matrix_vector_06.cu b/tests/cuda/matrix_free_matrix_vector_06.cu
new file mode 100644 (file)
index 0000000..55cf8de
--- /dev/null
@@ -0,0 +1,78 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// this function tests the correctness of the implementation of matrix free
+// matrix-vector products by comparing with the result of deal.II sparse
+// matrix. The mesh uses a mesh consisting of several different cell types
+// according to the create_mesh helper function. Quite large mesh
+// mesh so that the thread parallelization is actually used
+
+#include <deal.II/base/function.h>
+
+#include "../tests.h"
+#include "create_mesh.h"
+
+std::ofstream logfile("output");
+
+#include "matrix_vector_common.h"
+
+template <int dim, int fe_degree>
+void
+test()
+{
+  if (fe_degree > 1)
+    return;
+  Triangulation<dim> tria;
+  create_mesh(tria);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  typename Triangulation<dim>::active_cell_iterator cell, endc;
+  cell = tria.begin_active();
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+    if (cell->center().norm() < 0.5)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  tria.begin(tria.n_levels() - 1)->set_refine_flag();
+  tria.last()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  tria.refine_global(1);
+  cell = tria.begin_active();
+  for (unsigned int i = 0; i < 10 - 3 * dim; ++i)
+    {
+      cell                 = tria.begin_active();
+      endc                 = tria.end();
+      unsigned int counter = 0;
+      for (; cell != endc; ++cell, ++counter)
+        if (counter % (7 - i) == 0)
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement();
+    }
+
+  FE_Q<dim>       fe(fe_degree);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe);
+  AffineConstraints<double> constraints;
+  DoFTools::make_hanging_node_constraints(dof, constraints);
+  VectorTools::interpolate_boundary_values(dof,
+                                           0,
+                                           Functions::ZeroFunction<dim>(),
+                                           constraints);
+  constraints.close();
+
+  do_test<dim, fe_degree, double, fe_degree + 1>(dof, constraints);
+}
diff --git a/tests/cuda/matrix_free_matrix_vector_06.output b/tests/cuda/matrix_free_matrix_vector_06.output
new file mode 100644 (file)
index 0000000..e9e176a
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL:2d::Testing FE_Q<2>(1)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:3d::Testing FE_Q<3>(1)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
diff --git a/tests/cuda/matrix_free_matrix_vector_09.cu b/tests/cuda/matrix_free_matrix_vector_09.cu
new file mode 100644 (file)
index 0000000..5261136
--- /dev/null
@@ -0,0 +1,76 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+
+// same test as matrix_vector_06 (quite large mesh, hanging nodes, different
+// cell types), but tiny domain of size 1e-20 to test correctness of
+// relative scaling in mapping info
+
+#include <deal.II/base/function.h>
+
+#include "../tests.h"
+#include "create_mesh.h"
+
+std::ofstream logfile("output");
+
+#include "matrix_vector_common.h"
+
+template <int dim, int fe_degree>
+void
+test()
+{
+  if (fe_degree > 1)
+    return;
+
+  Triangulation<dim> tria;
+  create_mesh(tria, 1e-20);
+  tria.begin_active()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  typename Triangulation<dim>::active_cell_iterator cell, endc;
+  cell = tria.begin_active();
+  endc = tria.end();
+  for (; cell != endc; ++cell)
+    if (cell->center().norm() < 0.5 * 1e-20)
+      cell->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  tria.begin(tria.n_levels() - 1)->set_refine_flag();
+  tria.last()->set_refine_flag();
+  tria.execute_coarsening_and_refinement();
+  cell = tria.begin_active();
+  for (unsigned int i = 0; i < 10 - 3 * dim; ++i)
+    {
+      cell                 = tria.begin_active();
+      endc                 = tria.end();
+      unsigned int counter = 0;
+      for (; cell != endc; ++cell, ++counter)
+        if (counter % (7 - i) == 0)
+          cell->set_refine_flag();
+      tria.execute_coarsening_and_refinement();
+    }
+
+  FE_Q<dim>       fe(fe_degree);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe);
+  AffineConstraints<double> constraints;
+  DoFTools::make_hanging_node_constraints(dof, constraints);
+  VectorTools::interpolate_boundary_values(dof,
+                                           0,
+                                           Functions::ZeroFunction<dim>(),
+                                           constraints);
+  constraints.close();
+
+  do_test<dim, fe_degree, double, fe_degree + 1>(dof, constraints);
+}
diff --git a/tests/cuda/matrix_free_matrix_vector_09.output b/tests/cuda/matrix_free_matrix_vector_09.output
new file mode 100644 (file)
index 0000000..e9e176a
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL:2d::Testing FE_Q<2>(1)
+DEAL:2d::Norm of difference: 0
+DEAL:2d::
+DEAL:3d::Testing FE_Q<3>(1)
+DEAL:3d::Norm of difference: 0
+DEAL:3d::
index 1fe3f76aacbbc315188de1d267d139554a3ca782..a0f673a844628c820ccbb1611df5a2acc7cb19fd 100644 (file)
@@ -46,7 +46,7 @@
 #include "matrix_vector_mf.h"
 
 
-// forward declare this function. will be implemented in .cc files
+// forward declare this function. will be implemented in .cu files
 template <int dim, int fe_degree>
 void
 test();
@@ -90,7 +90,6 @@ do_test(const DoFHandler<dim> &          dof,
   cudaDeviceSynchronize();
   out.import(out_device, VectorOperation::insert);
 
-
   // assemble sparse matrix with (\nabla v, \nabla u) + (v, 10 * u)
   SparsityPattern sparsity;
   {
@@ -142,14 +141,11 @@ do_test(const DoFHandler<dim> &          dof,
     in_host[i] = in[i];
   sparse_matrix.vmult(out_host, in_host);
 
-  Number out_dist_cpu_norm = 0.;
-  Number out_norm          = 0.;
+  Number out_norm = 0.;
   for (unsigned i = 0; i < n_dofs; ++i)
-    {
-      out_norm += std::pow(out[i] - out_host[i], 2);
-      out_dist_cpu_norm += std::pow(out_host[i], 2);
-    }
-  const double diff_norm = out_norm / out_dist_cpu_norm;
+    out_norm += std::pow(out[i] - out_host[i], 2);
+  const double diff_norm = std::sqrt(out_norm) / out_host.linfty_norm();
+
   deallog << "Norm of difference: " << diff_norm << std::endl << std::endl;
 }
 
@@ -166,7 +162,7 @@ main()
   {
     deallog.push("2d");
     test<2, 1>();
-    test<2, 2>();
+    //  test<2, 2>();
     test<2, 3>();
     deallog.pop();
     deallog.push("3d");
index 39f7ac0fce7dcf05435528ef744ac152b058a685..e175d3791b00afb98737d81a536e417452426a12 100644 (file)
@@ -120,4 +120,5 @@ MatrixFreeTest<dim, fe_degree, Number, n_q_points_1d>::vmult(
   dst = static_cast<Number>(0.);
   HelmholtzOperator<dim, fe_degree, Number, n_q_points_1d> helmholtz_operator;
   data.cell_loop(helmholtz_operator, src, dst);
+  data.copy_constrained_values(src, dst);
 }

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.