From: Jiaqi Zhang <jiaqi2@clemson.edu>
Date: Mon, 11 Apr 2022 15:27:14 +0000 (-0400)
Subject: add a test
X-Git-Tag: v9.4.0-rc1~307^2~4
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=775b96d083d379b7d56fbd4690958e2c4c083cd0;p=dealii.git

add a test
---

diff --git a/tests/feinterface/fe_interface_values_10.cc b/tests/feinterface/fe_interface_values_10.cc
new file mode 100644
index 0000000000..a41c770636
--- /dev/null
+++ b/tests/feinterface/fe_interface_values_10.cc
@@ -0,0 +1,263 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 - 2021 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 FEInterfaceViews::get_*_function_values/gradients...
+
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_renumbering.h>
+
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_interface_values.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_values.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <fstream>
+#include <type_traits>
+
+#include "../tests.h"
+
+template <int dim = 2>
+void
+test_jump_function()
+{
+  FE_DGP<dim> fe(0);
+
+  Triangulation<dim> tria;
+  DoFHandler<dim>    dof_handler(tria);
+  Vector<double>     solution;
+
+  GridGenerator::hyper_cube(tria, -1, 1);
+  tria.refine_global(1);
+  dof_handler.distribute_dofs(fe);
+
+  solution.reinit(dof_handler.n_dofs());
+
+  std::vector<types::global_dof_index> local_dof_indices(1);
+  // Populate with non-trivial values
+  for (const auto &cell : dof_handler.active_cell_iterators())
+    {
+      cell->get_dof_indices(local_dof_indices);
+      solution(local_dof_indices[0]) =
+        static_cast<double>(cell->active_cell_index());
+    }
+
+  const QGauss<dim - 1>  face_quadrature(2);
+  FEInterfaceValues<dim> fe_iv(fe,
+                               face_quadrature,
+                               update_values | update_gradients |
+                                 update_quadrature_points |
+                                 update_3rd_derivatives);
+
+  std::vector<double> qp_jumps_global(face_quadrature.size());
+
+  typename DoFHandler<dim>::active_cell_iterator cell =
+    dof_handler.begin_active();
+  for (; cell != dof_handler.end(); ++cell)
+    for (const auto f : cell->face_indices())
+      if (!cell->face(f)->at_boundary())
+        {
+          {
+            fe_iv.reinit(cell,
+                         f,
+                         numbers::invalid_unsigned_int,
+                         cell->neighbor(f),
+                         cell->neighbor_of_neighbor(f),
+                         numbers::invalid_unsigned_int);
+
+            fe_iv.get_jump_in_function_values(solution, qp_jumps_global);
+
+            double exact = cell->active_cell_index() * 1.0 -
+                           cell->neighbor(f)->active_cell_index() * 1.0;
+
+            for (unsigned int q = 0; q < face_quadrature.size(); ++q)
+              Assert(std::fabs(qp_jumps_global[q] - exact) < 1e-15,
+                     ExcNotImplemented());
+          }
+        }
+
+  deallog << "OK" << std::endl;
+}
+
+template <int dim = 2>
+void
+test()
+{
+  FESystem<dim>              fe(FE_DGP<dim>(4), 1);
+  FEValuesExtractors::Scalar extractor(0);
+
+  Triangulation<dim> tria;
+  DoFHandler<dim>    dof_handler(tria);
+  Vector<double>     solution;
+
+  GridGenerator::hyper_cube(tria, -1, 1);
+  tria.refine_global(1);
+  dof_handler.distribute_dofs(fe);
+
+  solution.reinit(dof_handler.n_dofs());
+
+  const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
+
+  std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
+  // Populate with non-trivial values
+  for (const auto &cell : dof_handler.active_cell_iterators())
+    {
+      cell->get_dof_indices(local_dof_indices);
+      for (unsigned int i = 0; i < dofs_per_cell; ++i)
+        solution(local_dof_indices[i]) =
+          static_cast<double>(cell->active_cell_index() + 1);
+    }
+
+  const QGauss<dim - 1>  face_quadrature(2);
+  FEInterfaceValues<dim> fe_iv(fe,
+                               face_quadrature,
+                               update_values | update_gradients |
+                                 update_quadrature_points |
+                                 update_3rd_derivatives);
+  const unsigned int     n_face_q_points = face_quadrature.size();
+
+  std::vector<double> jump_in_values(n_face_q_points);
+  std::vector<double> exact_jump_in_values(n_face_q_points);
+
+  std::vector<Tensor<1, dim>> jump_in_grads(n_face_q_points);
+  std::vector<Tensor<1, dim>> exact_jump_in_grads(n_face_q_points);
+
+  std::vector<Tensor<2, dim>> jump_in_hessians(n_face_q_points);
+  std::vector<Tensor<2, dim>> exact_jump_in_hessians(n_face_q_points);
+
+  std::vector<Tensor<3, dim>> jump_in_3rd_derivatives(n_face_q_points);
+  std::vector<Tensor<3, dim>> exact_jump_in_3rd_derivatives(n_face_q_points);
+
+  std::vector<double> average_of_values(n_face_q_points);
+  std::vector<double> exact_average_of_values(n_face_q_points);
+
+  std::vector<Tensor<1, dim>> average_of_grads(n_face_q_points);
+  std::vector<Tensor<1, dim>> exact_average_of_grads(n_face_q_points);
+
+  std::vector<Tensor<2, dim>> average_of_hessians(n_face_q_points);
+  std::vector<Tensor<2, dim>> exact_average_of_hessians(n_face_q_points);
+
+  std::vector<Tensor<3, dim>> average_of_3rd_derivatives(n_face_q_points);
+  std::vector<Tensor<3, dim>> exact_average_of_3rd_derivatives(n_face_q_points);
+
+  const double tol = 1e-11;
+
+  typename DoFHandler<dim>::active_cell_iterator cell =
+    dof_handler.begin_active();
+  for (; cell != dof_handler.end(); ++cell)
+    for (const auto f : cell->face_indices())
+      if (!cell->face(f)->at_boundary())
+        {
+          {
+            fe_iv.reinit(cell,
+                         f,
+                         numbers::invalid_unsigned_int,
+                         cell->neighbor(f),
+                         cell->neighbor_of_neighbor(f),
+                         numbers::invalid_unsigned_int);
+
+            fe_iv.get_jump_in_function_values(solution, jump_in_values);
+            fe_iv[extractor].get_jump_in_function_values(solution,
+                                                         exact_jump_in_values);
+
+            fe_iv.get_jump_in_function_gradients(solution, jump_in_grads);
+            fe_iv[extractor].get_jump_in_function_gradients(
+              solution, exact_jump_in_grads);
+
+            fe_iv.get_jump_in_function_hessians(solution, jump_in_hessians);
+            fe_iv[extractor].get_jump_in_function_hessians(
+              solution, exact_jump_in_hessians);
+
+            fe_iv.get_jump_in_function_third_derivatives(
+              solution, jump_in_3rd_derivatives);
+            fe_iv[extractor].get_jump_in_function_third_derivatives(
+              solution, exact_jump_in_3rd_derivatives);
+
+            fe_iv.get_average_of_function_values(solution, average_of_values);
+            fe_iv[extractor].get_average_of_function_values(
+              solution, exact_average_of_values);
+
+            fe_iv.get_average_of_function_gradients(solution, average_of_grads);
+            fe_iv[extractor].get_average_of_function_gradients(
+              solution, exact_average_of_grads);
+
+            fe_iv.get_average_of_function_hessians(solution,
+                                                   average_of_hessians);
+            fe_iv[extractor].get_average_of_function_hessians(
+              solution, exact_average_of_hessians);
+
+            for (unsigned int q = 0; q < face_quadrature.size(); ++q)
+              {
+                Assert(std::fabs(jump_in_values[q] - exact_jump_in_values[q]) <
+                         tol,
+                       ExcNotImplemented());
+                Assert((jump_in_grads[q] - exact_jump_in_grads[q]).norm() < tol,
+                       ExcNotImplemented());
+                Assert((jump_in_hessians[q] - exact_jump_in_hessians[q])
+                           .norm() < tol,
+                       ExcNotImplemented());
+                Assert((jump_in_3rd_derivatives[q] -
+                        exact_jump_in_3rd_derivatives[q])
+                           .norm() < tol,
+                       ExcNotImplemented());
+
+                Assert(std::fabs(average_of_values[q] -
+                                 exact_average_of_values[q]) < tol,
+                       ExcNotImplemented());
+                Assert((average_of_grads[q] - exact_average_of_grads[q])
+                           .norm() < tol,
+                       ExcNotImplemented());
+                Assert((average_of_hessians[q] - exact_average_of_hessians[q])
+                           .norm() < tol,
+                       ExcNotImplemented());
+              }
+          }
+        }
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(
+    argc, argv, testing_max_num_threads());
+
+  deallog.push("Test jump in function values");
+  {
+    test_jump_function();
+  }
+  deallog.pop();
+  deallog << "OK" << std::endl;
+  deallog.push("Test all functions");
+  {
+    test();
+  }
+  deallog.pop();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/feinterface/fe_interface_values_10.output b/tests/feinterface/fe_interface_values_10.output
new file mode 100644
index 0000000000..62345a50d4
--- /dev/null
+++ b/tests/feinterface/fe_interface_values_10.output
@@ -0,0 +1,5 @@
+
+DEAL:Test jump in function values::OK
+DEAL::OK
+DEAL:Test all functions::OK
+DEAL::OK