From: Sebastian Kinnewig <kinnewig@ifam.uni-hannover.de>
Date: Mon, 15 Jan 2024 10:43:30 +0000 (+0100)
Subject: Added new test case based on fe_nedelec_sz_gradient_divergence_theorem.
X-Git-Tag: relicensing~125^2~1
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9d7a3e09378117f1d37dff50bc23bff128d12792;p=dealii.git

Added new test case based on fe_nedelec_sz_gradient_divergence_theorem.
---

diff --git a/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.cc b/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.cc
new file mode 100644
index 0000000000..d4e79a5e67
--- /dev/null
+++ b/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.cc
@@ -0,0 +1,195 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2023 - 2024 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 test is a variation from
+// fe_nedelec_sz_gradient_divergence_theorem.cc.
+// Here, we consider hanging nodes additionally.
+
+#include <deal.II/base/function.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_nedelec_sz.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_q1.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/manifold_lib.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <sstream>
+
+#include "../tests.h"
+
+template <int dim>
+Tensor<1, dim>
+ones()
+{
+  Tensor<1, dim> result;
+  for (unsigned int i = 0; i < dim; ++i)
+    result[i] = 1.0;
+  return result;
+}
+
+template <int dim>
+void
+test(const Triangulation<dim> &tr,
+     const FiniteElement<dim> &fe,
+     const double              tolerance)
+{
+  DoFHandler<dim> dof(tr);
+  dof.distribute_dofs(fe);
+
+  std::stringstream ss;
+
+  deallog << "FE=" << fe.get_name() << std::endl;
+
+  const QGauss<dim> quadrature(4);
+  FEValues<dim>     fe_values(fe,
+                          quadrature,
+                          update_gradients | update_quadrature_points |
+                            update_JxW_values);
+
+  const QGauss<dim - 1> face_quadrature(4);
+  FEFaceValues<dim>     fe_face_values(fe,
+                                   face_quadrature,
+                                   update_values | update_quadrature_points |
+                                     update_normal_vectors | update_JxW_values);
+
+  for (typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active();
+       cell != dof.end();
+       ++cell)
+    {
+      fe_values.reinit(cell);
+
+      deallog << "Cell nodes:" << std::endl;
+      for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
+        {
+          deallog << i << ": ( ";
+          for (unsigned int d = 0; d < dim; ++d)
+            deallog << cell->vertex(i)[d] << ' ';
+          deallog << ')' << std::endl;
+        }
+
+      bool cell_ok = true;
+
+      for (unsigned int c = 0; c < fe.n_components(); ++c)
+        {
+          FEValuesExtractors::Scalar single_component(c);
+
+          for (unsigned int i = 0; i < fe_values.dofs_per_cell; ++i)
+            {
+              ss << "component=" << c << ", dof=" << i << std::endl;
+
+              Tensor<1, dim> bulk_integral;
+              for (const auto q : fe_values.quadrature_point_indices())
+                {
+                  bulk_integral += fe_values[single_component].gradient(i, q) *
+                                   fe_values.JxW(q);
+                }
+
+              Tensor<1, dim> boundary_integral;
+              for (const unsigned int face : GeometryInfo<dim>::face_indices())
+                {
+                  fe_face_values.reinit(cell, face);
+                  for (const auto q : fe_face_values.quadrature_point_indices())
+                    {
+                      boundary_integral +=
+                        fe_face_values[single_component].value(i, q) *
+                        fe_face_values.normal_vector(q) * fe_face_values.JxW(q);
+                    }
+                }
+
+              if ((bulk_integral - boundary_integral).norm_square() >
+                  tolerance * (bulk_integral.norm() + boundary_integral.norm()))
+                {
+                  deallog << "Failed:" << std::endl;
+                  deallog << ss.str() << std::endl;
+                  deallog << "    bulk integral=" << bulk_integral << std::endl;
+                  deallog << "boundary integral=" << boundary_integral
+                          << std::endl;
+                  deallog
+                    << "Error! difference between bulk and surface integrals is greater than "
+                    << tolerance << "!\n\n"
+                    << std::endl;
+                  cell_ok = false;
+                }
+
+              ss.str("");
+            }
+        }
+
+      deallog << (cell_ok ? "OK: cell bulk and boundary integrals match...\n" :
+                            "Failed divergence test...\n")
+              << std::endl;
+    }
+}
+
+
+
+template <int dim>
+void
+test_hyper_ball(const double tolerance)
+{
+  Triangulation<dim> tr;
+  GridGenerator::hyper_ball(tr);
+
+  static const SphericalManifold<dim> boundary;
+  tr.set_manifold(0, boundary);
+
+  // Refine one cell at the boundary,
+  // such that the resulting grid contains
+  // hanging nodes.
+  for (auto &cell : tr.cell_iterators())
+    {
+      if (!cell->is_active())
+        continue;
+
+      if (cell->at_boundary())
+        {
+          cell->set_refine_flag();
+          break;
+        }
+    }
+
+  tr.execute_coarsening_and_refinement();
+
+
+  //  tr.refine_global(1); // taking too long,
+
+  FE_NedelecSZ<dim> fe(1);
+  test(tr, fe, tolerance);
+}
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(8);
+
+  test_hyper_ball<2>(1e-6);
+  test_hyper_ball<3>(1e-6);
+
+  deallog << "done..." << std::endl;
+}
diff --git a/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.output b/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.output
new file mode 100644
index 0000000000..23d544af09
--- /dev/null
+++ b/tests/fe/fe_nedelec_sz_divergence_theorem_hanging_nodes.output
@@ -0,0 +1,214 @@
+
+DEAL::FE=FE_NedelecSZ<2>(1)
+DEAL::Cell nodes:
+DEAL::0: ( -0.70710678 -0.70710678 )
+DEAL::1: ( -0.29289322 -0.29289322 )
+DEAL::2: ( -0.70710678 0.70710678 )
+DEAL::3: ( -0.29289322 0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.29289322 -0.29289322 )
+DEAL::1: ( 0.29289322 -0.29289322 )
+DEAL::2: ( -0.29289322 0.29289322 )
+DEAL::3: ( 0.29289322 0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.70710678 -0.70710678 )
+DEAL::1: ( 0.70710678 0.70710678 )
+DEAL::2: ( 0.29289322 -0.29289322 )
+DEAL::3: ( 0.29289322 0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.70710678 0.70710678 )
+DEAL::1: ( -0.29289322 0.29289322 )
+DEAL::2: ( 0.70710678 0.70710678 )
+DEAL::3: ( 0.29289322 0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.70710678 -0.70710678 )
+DEAL::1: ( -1.1102230e-16 -1.0000000 )
+DEAL::2: ( -0.50000000 -0.50000000 )
+DEAL::3: ( -5.5511151e-17 -0.64644661 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -1.1102230e-16 -1.0000000 )
+DEAL::1: ( 0.70710678 -0.70710678 )
+DEAL::2: ( -5.5511151e-17 -0.64644661 )
+DEAL::3: ( 0.50000000 -0.50000000 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.50000000 -0.50000000 )
+DEAL::1: ( -5.5511151e-17 -0.64644661 )
+DEAL::2: ( -0.29289322 -0.29289322 )
+DEAL::3: ( 0.0000000 -0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -5.5511151e-17 -0.64644661 )
+DEAL::1: ( 0.50000000 -0.50000000 )
+DEAL::2: ( 0.0000000 -0.29289322 )
+DEAL::3: ( 0.29289322 -0.29289322 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::FE=FE_NedelecSZ<3>(1)
+DEAL::Cell nodes:
+DEAL::0: ( -0.21132487 -0.21132487 -0.21132487 )
+DEAL::1: ( 0.21132487 -0.21132487 -0.21132487 )
+DEAL::2: ( -0.21132487 0.21132487 -0.21132487 )
+DEAL::3: ( 0.21132487 0.21132487 -0.21132487 )
+DEAL::4: ( -0.21132487 -0.21132487 0.21132487 )
+DEAL::5: ( 0.21132487 -0.21132487 0.21132487 )
+DEAL::6: ( -0.21132487 0.21132487 0.21132487 )
+DEAL::7: ( 0.21132487 0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.57735027 -0.57735027 -0.57735027 )
+DEAL::1: ( 0.57735027 0.57735027 -0.57735027 )
+DEAL::2: ( 0.21132487 -0.21132487 -0.21132487 )
+DEAL::3: ( 0.21132487 0.21132487 -0.21132487 )
+DEAL::4: ( 0.57735027 -0.57735027 0.57735027 )
+DEAL::5: ( 0.57735027 0.57735027 0.57735027 )
+DEAL::6: ( 0.21132487 -0.21132487 0.21132487 )
+DEAL::7: ( 0.21132487 0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.57735027 -0.57735027 0.57735027 )
+DEAL::1: ( 0.57735027 -0.57735027 0.57735027 )
+DEAL::2: ( -0.21132487 -0.21132487 0.21132487 )
+DEAL::3: ( 0.21132487 -0.21132487 0.21132487 )
+DEAL::4: ( -0.57735027 0.57735027 0.57735027 )
+DEAL::5: ( 0.57735027 0.57735027 0.57735027 )
+DEAL::6: ( -0.21132487 0.21132487 0.21132487 )
+DEAL::7: ( 0.21132487 0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 )
+DEAL::1: ( -0.21132487 -0.21132487 -0.21132487 )
+DEAL::2: ( -0.57735027 0.57735027 -0.57735027 )
+DEAL::3: ( -0.21132487 0.21132487 -0.21132487 )
+DEAL::4: ( -0.57735027 -0.57735027 0.57735027 )
+DEAL::5: ( -0.21132487 -0.21132487 0.21132487 )
+DEAL::6: ( -0.57735027 0.57735027 0.57735027 )
+DEAL::7: ( -0.21132487 0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 )
+DEAL::1: ( 0.57735027 -0.57735027 -0.57735027 )
+DEAL::2: ( -0.21132487 -0.21132487 -0.21132487 )
+DEAL::3: ( 0.21132487 -0.21132487 -0.21132487 )
+DEAL::4: ( -0.57735027 -0.57735027 0.57735027 )
+DEAL::5: ( 0.57735027 -0.57735027 0.57735027 )
+DEAL::6: ( -0.21132487 -0.21132487 0.21132487 )
+DEAL::7: ( 0.21132487 -0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.57735027 0.57735027 -0.57735027 )
+DEAL::1: ( -0.21132487 0.21132487 -0.21132487 )
+DEAL::2: ( 0.57735027 0.57735027 -0.57735027 )
+DEAL::3: ( 0.21132487 0.21132487 -0.21132487 )
+DEAL::4: ( -0.57735027 0.57735027 0.57735027 )
+DEAL::5: ( -0.21132487 0.21132487 0.21132487 )
+DEAL::6: ( 0.57735027 0.57735027 0.57735027 )
+DEAL::7: ( 0.21132487 0.21132487 0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.57735027 -0.57735027 -0.57735027 )
+DEAL::1: ( 0.0000000 -0.70710678 -0.70710678 )
+DEAL::2: ( -0.70710678 0.0000000 -0.70710678 )
+DEAL::3: ( 0.0000000 0.0000000 -1.0000000 )
+DEAL::4: ( -0.39433757 -0.39433757 -0.39433757 )
+DEAL::5: ( 0.0000000 -0.45921582 -0.45921582 )
+DEAL::6: ( -0.45921582 0.0000000 -0.45921582 )
+DEAL::7: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.0000000 -0.70710678 -0.70710678 )
+DEAL::1: ( 0.57735027 -0.57735027 -0.57735027 )
+DEAL::2: ( 0.0000000 0.0000000 -1.0000000 )
+DEAL::3: ( 0.70710678 0.0000000 -0.70710678 )
+DEAL::4: ( 0.0000000 -0.45921582 -0.45921582 )
+DEAL::5: ( 0.39433757 -0.39433757 -0.39433757 )
+DEAL::6: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::7: ( 0.45921582 0.0000000 -0.45921582 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.70710678 0.0000000 -0.70710678 )
+DEAL::1: ( 0.0000000 0.0000000 -1.0000000 )
+DEAL::2: ( -0.57735027 0.57735027 -0.57735027 )
+DEAL::3: ( 0.0000000 0.70710678 -0.70710678 )
+DEAL::4: ( -0.45921582 0.0000000 -0.45921582 )
+DEAL::5: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::6: ( -0.39433757 0.39433757 -0.39433757 )
+DEAL::7: ( 0.0000000 0.45921582 -0.45921582 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.0000000 0.0000000 -1.0000000 )
+DEAL::1: ( 0.70710678 0.0000000 -0.70710678 )
+DEAL::2: ( 0.0000000 0.70710678 -0.70710678 )
+DEAL::3: ( 0.57735027 0.57735027 -0.57735027 )
+DEAL::4: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::5: ( 0.45921582 0.0000000 -0.45921582 )
+DEAL::6: ( 0.0000000 0.45921582 -0.45921582 )
+DEAL::7: ( 0.39433757 0.39433757 -0.39433757 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.39433757 -0.39433757 -0.39433757 )
+DEAL::1: ( 0.0000000 -0.45921582 -0.45921582 )
+DEAL::2: ( -0.45921582 0.0000000 -0.45921582 )
+DEAL::3: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::4: ( -0.21132487 -0.21132487 -0.21132487 )
+DEAL::5: ( 0.0000000 -0.21132487 -0.21132487 )
+DEAL::6: ( -0.21132487 0.0000000 -0.21132487 )
+DEAL::7: ( 0.0000000 0.0000000 -0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.0000000 -0.45921582 -0.45921582 )
+DEAL::1: ( 0.39433757 -0.39433757 -0.39433757 )
+DEAL::2: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::3: ( 0.45921582 0.0000000 -0.45921582 )
+DEAL::4: ( 0.0000000 -0.21132487 -0.21132487 )
+DEAL::5: ( 0.21132487 -0.21132487 -0.21132487 )
+DEAL::6: ( 0.0000000 0.0000000 -0.21132487 )
+DEAL::7: ( 0.21132487 0.0000000 -0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( -0.45921582 0.0000000 -0.45921582 )
+DEAL::1: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::2: ( -0.39433757 0.39433757 -0.39433757 )
+DEAL::3: ( 0.0000000 0.45921582 -0.45921582 )
+DEAL::4: ( -0.21132487 0.0000000 -0.21132487 )
+DEAL::5: ( 0.0000000 0.0000000 -0.21132487 )
+DEAL::6: ( -0.21132487 0.21132487 -0.21132487 )
+DEAL::7: ( 0.0000000 0.21132487 -0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::Cell nodes:
+DEAL::0: ( 0.0000000 0.0000000 -0.60566243 )
+DEAL::1: ( 0.45921582 0.0000000 -0.45921582 )
+DEAL::2: ( 0.0000000 0.45921582 -0.45921582 )
+DEAL::3: ( 0.39433757 0.39433757 -0.39433757 )
+DEAL::4: ( 0.0000000 0.0000000 -0.21132487 )
+DEAL::5: ( 0.21132487 0.0000000 -0.21132487 )
+DEAL::6: ( 0.0000000 0.21132487 -0.21132487 )
+DEAL::7: ( 0.21132487 0.21132487 -0.21132487 )
+DEAL::OK: cell bulk and boundary integrals match...
+
+DEAL::done...