]> https://gitweb.dealii.org/ - dealii.git/commitdiff
symmetrize fourth order tensors 17408/head
authorVinayak Vijay <vinayak.dtu@gmail.com>
Wed, 31 Jul 2024 14:04:54 +0000 (19:34 +0530)
committerVinayak Vijay <vinayak.dtu@gmail.com>
Mon, 16 Sep 2024 04:52:53 +0000 (10:22 +0530)
doc/news/changes/minor/20240731Vinayak [new file with mode: 0644]
include/deal.II/base/symmetric_tensor.h
tests/tensors/symmetric_tensor_47.cc [new file with mode: 0644]
tests/tensors/symmetric_tensor_47.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20240731Vinayak b/doc/news/changes/minor/20240731Vinayak
new file mode 100644 (file)
index 0000000..d3be688
--- /dev/null
@@ -0,0 +1,4 @@
+New: Symmetrize fourth order tensors based
+on required symmetry type - major or minor
+<br>
+(Vinayak, 2024/07/31)
index 3fcccf0132c61c051487a425ddb3b8fa805e6a63..15545ccaf5bc42d8da48b6a01d44545089bfe592 100644 (file)
@@ -3546,6 +3546,71 @@ DEAL_II_HOST constexpr inline DEAL_II_ALWAYS_INLINE
 
 
 
+/**
+ * Return the symmetrized version of a full rank-4 tensor, i.e.
+ * as a symmetric rank-4 tensor. The symmetry could be of the following types:
+ * Only minor: $A_{ijkl}=A_{jikl}=A_{ijlk}=A_{jilk}$.
+ * Both minor and major: $A_{ijkl}=A_{jikl}=A_{ijlk}=A_{jilk}$ and
+ * $A_{ijkl}=A_{klij}$. This is the version for general dimensions.
+ * @param t The tensor to be symmetrized.
+ * @param major_symmetry  This argument decides the presence of major symmetry.
+ *
+ * @relatesalso SymmetricTensor
+ */
+template <int dim, typename Number>
+DEAL_II_HOST constexpr inline DEAL_II_ALWAYS_INLINE
+  SymmetricTensor<4, dim, Number>
+  symmetrize(const Tensor<4, dim, Number> &t, const bool major_symmetry)
+{
+  SymmetricTensor<4, dim, Number> result;
+
+  const Number half = internal::NumberType<Number>::value(0.5);
+
+  // minor symmetry - A_{ijkl}=A_{jikl}=A_{ijlk}=A_{jilk}
+  for (unsigned int i = 0; i < dim; ++i)
+    for (unsigned int j = 0; j < dim; ++j)
+      for (unsigned int k = 0; k < dim; ++k)
+        for (unsigned int l = 0; l < dim; ++l)
+          {
+            if (i != j && k == l)
+              {
+                // A_{ijkk}=A_{jikk}
+                result[i][j][k][k] = (t[i][j][k][k] + t[j][i][k][k]) * half;
+              }
+            else if (i == j && k != l)
+              {
+                // A_{iikl}=A_{iilk}
+                result[i][i][k][l] = (t[i][i][k][l] + t[i][i][l][k]) * half;
+              }
+            else if (i != j && k != l)
+              {
+                // A_{ijkl}=A_{jilk}
+                result[i][j][k][l] = (t[i][j][k][l] + t[j][i][k][l] +
+                                      t[i][j][l][k] + t[j][i][l][k]) *
+                                     half * half;
+              }
+            else
+              {
+                // A_{iijj} and A_{iiii} unchanged
+                result[i][j][k][l] = t[i][j][k][l];
+              }
+          }
+
+  // in case major symmetry is also required
+  if (major_symmetry)
+    {
+      // major symmetry - A_{ijkl}=A_{klij}
+      for (unsigned int i = 0; i < dim; ++i)
+        for (unsigned int j = i; j < dim; ++j)
+          for (unsigned int k = 0; k < dim; ++k)
+            for (unsigned int l = k; l < dim; ++l)
+              result[i][j][k][l] = (t[i][j][k][l] + t[k][l][i][j]) * half;
+    }
+  return result;
+}
+
+
+
 /**
  * Multiplication of a symmetric tensor of general rank with a scalar from the
  * right. This version of the operator is used if the scalar has the same data
diff --git a/tests/tensors/symmetric_tensor_47.cc b/tests/tensors/symmetric_tensor_47.cc
new file mode 100644 (file)
index 0000000..3cf5b90
--- /dev/null
@@ -0,0 +1,95 @@
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2023 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  SymmetricTensor::symmetrize on fourth order Tensors
+
+
+#include <deal.II/base/symmetric_tensor.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/vectorization.h>
+
+#include "../tests.h"
+
+
+template <int dim, typename Number>
+void
+make_fourth_order_tensor(Tensor<4, dim, Number> &t)
+{
+  for (unsigned int i = 0; i < dim; ++i)
+    for (unsigned int j = 0; j < dim; ++j)
+      for (unsigned int k = 0; k < dim; ++k)
+        for (unsigned int l = 0; l < dim; ++l)
+          t[i][j][k][l] = (i + 1) * (j + 2) * (k + 3) * (l + 4);
+}
+
+template <int dim, typename Number>
+void
+print_fourth_order_tensor(const SymmetricTensor<4, dim, Number> &t)
+{
+  const Number tolerance = 100. * std::numeric_limits<Number>::epsilon();
+  for (unsigned int i = 0; i < dim; ++i)
+    for (unsigned int j = 0; j < dim; ++j)
+      for (unsigned int k = 0; k < dim; ++k)
+        for (unsigned int l = 0; l < dim; ++l)
+          deallog << i << ' ' << j << ' ' << k << ' ' << l << ' '
+                  << filter_out_small_numbers(t[i][j][k][l], tolerance)
+                  << std::endl;
+}
+
+
+int
+main()
+{
+  initlog();
+  const unsigned int     dim = 3;
+  Tensor<4, dim, double> t;
+  make_fourth_order_tensor(t);
+
+  SymmetricTensor<4, dim, double> st_minor;
+
+  st_minor = symmetrize(t, false);
+
+  // test minor symmetry
+  for (unsigned int i = 0; i < dim; ++i)
+    for (unsigned int j = 0; j < dim; ++j)
+      for (unsigned int k = 0; k < dim; ++k)
+        for (unsigned int l = 0; l < dim; ++l)
+          AssertThrow(st_minor[i][j][k][l] == st_minor[i][j][l][k] &&
+                        st_minor[j][i][k][l] == st_minor[i][j][k][l] &&
+                        st_minor[i][j][k][l] == st_minor[j][i][l][k],
+                      ExcInternalError());
+
+  print_fourth_order_tensor(st_minor);
+
+  deallog << "MINOR SYMMETRY OK" << std::endl;
+
+
+  SymmetricTensor<4, dim, double> st_major;
+
+  st_major = symmetrize(t, true);
+
+  // test major symmetry
+  for (unsigned int i = 0; i < dim; ++i)
+    for (unsigned int j = 0; j < dim; ++j)
+      for (unsigned int k = 0; k < dim; ++k)
+        for (unsigned int l = 0; l < dim; ++l)
+          AssertThrow(st_major[i][j][k][l] == st_major[k][l][i][j],
+                      ExcInternalError());
+
+  print_fourth_order_tensor(st_major);
+
+  deallog << "MAJOR SYMMETRY OK" << std::endl;
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/tensors/symmetric_tensor_47.output b/tests/tensors/symmetric_tensor_47.output
new file mode 100644 (file)
index 0000000..360f4ea
--- /dev/null
@@ -0,0 +1,166 @@
+
+DEAL::0 0 0 0 24.0000
+DEAL::0 0 0 1 31.0000
+DEAL::0 0 0 2 38.0000
+DEAL::0 0 1 0 31.0000
+DEAL::0 0 1 1 40.0000
+DEAL::0 0 1 2 49.0000
+DEAL::0 0 2 0 38.0000
+DEAL::0 0 2 1 49.0000
+DEAL::0 0 2 2 60.0000
+DEAL::0 1 0 0 42.0000
+DEAL::0 1 0 1 54.2500
+DEAL::0 1 0 2 66.5000
+DEAL::0 1 1 0 54.2500
+DEAL::0 1 1 1 70.0000
+DEAL::0 1 1 2 85.7500
+DEAL::0 1 2 0 66.5000
+DEAL::0 1 2 1 85.7500
+DEAL::0 1 2 2 105.000
+DEAL::0 2 0 0 60.0000
+DEAL::0 2 0 1 77.5000
+DEAL::0 2 0 2 95.0000
+DEAL::0 2 1 0 77.5000
+DEAL::0 2 1 1 100.000
+DEAL::0 2 1 2 122.500
+DEAL::0 2 2 0 95.0000
+DEAL::0 2 2 1 122.500
+DEAL::0 2 2 2 150.000
+DEAL::1 0 0 0 42.0000
+DEAL::1 0 0 1 54.2500
+DEAL::1 0 0 2 66.5000
+DEAL::1 0 1 0 54.2500
+DEAL::1 0 1 1 70.0000
+DEAL::1 0 1 2 85.7500
+DEAL::1 0 2 0 66.5000
+DEAL::1 0 2 1 85.7500
+DEAL::1 0 2 2 105.000
+DEAL::1 1 0 0 72.0000
+DEAL::1 1 0 1 93.0000
+DEAL::1 1 0 2 114.000
+DEAL::1 1 1 0 93.0000
+DEAL::1 1 1 1 120.000
+DEAL::1 1 1 2 147.000
+DEAL::1 1 2 0 114.000
+DEAL::1 1 2 1 147.000
+DEAL::1 1 2 2 180.000
+DEAL::1 2 0 0 102.000
+DEAL::1 2 0 1 131.750
+DEAL::1 2 0 2 161.500
+DEAL::1 2 1 0 131.750
+DEAL::1 2 1 1 170.000
+DEAL::1 2 1 2 208.250
+DEAL::1 2 2 0 161.500
+DEAL::1 2 2 1 208.250
+DEAL::1 2 2 2 255.000
+DEAL::2 0 0 0 60.0000
+DEAL::2 0 0 1 77.5000
+DEAL::2 0 0 2 95.0000
+DEAL::2 0 1 0 77.5000
+DEAL::2 0 1 1 100.000
+DEAL::2 0 1 2 122.500
+DEAL::2 0 2 0 95.0000
+DEAL::2 0 2 1 122.500
+DEAL::2 0 2 2 150.000
+DEAL::2 1 0 0 102.000
+DEAL::2 1 0 1 131.750
+DEAL::2 1 0 2 161.500
+DEAL::2 1 1 0 131.750
+DEAL::2 1 1 1 170.000
+DEAL::2 1 1 2 208.250
+DEAL::2 1 2 0 161.500
+DEAL::2 1 2 1 208.250
+DEAL::2 1 2 2 255.000
+DEAL::2 2 0 0 144.000
+DEAL::2 2 0 1 186.000
+DEAL::2 2 0 2 228.000
+DEAL::2 2 1 0 186.000
+DEAL::2 2 1 1 240.000
+DEAL::2 2 1 2 294.000
+DEAL::2 2 2 0 228.000
+DEAL::2 2 2 1 294.000
+DEAL::2 2 2 2 360.000
+DEAL::MINOR SYMMETRY OK
+DEAL::0 0 0 0 24.0000
+DEAL::0 0 0 1 33.0000
+DEAL::0 0 0 2 42.0000
+DEAL::0 0 1 0 33.0000
+DEAL::0 0 1 1 56.0000
+DEAL::0 0 1 2 72.0000
+DEAL::0 0 2 0 42.0000
+DEAL::0 0 2 1 72.0000
+DEAL::0 0 2 2 102.000
+DEAL::0 1 0 0 33.0000
+DEAL::0 1 0 1 45.0000
+DEAL::0 1 0 2 57.0000
+DEAL::0 1 1 0 45.0000
+DEAL::0 1 1 1 75.0000
+DEAL::0 1 1 2 96.0000
+DEAL::0 1 2 0 57.0000
+DEAL::0 1 2 1 96.0000
+DEAL::0 1 2 2 135.000
+DEAL::0 2 0 0 42.0000
+DEAL::0 2 0 1 57.0000
+DEAL::0 2 0 2 72.0000
+DEAL::0 2 1 0 57.0000
+DEAL::0 2 1 1 94.0000
+DEAL::0 2 1 2 120.000
+DEAL::0 2 2 0 72.0000
+DEAL::0 2 2 1 120.000
+DEAL::0 2 2 2 168.000
+DEAL::1 0 0 0 33.0000
+DEAL::1 0 0 1 45.0000
+DEAL::1 0 0 2 57.0000
+DEAL::1 0 1 0 45.0000
+DEAL::1 0 1 1 75.0000
+DEAL::1 0 1 2 96.0000
+DEAL::1 0 2 0 57.0000
+DEAL::1 0 2 1 96.0000
+DEAL::1 0 2 2 135.000
+DEAL::1 1 0 0 56.0000
+DEAL::1 1 0 1 75.0000
+DEAL::1 1 0 2 94.0000
+DEAL::1 1 1 0 75.0000
+DEAL::1 1 1 1 120.000
+DEAL::1 1 1 2 152.000
+DEAL::1 1 2 0 94.0000
+DEAL::1 1 2 1 152.000
+DEAL::1 1 2 2 210.000
+DEAL::1 2 0 0 72.0000
+DEAL::1 2 0 1 96.0000
+DEAL::1 2 0 2 120.000
+DEAL::1 2 1 0 96.0000
+DEAL::1 2 1 1 152.000
+DEAL::1 2 1 2 192.000
+DEAL::1 2 2 0 120.000
+DEAL::1 2 2 1 192.000
+DEAL::1 2 2 2 264.000
+DEAL::2 0 0 0 42.0000
+DEAL::2 0 0 1 57.0000
+DEAL::2 0 0 2 72.0000
+DEAL::2 0 1 0 57.0000
+DEAL::2 0 1 1 94.0000
+DEAL::2 0 1 2 120.000
+DEAL::2 0 2 0 72.0000
+DEAL::2 0 2 1 120.000
+DEAL::2 0 2 2 168.000
+DEAL::2 1 0 0 72.0000
+DEAL::2 1 0 1 96.0000
+DEAL::2 1 0 2 120.000
+DEAL::2 1 1 0 96.0000
+DEAL::2 1 1 1 152.000
+DEAL::2 1 1 2 192.000
+DEAL::2 1 2 0 120.000
+DEAL::2 1 2 1 192.000
+DEAL::2 1 2 2 264.000
+DEAL::2 2 0 0 102.000
+DEAL::2 2 0 1 135.000
+DEAL::2 2 0 2 168.000
+DEAL::2 2 1 0 135.000
+DEAL::2 2 1 1 210.000
+DEAL::2 2 1 2 264.000
+DEAL::2 2 2 0 168.000
+DEAL::2 2 2 1 264.000
+DEAL::2 2 2 2 360.000
+DEAL::MAJOR SYMMETRY OK
+DEAL::OK

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.