]> https://gitweb.dealii.org/ - dealii.git/commitdiff
tests
authorMichał Wichrowski <mtwichrowski@gmail.com>
Thu, 1 May 2025 17:40:42 +0000 (19:40 +0200)
committerMichał Wichrowski <mtwichrowski@gmail.com>
Thu, 1 May 2025 17:40:42 +0000 (19:40 +0200)
include/deal.II/numerics/tensor_product_matrix_creator.h
tests/numerics/tensor_product_creator_01.cc [new file with mode: 0644]
tests/numerics/tensor_product_creator_01.output [new file with mode: 0644]
tests/numerics/tensor_product_creator_02.cc [new file with mode: 0644]
tests/numerics/tensor_product_creator_02.output [new file with mode: 0644]

index 5ea67a6ab672af8fc4eaf4225ec3281ce6171db5..c1caa78c6a043fc89bc3d907095df739e674be01 100644 (file)
@@ -23,6 +23,8 @@
 #include <deal.II/dofs/dof_handler.h>
 
 #include <deal.II/fe/fe.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_tools.h>
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/fe/mapping_q1.h>
@@ -121,7 +123,7 @@ namespace TensorProductMatrixCreator
 
 
   /**
-   * Compute a 1D cell derivative matrix for a given finite element.
+   * Compute a 1D cell laplace matrix for a given finite element.
    *
    * @param fe The finite element object defining the shape functions.
    *
@@ -139,7 +141,7 @@ namespace TensorProductMatrixCreator
    * numbering of the finite element is used.
    */
   FullMatrix<double>
-  create_1d_cell_derivative_matrix(
+  create_1d_cell_laplace_matrix(
     const FiniteElement<1>     &fe,
     const double               &h,
     const std::pair<bool, bool> include_endpoints = {true, true},
@@ -525,11 +527,10 @@ namespace TensorProductMatrixCreator
   }
 
   FullMatrix<double>
-  create_1d_cell_derivative_matrix(
-    const FiniteElement<1>     &fe,
-    const double               &h,
-    const std::pair<bool, bool> include_endpoints,
-    std::vector<unsigned int>   numbering)
+  create_1d_cell_laplace_matrix(const FiniteElement<1>     &fe,
+                                const double               &h,
+                                const std::pair<bool, bool> include_endpoints,
+                                std::vector<unsigned int>   numbering)
   {
     if (dynamic_cast<const FE_DGQ<1> *>(&fe) == nullptr &&
         numbering.size() == 0)
diff --git a/tests/numerics/tensor_product_creator_01.cc b/tests/numerics/tensor_product_creator_01.cc
new file mode 100644 (file)
index 0000000..a7eccb2
--- /dev/null
@@ -0,0 +1,96 @@
+
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2021 - 2024 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+// chceck creation of 1D mass and laplace matrices.
+
+#include <deal.II/matrix_free/fe_evaluation.h>
+
+#include <deal.II/numerics/tensor_product_matrix_creator.h>
+
+#include "../tests.h"
+
+
+void
+print_matrix(FullMatrix<double> &matrix)
+{
+  for (unsigned int i = 0; i < matrix.m(); ++i)
+    {
+      for (unsigned int j = 0; j < matrix.n(); ++j)
+        deallog << std::setw(10) << std::fixed << std::setprecision(4)
+                << matrix(i, j) << " ";
+
+      deallog << std::endl;
+    }
+  deallog << std::endl;
+}
+
+template <int fe_degree>
+void
+test()
+{
+  FE_Q<1> fe_q(fe_degree);
+  FE_Q<1> fe_dgq(fe_degree);
+  auto    mass_matrix =
+    TensorProductMatrixCreator::create_1d_cell_mass_matrix(fe_q,
+                                                           1,
+                                                           {true, true});
+  mass_matrix.print_formatted(std::cout, 4, true, 10, "0.");
+
+  std::cout << std::endl;
+
+
+  FEEvaluation<1, fe_degree, fe_degree + 1, 1, double> fe_eval(
+    fe_q, QGauss<1>(fe_degree + 1), update_values);
+
+
+  auto mass_matrix_renumbered =
+    TensorProductMatrixCreator::create_1d_cell_mass_matrix(
+      fe_q, 1., {true, true}, fe_eval.get_internal_dof_numbering());
+  mass_matrix_renumbered.print_formatted(std::cout, 3, true, 10, "0.");
+  std::cout << std::endl;
+
+  auto mass_matrix_dgq =
+    TensorProductMatrixCreator::create_1d_cell_mass_matrix(fe_dgq,
+                                                           1.,
+                                                           {true, true});
+
+  for (unsigned int i = 0; i < mass_matrix_dgq.m(); ++i)
+    for (unsigned int j = 0; j < mass_matrix_dgq.n(); ++j)
+      if (std::abs(mass_matrix_dgq(i, j) - mass_matrix_renumbered(i, j)) >
+          1e-12)
+        std::cout << "Different at " << i << " " << j << std::endl;
+
+  std::cout << std::endl;
+
+
+  auto mass_matrix_patch =
+    TensorProductMatrixCreator::create_1D_discretization_matrix(
+      mass_matrix_renumbered, 4, 1, {true, true});
+
+  mass_matrix_patch.print_formatted(std::cout, 4, true, 10, "0.");
+  std::cout << std::endl;
+}
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(2);
+  deallog << std::fixed;
+
+  test<1>();
+  test<3>();
+}
\ No newline at end of file
diff --git a/tests/numerics/tensor_product_creator_01.output b/tests/numerics/tensor_product_creator_01.output
new file mode 100644 (file)
index 0000000..f7c81d7
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL::1.0000     -2.0000    1.0000     
+DEAL::-2.0000    4.0000     -2.0000    
+DEAL::1.0000     -2.0000    1.0000     
+DEAL::
+DEAL::51.0000    -130.3444  160.3444   -62.0000   -32.3607   12.3607    1.0000     
+DEAL::-130.3444  338.1966   -425.0000  148.8854   125.0000   -69.0983   12.3607    
+DEAL::160.3444   -425.0000  561.8034   -208.8854  -180.9017  125.0000   -32.3607   
+DEAL::-62.0000   148.8854   -208.8854  244.0000   -208.8854  148.8854   -62.0000   
+DEAL::-32.3607   125.0000   -180.9017  -208.8854  561.8034   -425.0000  160.3444   
+DEAL::12.3607    -69.0983   125.0000   148.8854   -425.0000  338.1966   -130.3444  
+DEAL::1.0000     12.3607    -32.3607   -62.0000   160.3444   -130.3444  51.0000    
+DEAL::
diff --git a/tests/numerics/tensor_product_creator_02.cc b/tests/numerics/tensor_product_creator_02.cc
new file mode 100644 (file)
index 0000000..a19a6ef
--- /dev/null
@@ -0,0 +1,52 @@
+
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2021 - 2024 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+// check creation of 1d ghost penalty matrix
+
+
+#include <deal.II/numerics/tensor_product_matrix_creator.h>
+
+#include "../tests.h"
+
+void
+test(const unsigned int fe_degree)
+{
+  FullMatrix<double> ghost_penalty_1d =
+    TensorProductMatrixCreator::create_1d_ghost_penalty_matrix(
+      FE_Q<1>(fe_degree), 1.);
+
+
+  for (unsigned int row = 0; row < ghost_penalty_1d.n_rows(); ++row)
+    {
+      for (unsigned int col = 0; col < ghost_penalty_1d.n_cols(); ++col)
+        deallog << std::setw(10) << std::fixed << std::setprecision(4)
+                << ghost_penalty_1d(row, col) << " ";
+
+      deallog << std::endl;
+    }
+  deallog << std::endl;
+}
+
+
+int
+main()
+{
+  initlog();
+  deallog << std::setprecision(2);
+  deallog << std::fixed;
+
+  test(1);
+  test(3);
+}
\ No newline at end of file
diff --git a/tests/numerics/tensor_product_creator_02.output b/tests/numerics/tensor_product_creator_02.output
new file mode 100644 (file)
index 0000000..f7c81d7
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL::1.0000     -2.0000    1.0000     
+DEAL::-2.0000    4.0000     -2.0000    
+DEAL::1.0000     -2.0000    1.0000     
+DEAL::
+DEAL::51.0000    -130.3444  160.3444   -62.0000   -32.3607   12.3607    1.0000     
+DEAL::-130.3444  338.1966   -425.0000  148.8854   125.0000   -69.0983   12.3607    
+DEAL::160.3444   -425.0000  561.8034   -208.8854  -180.9017  125.0000   -32.3607   
+DEAL::-62.0000   148.8854   -208.8854  244.0000   -208.8854  148.8854   -62.0000   
+DEAL::-32.3607   125.0000   -180.9017  -208.8854  561.8034   -425.0000  160.3444   
+DEAL::12.3607    -69.0983   125.0000   148.8854   -425.0000  338.1966   -130.3444  
+DEAL::1.0000     12.3607    -32.3607   -62.0000   160.3444   -130.3444  51.0000    
+DEAL::

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.