Add tests for std::complex make_periodicity_constraints()

diff --git a/tests/bits/periodicity_08.cc b/tests/bits/periodicity_08.cc
new file mode 100644 (file)
index 0000000..0101fed
--- /dev/null
+++ b/tests/bits/periodicity_08.cc
@@ -0,0 +1,74 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check periodic boundary conditions for a simple enough case where we know
+// the exact set of constraints
+//
+// this test simply uses two hypercubes and matches the faces at the far ends
+// and the AffineConstraints is std::complex<double>
+
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include "../tests.h"
+
+
+
+template <int dim>
+void
+test()
+{
+  // create a 2x1 (or 2x1x1) mesh
+  Triangulation<dim>        triangulation;
+  std::vector<unsigned int> repetitions(dim, 1);
+  repetitions[0] = 2;
+  GridGenerator::subdivided_hyper_rectangle(triangulation,
+                                            repetitions,
+                                            Point<dim>(),
+                                            (dim == 2 ? Point<dim>(2, 1) :
+                                                        Point<dim>(2, 1, 1)));
+
+  FE_Q<dim>       fe(1);
+  DoFHandler<dim> dof_handler(triangulation);
+  dof_handler.distribute_dofs(fe);
+
+  AffineConstraints<std::complex<double>> cm;
+  DoFTools::make_periodicity_constraints(dof_handler.begin(0)->face(0),
+                                         (++dof_handler.begin(0))->face(1),
+                                         cm);
+  cm.print(deallog.get_file_stream());
+}
+
+
+
+int
+main()
+{
+  initlog();
+
+  test<2>();
+  test<3>();
+  return 0;
+}

diff --git a/tests/bits/periodicity_08.output b/tests/bits/periodicity_08.output
new file mode 100644 (file)
index 0000000..f4ef592
--- /dev/null
+++ b/tests/bits/periodicity_08.output
@@ -0,0 +1,7 @@
+
+    4 0:  (1.00000,0.00000)
+    5 2:  (1.00000,0.00000)
+    8 0:  (1.00000,0.00000)
+    9 2:  (1.00000,0.00000)
+    10 4:  (1.00000,0.00000)
+    11 6:  (1.00000,0.00000)

diff --git a/tests/dofs/dof_tools_23.cc b/tests/dofs/dof_tools_23.cc
new file mode 100644 (file)
index 0000000..767b62e
--- /dev/null
+++ b/tests/dofs/dof_tools_23.cc
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/function_lib.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/lac/affine_constraints.h>
+#include <deal.II/lac/vector.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include <sstream>
+
+#include "../tests.h"
+#include "dof_tools_periodic.h"
+
+// A simple test for
+//   DoFTools::
+//   make_periodicity_constraints (const FaceIterator       &,
+//                                 const FaceIterator       &,
+//                                 dealii::AffineConstraints<double> &,
+//                                 const std::vector<bool>  &)
+//
+// We project an already periodic function onto the FE space of
+// periodic functions and store the resulting L2 difference
+// between the projected and unprojected function.
+// This should reveal any errors in the constraint_matrix.
+// AffineConstraints is std::complex.
+
+
+
+template <typename DoFHandlerType>
+void
+check_this(const DoFHandlerType &dof_handler)
+{
+  Functions::CosineFunction<DoFHandlerType::dimension> test_func(
+    dof_handler.get_fe().n_components());
+
+  AffineConstraints<std::complex<double>> cm;
+
+  // Apply periodic boundary conditions only in the one direction where
+  // we can match the (locally refined) faces:
+  DoFTools::make_periodicity_constraints(dof_handler.begin(0)->face(0),
+                                         dof_handler.begin(0)->face(1),
+                                         cm);
+  cm.close();
+
+  deallog << cm.n_constraints() << std::endl;
+  deallog << cm.max_constraint_indirections() << std::endl;
+
+  QGauss<DoFHandlerType::dimension> quadrature(6);
+
+  Vector<double>               unconstrained(dof_handler.n_dofs());
+  Vector<std::complex<double>> unconstrained_complex(dof_handler.n_dofs());
+  Vector<std::complex<double>> constrained(dof_handler.n_dofs());
+
+  // Ensure that we can interpolate:
+  if (dof_handler.get_fe().get_unit_support_points().size() == 0)
+    return;
+
+  VectorTools::interpolate(dof_handler, test_func, unconstrained);
+
+  constrained = unconstrained;
+  cm.distribute(constrained);
+
+  // Cast Vector<double> to Vector<std::complex<double>>
+  for (unsigned int index = 0; index < unconstrained.size(); index++)
+    {
+      unconstrained_complex[index] = unconstrained[index];
+    }
+
+  constrained -= unconstrained_complex;
+
+  const double p_l2_error = constrained.l2_norm();
+
+  Assert(p_l2_error < 1e-11, ExcInternalError());
+
+  deallog << "L2_Error : " << p_l2_error << std::endl;
+}

diff --git a/tests/dofs/dof_tools_23.output b/tests/dofs/dof_tools_23.output
new file mode 100644 (file)
index 0000000..7ac94c5
--- /dev/null
+++ b/tests/dofs/dof_tools_23.output
@@ -0,0 +1,159 @@
+
+DEAL::Checking Q1 in 2d:
+DEAL::7
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking Q1 in 3d:
+DEAL::41
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking Q2 in 2d:
+DEAL::13
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking Q2 in 3d:
+DEAL::141
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking Q3 in 2d:
+DEAL::19
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking Q3 in 3d:
+DEAL::297
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking DGQ0 in 2d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGQ0 in 3d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGQ1 in 2d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGQ1 in 3d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGQ3 in 2d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGQ3 in 3d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking DGP0 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking DGP0 in 3d:
+DEAL::0
+DEAL::0
+DEAL::Checking DGP1 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking DGP1 in 3d:
+DEAL::0
+DEAL::0
+DEAL::Checking DGP3 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking DGP3 in 3d:
+DEAL::0
+DEAL::0
+DEAL::Checking Nedelec0 in 2d:
+DEAL::6
+DEAL::1
+DEAL::Checking Nedelec0 in 3d:
+DEAL::72
+DEAL::1
+DEAL::Checking RaviartThomas0 in 2d:
+DEAL::6
+DEAL::1
+DEAL::Checking RaviartThomas1 in 2d:
+DEAL::12
+DEAL::1
+DEAL::Checking RaviartThomas2 in 2d:
+DEAL::18
+DEAL::1
+DEAL::Checking RaviartThomasNodal0 in 2d:
+DEAL::6
+DEAL::1
+DEAL::Checking RaviartThomasNodal1 in 2d:
+DEAL::12
+DEAL::1
+DEAL::Checking RaviartThomasNodal2 in 2d:
+DEAL::18
+DEAL::1
+DEAL::Checking RaviartThomasNodal0 in 3d:
+DEAL::28
+DEAL::1
+DEAL::Checking RaviartThomasNodal1 in 3d:
+DEAL::112
+DEAL::1
+DEAL::Checking RaviartThomasNodal2 in 3d:
+DEAL::252
+DEAL::1
+DEAL::Checking FE_Q<2>(1)3 in 2d:
+DEAL::21
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGQ<2>(2)2 in 2d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGP<2>(3)1 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking FE_Q<3>(1)3 in 3d:
+DEAL::123
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGQ<3>(2)2 in 3d:
+DEAL::0
+DEAL::0
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGP<3>(3)1 in 3d:
+DEAL::0
+DEAL::0
+DEAL::Checking FE_Q<2>(1)3FE_DGQ<2>(2)2 in 2d:
+DEAL::21
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGQ<2>(2)2FE_DGP<2>(3)1 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking FE_DGP<2>(3)1FE_DGQ<2>(2)2 in 2d:
+DEAL::0
+DEAL::0
+DEAL::Checking FE_Q<2>(1)3FE_DGP<2>(3)1FE_Q<2>(1)3 in 2d:
+DEAL::42
+DEAL::1
+DEAL::Checking FE_DGQ<2>(2)2FE_DGQ<2>(2)2FE_Q<2>(3)3 in 2d:
+DEAL::57
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGP<2>(3)1FE_DGP<2>(3)1FE_Q<2>(2)3 in 2d:
+DEAL::39
+DEAL::1
+DEAL::Checking FE_DGQ<3>(1)3FE_DGP<3>(3)1FE_Q<3>(1)3 in 3d:
+DEAL::123
+DEAL::1
+DEAL::Checking (FESystem<2>(FE_Q<2>(1),3))3FE_DGQ<2>(0)1FE_Q<2>(1)3 in 2d:
+DEAL::84
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGQ<2>(3)1FESystem<2>(FE_DGQ<2>(0),3)1FESystem<2>(FE_Q<2>(2),1, FE_DGQ<2>(0),1)2 in 2d:
+DEAL::26
+DEAL::1
+DEAL::L2_Error : 0
+DEAL::Checking FE_DGQ<2>(3)1FE_Nedelec<2>(0)2 in 2d:
+DEAL::12
+DEAL::1
+DEAL::Checking FE_Nedelec<2>(0)1FESystem<2>(FE_DGQ<2>(1),2)1FESystem<2>(FE_Q<2>(2),1, FE_Nedelec<2>(0),2)2 in 2d:
+DEAL::56
+DEAL::1