]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add some SD-related convenience functions for tensor operations.
authorJean-Paul Pelteret <jppelteret@gmail.com>
Wed, 17 Apr 2019 20:29:31 +0000 (22:29 +0200)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Mon, 22 Apr 2019 08:14:28 +0000 (10:14 +0200)
Specifically, these can be used to create a tensor symbol and symbolic
function.

include/deal.II/differentiation/sd.h
include/deal.II/differentiation/sd/symengine_number_types.h
include/deal.II/differentiation/sd/symengine_tensor_operations.h [new file with mode: 0644]
source/differentiation/sd/CMakeLists.txt
source/differentiation/sd/symengine_tensor_operations.cc [new file with mode: 0644]
source/differentiation/sd/symengine_tensor_operations.inst.in [new file with mode: 0644]
tests/symengine/symengine_tensor_operations_01.cc [new file with mode: 0644]
tests/symengine/symengine_tensor_operations_01.output [new file with mode: 0644]

index 5039e8a4a2b56db14a5186d4c68e226a71a88810..44dca06e9e2fefd2e77961176b43f45d5cd1c065 100644 (file)
@@ -25,6 +25,7 @@
 #  include <deal.II/differentiation/sd/symengine_number_types.h>
 #  include <deal.II/differentiation/sd/symengine_product_types.h>
 #  include <deal.II/differentiation/sd/symengine_scalar_operations.h>
+#  include <deal.II/differentiation/sd/symengine_tensor_operations.h>
 #  include <deal.II/differentiation/sd/symengine_types.h>
 #  include <deal.II/differentiation/sd/symengine_utilities.h>
 
index 23db60cd6ec96a00efba107e5f949cf594b411bb..2b4c65330874316b8f92638db440e80f439e68a2 100644 (file)
@@ -46,6 +46,7 @@
 
 // Differentiation
 #  include <deal.II/base/exceptions.h>
+#  include <deal.II/base/numbers.h>
 
 #  include <deal.II/differentiation/sd/symengine_number_traits.h>
 #  include <deal.II/differentiation/sd/symengine_types.h>
@@ -1296,6 +1297,66 @@ namespace Differentiation
 // Specializations from numbers.h
 // These are required in order to make the SD::Expression class a compatible
 // number for the Tensor and SymmetricTensor classes
+
+// Forward declarations:
+template <int rank_, int dim, typename Number>
+class Tensor;
+
+namespace internal
+{
+  template <>
+  struct NumberType<Differentiation::SD::Expression>
+  {
+    static const Differentiation::SD::Expression &
+    value(const Differentiation::SD::Expression &t)
+    {
+      return t;
+    }
+
+    template <
+      typename T,
+      typename = typename std::enable_if<std::is_arithmetic<T>::value>::type>
+    static Differentiation::SD::Expression
+    value(const T &t)
+    {
+      return Differentiation::SD::Expression(t);
+    }
+
+    template <
+      typename T,
+      typename = typename std::enable_if<std::is_arithmetic<T>::value>::type>
+    static Differentiation::SD::Expression
+    value(T &&t)
+    {
+      return Differentiation::SD::Expression(t);
+    }
+
+    template <
+      typename T,
+      typename = typename std::enable_if<std::is_arithmetic<T>::value>::type>
+    static Differentiation::SD::Expression
+    value(const std::complex<T> &t)
+    {
+      return Differentiation::SD::Expression(t);
+    }
+
+    template <typename T, int dim>
+    static Differentiation::SD::Expression
+    value(const Tensor<0, dim, T> &t)
+    {
+      return Differentiation::SD::Expression(static_cast<T>(t));
+    }
+
+    template <typename T, int dim>
+    static Differentiation::SD::Expression
+    value(const Tensor<0, dim, std::complex<T>> &t)
+    {
+      return Differentiation::SD::Expression(static_cast<std::complex<T>>(t));
+    }
+  };
+} // namespace internal
+
+
 namespace numbers
 {
   template <>
@@ -1329,11 +1390,11 @@ namespace numbers
   }
 } // namespace numbers
 
-DEAL_II_NAMESPACE_CLOSE
-
 
 #  endif // DOXYGEN
 
+DEAL_II_NAMESPACE_CLOSE
+
 #endif // DEAL_II_WITH_SYMENGINE
 
 #endif // dealii_differentiation_sd_symengine_number_types_h
diff --git a/include/deal.II/differentiation/sd/symengine_tensor_operations.h b/include/deal.II/differentiation/sd/symengine_tensor_operations.h
new file mode 100644 (file)
index 0000000..d8ba720
--- /dev/null
@@ -0,0 +1,196 @@
+// ---------------------------------------------------------------------
+//
+// 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_differentiation_sd_symengine_tensor_operations_h
+#define dealii_differentiation_sd_symengine_tensor_operations_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_SYMENGINE
+
+#  include <deal.II/base/symmetric_tensor.h>
+#  include <deal.II/base/tensor.h>
+
+#  include <deal.II/differentiation/sd/symengine_number_types.h>
+#  include <deal.II/differentiation/sd/symengine_types.h>
+
+#  include <utility>
+#  include <vector>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Differentiation
+{
+  namespace SD
+  {
+    /**
+     * @name Symbolic variable creation
+     */
+    //@{
+
+    /**
+     * Return a vector of Expressions representing a vectorial symbolic
+     * variable with the identifier specified by @p symbol.
+     *
+     * For example, if the @p symbol is the string `"v"` then the vectorial
+     * symbolic variable that is returned represents the vector $v$. Each
+     * component of $v$ is prefixed by the given @p symbol, and has a suffix that
+     * indicates its component index.
+     *
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the vector of returned
+     *            symbolic variables.
+     * @return A vector (a rank-1 tensor) of symbolic variables with the name of
+     *         each individual component prefixed by @p symbol.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int dim>
+    Tensor<1, dim, Expression>
+    make_vector_of_symbols(const std::string &symbol);
+
+    /**
+     * Return a tensor of Expressions representing a tensorial symbolic
+     * variable with the identifier specified by @p symbol.
+     *
+     * For example, if the @p symbol is the string `"T"` then the tensorial
+     * symbolic variable that is returned represents the vector $T$. Each
+     * component of $T$ is prefixed by the given @p symbol, and has a suffix that
+     * indicates its component indices.
+     *
+     * @tparam rank The rank of the returned tensor.
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the tensor of returned
+     *            symbolic variables.
+     * @return A tensor of symbolic variables with the name of each individual
+     *         component prefixed by @p symbol.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int rank, int dim>
+    Tensor<rank, dim, Expression>
+    make_tensor_of_symbols(const std::string &symbol);
+
+    /**
+     * Return a symmetric tensor of Expressions representing a tensorial
+     * symbolic variable with the identifier specified by @p symbol.
+     *
+     * For example, if the @p symbol is the string `"S"` then the tensorial
+     * symbolic variable that is returned represents the vector $S$. Each
+     * component of $S$ is prefixed by the given @p symbol, and has a suffix that
+     * indicates its component indices.
+     *
+     * @tparam rank The rank of the returned tensor.
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the tensor of returned
+     *            symbolic variables.
+     * @return A tensor of symbolic variables with the name of each individual
+     *         component prefixed by @p symbol.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, Expression>
+    make_symmetric_tensor_of_symbols(const std::string &symbol);
+
+    /**
+     * Return a vector of Expression representing a vectorial symbolic
+     * function with the identifier specified by @p symbol. The functions'
+     * symbolic dependencies are specified by the keys to the input
+     * @p arguments map; the values stored in the map are ignored.
+     *
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the vector of returned
+     *            symbolic functions.
+     * @param[in] arguments A map of input arguments to the returned
+     *            symbolic functions.
+     * @return A vector (a rank-1 tensor) of generic symbolic functions with the
+     *         name of each individual component prefixed by @p symbol, a suffix
+     *         that indicates its component index, and the number of input
+     *         arguments equal to the length of @p arguments.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int dim>
+    Tensor<1, dim, Expression>
+    make_vector_of_symbolic_functions(const std::string &            symbol,
+                                      const types::substitution_map &arguments);
+
+    /**
+     * Return a tensor of Expression representing a tensorial symbolic
+     * function with the identifier specified by @p symbol. The functions'
+     * symbolic dependencies are specified by the keys to the input
+     * @p arguments map; the values stored in the map are ignored.
+     *
+     * @tparam rank The rank of the returned tensor.
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the tensor of returned
+     *            symbolic functions.
+     * @param[in] arguments A map of input arguments to the returned
+     *            symbolic functions.
+     * @return A tensor of generic symbolic functions with the name of each
+     *         individual component prefixed by @p symbol, a suffix that
+     *         indicates its component indeices, and the number of input
+     *         arguments equal to  the length of @p arguments.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int rank, int dim>
+    Tensor<rank, dim, Expression>
+    make_tensor_of_symbolic_functions(const std::string &            symbol,
+                                      const types::substitution_map &arguments);
+
+    /**
+     * Return a symmetric tensor of Expression representing a tensorial
+     * symbolic function with the identifier specified by @p symbol. The
+     * functions' symbolic dependencies are specified by the keys to the input
+     * @p arguments map; the values stored in the map are ignored.
+     *
+     * @tparam rank The rank of the returned tensor.
+     * @tparam dim The dimension of the returned tensor.
+     * @param[in] symbol An indentifier (or name) for the tensor of returned
+     *            symbolic functions.
+     * @param[in] arguments A map of input arguments to the returned
+     *            symbolic functions.
+     * @return A symmetric tensor of generic symbolic functions with the name of
+     *         each individual component prefixed by @p symbol, a suffix that
+     *         indicates its component indeices, and the number of input
+     *         arguments equal to the length of @p arguments.
+     *
+     * @warning It is up to the user to ensure that there is no ambiguity in the
+     * symbols used within a section of code.
+     */
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, Expression>
+    make_symmetric_tensor_of_symbolic_functions(
+      const std::string &            symbol,
+      const types::substitution_map &arguments);
+
+    //@}
+
+  } // namespace SD
+} // namespace Differentiation
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_SYMENGINE
+
+#endif
index 0cb7b56a5e7e23aa5dc279c7721e4fb1443cee77..a0de87e090d64d0aab42ca1c5396a2fe4e3562ba 100644 (file)
@@ -19,11 +19,13 @@ SET(_src
   symengine_math.cc
   symengine_number_types.cc
   symengine_scalar_operations.cc
+  symengine_tensor_operations.cc
   symengine_types.cc
   symengine_utilities.cc
   )
 
 SET(_inst
+  symengine_tensor_operations.inst.in
   )
 
 FILE(GLOB _header
diff --git a/source/differentiation/sd/symengine_tensor_operations.cc b/source/differentiation/sd/symengine_tensor_operations.cc
new file mode 100644 (file)
index 0000000..bd2e743
--- /dev/null
@@ -0,0 +1,313 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 - 2017 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_SYMENGINE
+
+#  include <deal.II/base/symmetric_tensor.h>
+#  include <deal.II/base/tensor.h>
+#  include <deal.II/base/utilities.h>
+
+#  include <deal.II/differentiation/sd/symengine_product_types.h>
+#  include <deal.II/differentiation/sd/symengine_scalar_operations.h>
+#  include <deal.II/differentiation/sd/symengine_tensor_operations.h>
+#  include <deal.II/differentiation/sd/symengine_types.h>
+#  include <deal.II/differentiation/sd/symengine_utilities.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace Differentiation
+{
+  namespace SD
+  {
+    // --- Symbolic variable creation ---
+
+    namespace internal
+    {
+      namespace
+      {
+        template <int dim>
+        TableIndices<4>
+        make_rank_4_tensor_indices(const unsigned int &idx_i,
+                                   const unsigned int &idx_j)
+        {
+          const TableIndices<2> indices_i(
+            SymmetricTensor<2, dim>::unrolled_to_component_indices(idx_i));
+          const TableIndices<2> indices_j(
+            SymmetricTensor<2, dim>::unrolled_to_component_indices(idx_j));
+          return TableIndices<4>(indices_i[0],
+                                 indices_i[1],
+                                 indices_j[0],
+                                 indices_j[1]);
+        }
+
+
+        template <int rank>
+        std::string
+        make_index_string(const TableIndices<rank> &indices)
+        {
+          std::string out;
+          for (unsigned int i = 0; i < rank; ++i)
+            out += dealii::Utilities::to_string(indices[i]);
+          return out;
+        }
+
+
+        template <int rank, int dim>
+        struct Symbol_Tensor;
+
+
+        template <int rank, int dim>
+        struct Symbol_Tensor
+        {
+          static Tensor<rank, dim, Expression>
+          create(const std::string &sym)
+          {
+            Tensor<rank, dim, Expression> out;
+            for (unsigned int i = 0; i < out.n_independent_components; ++i)
+              {
+                const TableIndices<rank> indices(
+                  out.unrolled_to_component_indices(i));
+                out[indices] = Expression(sym + std::string("_") +
+                                          internal::make_index_string(indices));
+              }
+            return out;
+          }
+        };
+
+
+        template <int dim>
+        struct Symbol_Tensor<0, dim>
+        {
+          static Tensor<0, dim, Expression>
+          create(const std::string &sym)
+          {
+            return make_symbol(sym);
+          }
+        };
+
+
+        template <int rank, int dim>
+        struct Symbol_SymmetricTensor;
+
+
+        template <int rank, int dim>
+        struct Symbol_SymmetricTensor
+        {
+          static SymmetricTensor<2, dim, Expression>
+          create(const std::string &sym)
+          {
+            SymmetricTensor<rank, dim, Expression> out;
+            for (unsigned int i = 0; i < out.n_independent_components; ++i)
+              {
+                const TableIndices<rank> indices(
+                  out.unrolled_to_component_indices(i));
+                out[indices] = Expression(sym + std::string("_") +
+                                          internal::make_index_string(indices));
+              }
+            return out;
+          }
+        };
+
+
+        template <int dim>
+        struct Symbol_SymmetricTensor<4, dim>
+        {
+          static SymmetricTensor<4, dim, Expression>
+          create(const std::string &sym)
+          {
+            SymmetricTensor<4, dim, Expression> out;
+            for (unsigned int i = 0;
+                 i < SymmetricTensor<2, dim>::n_independent_components;
+                 ++i)
+              for (unsigned int j = 0;
+                   j < SymmetricTensor<2, dim>::n_independent_components;
+                   ++j)
+                {
+                  const TableIndices<4> indices =
+                    make_rank_4_tensor_indices<dim>(i, j);
+                  out[indices] =
+                    Expression(sym + std::string("_") +
+                               internal::make_index_string(indices));
+                }
+            return out;
+          }
+        };
+
+
+        template <int rank, int dim>
+        struct Symbol_Function_Tensor;
+
+
+        template <int rank, int dim>
+        struct Symbol_Function_Tensor
+        {
+          static Tensor<rank, dim, Expression>
+          create(const std::string &            sym,
+                 const types::substitution_map &arguments)
+          {
+            Tensor<rank, dim, Expression> out;
+            const types::symbol_vector    args =
+              Utilities::extract_symbols(arguments);
+            for (unsigned int i = 0; i < out.n_independent_components; ++i)
+              {
+                const TableIndices<rank> indices(
+                  out.unrolled_to_component_indices(i));
+                out[indices] =
+                  Expression(sym + std::string("_") +
+                               internal::make_index_string(indices),
+                             args);
+              }
+            return out;
+          }
+        };
+
+
+        template <int dim>
+        struct Symbol_Function_Tensor<0, dim>
+        {
+          static Tensor<0, dim, Expression>
+          create(const std::string &            sym,
+                 const types::substitution_map &arguments)
+          {
+            return make_symbolic_function(sym, arguments);
+          }
+        };
+
+
+        template <int rank, int dim>
+        struct Symbol_Function_SymmetricTensor;
+
+
+        template <int rank, int dim>
+        struct Symbol_Function_SymmetricTensor
+        {
+          static SymmetricTensor<2, dim, Expression>
+          create(const std::string &            sym,
+                 const types::substitution_map &arguments)
+          {
+            SymmetricTensor<rank, dim, Expression> out;
+            const types::symbol_vector             args =
+              Utilities::extract_symbols(arguments);
+            for (unsigned int i = 0; i < out.n_independent_components; ++i)
+              {
+                const TableIndices<rank> indices(
+                  out.unrolled_to_component_indices(i));
+                out[indices] =
+                  Expression(sym + std::string("_") +
+                               internal::make_index_string(indices),
+                             args);
+              }
+            return out;
+          }
+        };
+
+
+        template <int dim>
+        struct Symbol_Function_SymmetricTensor<4, dim>
+        {
+          static SymmetricTensor<4, dim, Expression>
+          create(const std::string &            sym,
+                 const types::substitution_map &arguments)
+          {
+            SymmetricTensor<4, dim, Expression> out;
+            const types::symbol_vector          args =
+              Utilities::extract_symbols(arguments);
+            for (unsigned int i = 0;
+                 i < SymmetricTensor<2, dim>::n_independent_components;
+                 ++i)
+              for (unsigned int j = 0;
+                   j < SymmetricTensor<2, dim>::n_independent_components;
+                   ++j)
+                {
+                  const TableIndices<4> indices =
+                    make_rank_4_tensor_indices<dim>(i, j);
+                  out[indices] =
+                    Expression(sym + std::string("_") +
+                                 internal::make_index_string(indices),
+                               args);
+                }
+            return out;
+          }
+        };
+      } // namespace
+    }   // namespace internal
+
+
+    template <int dim>
+    Tensor<1, dim, Expression>
+    make_vector_of_symbols(const std::string &sym)
+    {
+      return internal::Symbol_Tensor<1, dim>::create(sym);
+    }
+
+
+    template <int dim>
+    Tensor<1, dim, Expression>
+    make_vector_of_symbolic_functions(const std::string &            sym,
+                                      const types::substitution_map &arguments)
+    {
+      return internal::Symbol_Function_Tensor<1, dim>::create(sym, arguments);
+    }
+
+
+    template <int rank, int dim>
+    Tensor<rank, dim, Expression>
+    make_tensor_of_symbols(const std::string &sym)
+    {
+      return internal::Symbol_Tensor<rank, dim>::create(sym);
+    }
+
+
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, Expression>
+    make_symmetric_tensor_of_symbols(const std::string &sym)
+    {
+      return internal::Symbol_SymmetricTensor<rank, dim>::create(sym);
+    }
+
+
+    template <int rank, int dim>
+    Tensor<rank, dim, Expression>
+    make_tensor_of_symbolic_functions(const std::string &            sym,
+                                      const types::substitution_map &arguments)
+    {
+      return internal::Symbol_Function_Tensor<rank, dim>::create(sym,
+                                                                 arguments);
+    }
+
+
+    template <int rank, int dim>
+    SymmetricTensor<rank, dim, Expression>
+    make_symmetric_tensor_of_symbolic_functions(
+      const std::string &            sym,
+      const types::substitution_map &arguments)
+    {
+      return internal::Symbol_Function_SymmetricTensor<rank, dim>::create(
+        sym, arguments);
+    }
+
+  } // namespace SD
+} // namespace Differentiation
+
+/* --- Explicit instantiations --- */
+
+#  include "symengine_tensor_operations.inst"
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_WITH_SYMENGINE
diff --git a/source/differentiation/sd/symengine_tensor_operations.inst.in b/source/differentiation/sd/symengine_tensor_operations.inst.in
new file mode 100644 (file)
index 0000000..1a4a802
--- /dev/null
@@ -0,0 +1,85 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+for (deal_II_dimension : DIMENSIONS)
+  {
+    namespace Differentiation
+    \{
+      namespace SD
+      \{
+
+        template Tensor<1, deal_II_dimension, Expression>
+        make_vector_of_symbols<deal_II_dimension>(const std::string &symbol);
+
+        template Tensor<1, deal_II_dimension, Expression>
+        make_vector_of_symbolic_functions<deal_II_dimension>(
+          const std::string &            symbol,
+          const types::substitution_map &arguments);
+
+        /*--- Rank-0 definitions ---*/
+
+        template Tensor<0, deal_II_dimension, Expression>
+        make_tensor_of_symbols<0, deal_II_dimension>(const std::string &symbol);
+
+        template Tensor<0, deal_II_dimension, Expression>
+        make_tensor_of_symbolic_functions<0, deal_II_dimension>(
+          const std::string &            symbol,
+          const types::substitution_map &arguments);
+      \}
+    \}
+  }
+
+
+for (rank : RANKS; deal_II_dimension : DIMENSIONS)
+  {
+    namespace Differentiation
+    \{
+      namespace SD
+      \{
+
+        template Tensor<rank, deal_II_dimension, Expression>
+        make_tensor_of_symbols<rank, deal_II_dimension>(
+          const std::string &symbol);
+
+        template Tensor<rank, deal_II_dimension, Expression>
+        make_tensor_of_symbolic_functions<rank, deal_II_dimension>(
+          const std::string &            symbol,
+          const types::substitution_map &arguments);
+
+      \}
+    \}
+  }
+
+
+for (rank : SYM_RANKS; deal_II_dimension : DIMENSIONS)
+  {
+    namespace Differentiation
+    \{
+      namespace SD
+      \{
+
+        template SymmetricTensor<rank, deal_II_dimension, Expression>
+        make_symmetric_tensor_of_symbols<rank, deal_II_dimension>(
+          const std::string &symbol);
+
+        template SymmetricTensor<rank, deal_II_dimension, Expression>
+        make_symmetric_tensor_of_symbolic_functions<rank, deal_II_dimension>(
+          const std::string &            symbol,
+          const types::substitution_map &arguments);
+
+      \}
+    \}
+  }
diff --git a/tests/symengine/symengine_tensor_operations_01.cc b/tests/symengine/symengine_tensor_operations_01.cc
new file mode 100644 (file)
index 0000000..b750d40
--- /dev/null
@@ -0,0 +1,112 @@
+// ---------------------------------------------------------------------
+//
+// 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Check that symbol tensors and symbolic function tensors can be created
+// using the utility functions
+
+#include <deal.II/differentiation/sd.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+namespace SD = Differentiation::SD;
+
+
+template <int dim>
+void
+test()
+{
+  deallog.push("dim " + Utilities::to_string(dim));
+
+  using SD_number_t = SD::Expression;
+
+  const Tensor<1, dim, SD_number_t> symb_v =
+    SD::make_vector_of_symbols<dim>("v");
+  const Tensor<0, dim, SD_number_t> symb_T0 =
+    SD::make_tensor_of_symbols<0, dim>("T0");
+  const Tensor<1, dim, SD_number_t> symb_T1 =
+    SD::make_tensor_of_symbols<1, dim>("T1");
+  const Tensor<2, dim, SD_number_t> symb_T2 =
+    SD::make_tensor_of_symbols<2, dim>("T2");
+  const Tensor<3, dim, SD_number_t> symb_T3 =
+    SD::make_tensor_of_symbols<3, dim>("T3");
+  const Tensor<4, dim, SD_number_t> symb_T4 =
+    SD::make_tensor_of_symbols<4, dim>("T4");
+  const SymmetricTensor<2, dim, SD_number_t> symb_ST2 =
+    SD::make_symmetric_tensor_of_symbols<2, dim>("ST2");
+  const SymmetricTensor<4, dim, SD_number_t> symb_ST4 =
+    SD::make_symmetric_tensor_of_symbols<4, dim>("ST4");
+
+  deallog << "Symbol vector: " << symb_v << std::endl;
+  deallog << "Symbol tensor (rank-0): " << symb_T0 << std::endl;
+  deallog << "Symbol tensor (rank-1): " << symb_T1 << std::endl;
+  deallog << "Symbol tensor (rank-2): " << symb_T2 << std::endl;
+  deallog << "Symbol tensor (rank-3): " << symb_T3 << std::endl;
+  deallog << "Symbol tensor (rank-4): " << symb_T4 << std::endl;
+  deallog << "Symbol symmetric tensor (rank-2): " << symb_ST2 << std::endl;
+  deallog << "Symbol symmetric tensor (rank-4): " << symb_ST4 << std::endl;
+
+  SD::types::substitution_map sub_map;
+  sub_map[symb_T0] = 0.0;
+  for (unsigned int i = 0; i < dim; ++i)
+    {
+      sub_map[symb_v[i]] = 0.0;
+    }
+
+  const Tensor<1, dim, SD_number_t> symb_func_v =
+    SD::make_vector_of_symbolic_functions<dim>("f_v", sub_map);
+  const Tensor<0, dim, SD_number_t> symb_func_T0 =
+    SD::make_tensor_of_symbolic_functions<0, dim>("f_T0", sub_map);
+  const Tensor<1, dim, SD_number_t> symb_func_T1 =
+    SD::make_tensor_of_symbolic_functions<1, dim>("f_T1", sub_map);
+  const Tensor<2, dim, SD_number_t> symb_func_T2 =
+    SD::make_tensor_of_symbolic_functions<2, dim>("f_T2", sub_map);
+  const Tensor<3, dim, SD_number_t> symb_func_T3 =
+    SD::make_tensor_of_symbolic_functions<3, dim>("f_T3", sub_map);
+  const Tensor<4, dim, SD_number_t> symb_func_T4 =
+    SD::make_tensor_of_symbolic_functions<4, dim>("f_T4", sub_map);
+  const SymmetricTensor<2, dim, SD_number_t> symb_func_ST2 =
+    SD::make_symmetric_tensor_of_symbolic_functions<2, dim>("f_ST2", sub_map);
+  const SymmetricTensor<4, dim, SD_number_t> symb_func_ST4 =
+    SD::make_symmetric_tensor_of_symbolic_functions<4, dim>("f_ST4", sub_map);
+
+  deallog << "Symbol function vector: " << symb_func_v << std::endl;
+  deallog << "Symbol function tensor (rank-0): " << symb_func_T0 << std::endl;
+  deallog << "Symbol function tensor (rank-1): " << symb_func_T1 << std::endl;
+  deallog << "Symbol function tensor (rank-2): " << symb_func_T2 << std::endl;
+  deallog << "Symbol function tensor (rank-3): " << symb_func_T3 << std::endl;
+  deallog << "Symbol function tensor (rank-4): " << symb_func_T4 << std::endl;
+  deallog << "Symbol function symmetric tensor (rank-2): " << symb_func_ST2
+          << std::endl;
+  deallog << "Symbol function symmetric tensor (rank-4): " << symb_func_ST4
+          << std::endl;
+
+  deallog << "OK" << std::endl;
+  deallog.pop();
+}
+
+
+int
+main()
+{
+  initlog();
+
+  test<1>();
+  test<2>();
+  test<3>();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/symengine/symengine_tensor_operations_01.output b/tests/symengine/symengine_tensor_operations_01.output
new file mode 100644 (file)
index 0000000..ee6a977
--- /dev/null
@@ -0,0 +1,53 @@
+
+DEAL:dim 1::Symbol vector: v_0
+DEAL:dim 1::Symbol tensor (rank-0): T0
+DEAL:dim 1::Symbol tensor (rank-1): T1_0
+DEAL:dim 1::Symbol tensor (rank-2): T2_00
+DEAL:dim 1::Symbol tensor (rank-3): T3_000
+DEAL:dim 1::Symbol tensor (rank-4): T4_0000
+DEAL:dim 1::Symbol symmetric tensor (rank-2): ST2_00
+DEAL:dim 1::Symbol symmetric tensor (rank-4): ST4_0000
+DEAL:dim 1::Symbol function vector: f_v_0(T0, v_0)
+DEAL:dim 1::Symbol function tensor (rank-0): f_T0(T0, v_0)
+DEAL:dim 1::Symbol function tensor (rank-1): f_T1_0(T0, v_0)
+DEAL:dim 1::Symbol function tensor (rank-2): f_T2_00(T0, v_0)
+DEAL:dim 1::Symbol function tensor (rank-3): f_T3_000(T0, v_0)
+DEAL:dim 1::Symbol function tensor (rank-4): f_T4_0000(T0, v_0)
+DEAL:dim 1::Symbol function symmetric tensor (rank-2): f_ST2_00(T0, v_0)
+DEAL:dim 1::Symbol function symmetric tensor (rank-4): f_ST4_0000(T0, v_0)
+DEAL:dim 1::OK
+DEAL:dim 2::Symbol vector: v_0 v_1
+DEAL:dim 2::Symbol tensor (rank-0): T0
+DEAL:dim 2::Symbol tensor (rank-1): T1_0 T1_1
+DEAL:dim 2::Symbol tensor (rank-2): T2_00 T2_01 T2_10 T2_11
+DEAL:dim 2::Symbol tensor (rank-3): T3_000 T3_001 T3_010 T3_011 T3_100 T3_101 T3_110 T3_111
+DEAL:dim 2::Symbol tensor (rank-4): T4_0000 T4_0001 T4_0010 T4_0011 T4_0100 T4_0101 T4_0110 T4_0111 T4_1000 T4_1001 T4_1010 T4_1011 T4_1100 T4_1101 T4_1110 T4_1111
+DEAL:dim 2::Symbol symmetric tensor (rank-2): ST2_00 ST2_01 ST2_01 ST2_11
+DEAL:dim 2::Symbol symmetric tensor (rank-4): ST4_0000 ST4_0001 ST4_0001 ST4_0011 ST4_0100 ST4_0101 ST4_0101 ST4_0111 ST4_0100 ST4_0101 ST4_0101 ST4_0111 ST4_1100 ST4_1101 ST4_1101 ST4_1111
+DEAL:dim 2::Symbol function vector: f_v_0(T0, v_0, v_1) f_v_1(T0, v_0, v_1)
+DEAL:dim 2::Symbol function tensor (rank-0): f_T0(T0, v_0, v_1)
+DEAL:dim 2::Symbol function tensor (rank-1): f_T1_0(T0, v_0, v_1) f_T1_1(T0, v_0, v_1)
+DEAL:dim 2::Symbol function tensor (rank-2): f_T2_00(T0, v_0, v_1) f_T2_01(T0, v_0, v_1) f_T2_10(T0, v_0, v_1) f_T2_11(T0, v_0, v_1)
+DEAL:dim 2::Symbol function tensor (rank-3): f_T3_000(T0, v_0, v_1) f_T3_001(T0, v_0, v_1) f_T3_010(T0, v_0, v_1) f_T3_011(T0, v_0, v_1) f_T3_100(T0, v_0, v_1) f_T3_101(T0, v_0, v_1) f_T3_110(T0, v_0, v_1) f_T3_111(T0, v_0, v_1)
+DEAL:dim 2::Symbol function tensor (rank-4): f_T4_0000(T0, v_0, v_1) f_T4_0001(T0, v_0, v_1) f_T4_0010(T0, v_0, v_1) f_T4_0011(T0, v_0, v_1) f_T4_0100(T0, v_0, v_1) f_T4_0101(T0, v_0, v_1) f_T4_0110(T0, v_0, v_1) f_T4_0111(T0, v_0, v_1) f_T4_1000(T0, v_0, v_1) f_T4_1001(T0, v_0, v_1) f_T4_1010(T0, v_0, v_1) f_T4_1011(T0, v_0, v_1) f_T4_1100(T0, v_0, v_1) f_T4_1101(T0, v_0, v_1) f_T4_1110(T0, v_0, v_1) f_T4_1111(T0, v_0, v_1)
+DEAL:dim 2::Symbol function symmetric tensor (rank-2): f_ST2_00(T0, v_0, v_1) f_ST2_01(T0, v_0, v_1) f_ST2_01(T0, v_0, v_1) f_ST2_11(T0, v_0, v_1)
+DEAL:dim 2::Symbol function symmetric tensor (rank-4): f_ST4_0000(T0, v_0, v_1) f_ST4_0001(T0, v_0, v_1) f_ST4_0001(T0, v_0, v_1) f_ST4_0011(T0, v_0, v_1) f_ST4_0100(T0, v_0, v_1) f_ST4_0101(T0, v_0, v_1) f_ST4_0101(T0, v_0, v_1) f_ST4_0111(T0, v_0, v_1) f_ST4_0100(T0, v_0, v_1) f_ST4_0101(T0, v_0, v_1) f_ST4_0101(T0, v_0, v_1) f_ST4_0111(T0, v_0, v_1) f_ST4_1100(T0, v_0, v_1) f_ST4_1101(T0, v_0, v_1) f_ST4_1101(T0, v_0, v_1) f_ST4_1111(T0, v_0, v_1)
+DEAL:dim 2::OK
+DEAL:dim 3::Symbol vector: v_0 v_1 v_2
+DEAL:dim 3::Symbol tensor (rank-0): T0
+DEAL:dim 3::Symbol tensor (rank-1): T1_0 T1_1 T1_2
+DEAL:dim 3::Symbol tensor (rank-2): T2_00 T2_01 T2_02 T2_10 T2_11 T2_12 T2_20 T2_21 T2_22
+DEAL:dim 3::Symbol tensor (rank-3): T3_000 T3_001 T3_002 T3_010 T3_011 T3_012 T3_020 T3_021 T3_022 T3_100 T3_101 T3_102 T3_110 T3_111 T3_112 T3_120 T3_121 T3_122 T3_200 T3_201 T3_202 T3_210 T3_211 T3_212 T3_220 T3_221 T3_222
+DEAL:dim 3::Symbol tensor (rank-4): T4_0000 T4_0001 T4_0002 T4_0010 T4_0011 T4_0012 T4_0020 T4_0021 T4_0022 T4_0100 T4_0101 T4_0102 T4_0110 T4_0111 T4_0112 T4_0120 T4_0121 T4_0122 T4_0200 T4_0201 T4_0202 T4_0210 T4_0211 T4_0212 T4_0220 T4_0221 T4_0222 T4_1000 T4_1001 T4_1002 T4_1010 T4_1011 T4_1012 T4_1020 T4_1021 T4_1022 T4_1100 T4_1101 T4_1102 T4_1110 T4_1111 T4_1112 T4_1120 T4_1121 T4_1122 T4_1200 T4_1201 T4_1202 T4_1210 T4_1211 T4_1212 T4_1220 T4_1221 T4_1222 T4_2000 T4_2001 T4_2002 T4_2010 T4_2011 T4_2012 T4_2020 T4_2021 T4_2022 T4_2100 T4_2101 T4_2102 T4_2110 T4_2111 T4_2112 T4_2120 T4_2121 T4_2122 T4_2200 T4_2201 T4_2202 T4_2210 T4_2211 T4_2212 T4_2220 T4_2221 T4_2222
+DEAL:dim 3::Symbol symmetric tensor (rank-2): ST2_00 ST2_01 ST2_02 ST2_01 ST2_11 ST2_12 ST2_02 ST2_12 ST2_22
+DEAL:dim 3::Symbol symmetric tensor (rank-4): ST4_0000 ST4_0001 ST4_0002 ST4_0001 ST4_0011 ST4_0012 ST4_0002 ST4_0012 ST4_0022 ST4_0100 ST4_0101 ST4_0102 ST4_0101 ST4_0111 ST4_0112 ST4_0102 ST4_0112 ST4_0122 ST4_0200 ST4_0201 ST4_0202 ST4_0201 ST4_0211 ST4_0212 ST4_0202 ST4_0212 ST4_0222 ST4_0100 ST4_0101 ST4_0102 ST4_0101 ST4_0111 ST4_0112 ST4_0102 ST4_0112 ST4_0122 ST4_1100 ST4_1101 ST4_1102 ST4_1101 ST4_1111 ST4_1112 ST4_1102 ST4_1112 ST4_1122 ST4_1200 ST4_1201 ST4_1202 ST4_1201 ST4_1211 ST4_1212 ST4_1202 ST4_1212 ST4_1222 ST4_0200 ST4_0201 ST4_0202 ST4_0201 ST4_0211 ST4_0212 ST4_0202 ST4_0212 ST4_0222 ST4_1200 ST4_1201 ST4_1202 ST4_1201 ST4_1211 ST4_1212 ST4_1202 ST4_1212 ST4_1222 ST4_2200 ST4_2201 ST4_2202 ST4_2201 ST4_2211 ST4_2212 ST4_2202 ST4_2212 ST4_2222
+DEAL:dim 3::Symbol function vector: f_v_0(T0, v_0, v_1, v_2) f_v_1(T0, v_0, v_1, v_2) f_v_2(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function tensor (rank-0): f_T0(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function tensor (rank-1): f_T1_0(T0, v_0, v_1, v_2) f_T1_1(T0, v_0, v_1, v_2) f_T1_2(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function tensor (rank-2): f_T2_00(T0, v_0, v_1, v_2) f_T2_01(T0, v_0, v_1, v_2) f_T2_02(T0, v_0, v_1, v_2) f_T2_10(T0, v_0, v_1, v_2) f_T2_11(T0, v_0, v_1, v_2) f_T2_12(T0, v_0, v_1, v_2) f_T2_20(T0, v_0, v_1, v_2) f_T2_21(T0, v_0, v_1, v_2) f_T2_22(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function tensor (rank-3): f_T3_000(T0, v_0, v_1, v_2) f_T3_001(T0, v_0, v_1, v_2) f_T3_002(T0, v_0, v_1, v_2) f_T3_010(T0, v_0, v_1, v_2) f_T3_011(T0, v_0, v_1, v_2) f_T3_012(T0, v_0, v_1, v_2) f_T3_020(T0, v_0, v_1, v_2) f_T3_021(T0, v_0, v_1, v_2) f_T3_022(T0, v_0, v_1, v_2) f_T3_100(T0, v_0, v_1, v_2) f_T3_101(T0, v_0, v_1, v_2) f_T3_102(T0, v_0, v_1, v_2) f_T3_110(T0, v_0, v_1, v_2) f_T3_111(T0, v_0, v_1, v_2) f_T3_112(T0, v_0, v_1, v_2) f_T3_120(T0, v_0, v_1, v_2) f_T3_121(T0, v_0, v_1, v_2) f_T3_122(T0, v_0, v_1, v_2) f_T3_200(T0, v_0, v_1, v_2) f_T3_201(T0, v_0, v_1, v_2) f_T3_202(T0, v_0, v_1, v_2) f_T3_210(T0, v_0, v_1, v_2) f_T3_211(T0, v_0, v_1, v_2) f_T3_212(T0, v_0, v_1, v_2) f_T3_220(T0, v_0, v_1, v_2) f_T3_221(T0, v_0, v_1, v_2) f_T3_222(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function tensor (rank-4): f_T4_0000(T0, v_0, v_1, v_2) f_T4_0001(T0, v_0, v_1, v_2) f_T4_0002(T0, v_0, v_1, v_2) f_T4_0010(T0, v_0, v_1, v_2) f_T4_0011(T0, v_0, v_1, v_2) f_T4_0012(T0, v_0, v_1, v_2) f_T4_0020(T0, v_0, v_1, v_2) f_T4_0021(T0, v_0, v_1, v_2) f_T4_0022(T0, v_0, v_1, v_2) f_T4_0100(T0, v_0, v_1, v_2) f_T4_0101(T0, v_0, v_1, v_2) f_T4_0102(T0, v_0, v_1, v_2) f_T4_0110(T0, v_0, v_1, v_2) f_T4_0111(T0, v_0, v_1, v_2) f_T4_0112(T0, v_0, v_1, v_2) f_T4_0120(T0, v_0, v_1, v_2) f_T4_0121(T0, v_0, v_1, v_2) f_T4_0122(T0, v_0, v_1, v_2) f_T4_0200(T0, v_0, v_1, v_2) f_T4_0201(T0, v_0, v_1, v_2) f_T4_0202(T0, v_0, v_1, v_2) f_T4_0210(T0, v_0, v_1, v_2) f_T4_0211(T0, v_0, v_1, v_2) f_T4_0212(T0, v_0, v_1, v_2) f_T4_0220(T0, v_0, v_1, v_2) f_T4_0221(T0, v_0, v_1, v_2) f_T4_0222(T0, v_0, v_1, v_2) f_T4_1000(T0, v_0, v_1, v_2) f_T4_1001(T0, v_0, v_1, v_2) f_T4_1002(T0, v_0, v_1, v_2) f_T4_1010(T0, v_0, v_1, v_2) f_T4_1011(T0, v_0, v_1, v_2) f_T4_1012(T0, v_0, v_1, v_2) f_T4_1020(T0, v_0, v_1, v_2) f_T4_1021(T0, v_0, v_1, v_2) f_T4_1022(T0, v_0, v_1, v_2) f_T4_1100(T0, v_0, v_1, v_2) f_T4_1101(T0, v_0, v_1, v_2) f_T4_1102(T0, v_0, v_1, v_2) f_T4_1110(T0, v_0, v_1, v_2) f_T4_1111(T0, v_0, v_1, v_2) f_T4_1112(T0, v_0, v_1, v_2) f_T4_1120(T0, v_0, v_1, v_2) f_T4_1121(T0, v_0, v_1, v_2) f_T4_1122(T0, v_0, v_1, v_2) f_T4_1200(T0, v_0, v_1, v_2) f_T4_1201(T0, v_0, v_1, v_2) f_T4_1202(T0, v_0, v_1, v_2) f_T4_1210(T0, v_0, v_1, v_2) f_T4_1211(T0, v_0, v_1, v_2) f_T4_1212(T0, v_0, v_1, v_2) f_T4_1220(T0, v_0, v_1, v_2) f_T4_1221(T0, v_0, v_1, v_2) f_T4_1222(T0, v_0, v_1, v_2) f_T4_2000(T0, v_0, v_1, v_2) f_T4_2001(T0, v_0, v_1, v_2) f_T4_2002(T0, v_0, v_1, v_2) f_T4_2010(T0, v_0, v_1, v_2) f_T4_2011(T0, v_0, v_1, v_2) f_T4_2012(T0, v_0, v_1, v_2) f_T4_2020(T0, v_0, v_1, v_2) f_T4_2021(T0, v_0, v_1, v_2) f_T4_2022(T0, v_0, v_1, v_2) f_T4_2100(T0, v_0, v_1, v_2) f_T4_2101(T0, v_0, v_1, v_2) f_T4_2102(T0, v_0, v_1, v_2) f_T4_2110(T0, v_0, v_1, v_2) f_T4_2111(T0, v_0, v_1, v_2) f_T4_2112(T0, v_0, v_1, v_2) f_T4_2120(T0, v_0, v_1, v_2) f_T4_2121(T0, v_0, v_1, v_2) f_T4_2122(T0, v_0, v_1, v_2) f_T4_2200(T0, v_0, v_1, v_2) f_T4_2201(T0, v_0, v_1, v_2) f_T4_2202(T0, v_0, v_1, v_2) f_T4_2210(T0, v_0, v_1, v_2) f_T4_2211(T0, v_0, v_1, v_2) f_T4_2212(T0, v_0, v_1, v_2) f_T4_2220(T0, v_0, v_1, v_2) f_T4_2221(T0, v_0, v_1, v_2) f_T4_2222(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function symmetric tensor (rank-2): f_ST2_00(T0, v_0, v_1, v_2) f_ST2_01(T0, v_0, v_1, v_2) f_ST2_02(T0, v_0, v_1, v_2) f_ST2_01(T0, v_0, v_1, v_2) f_ST2_11(T0, v_0, v_1, v_2) f_ST2_12(T0, v_0, v_1, v_2) f_ST2_02(T0, v_0, v_1, v_2) f_ST2_12(T0, v_0, v_1, v_2) f_ST2_22(T0, v_0, v_1, v_2)
+DEAL:dim 3::Symbol function symmetric tensor (rank-4): f_ST4_0000(T0, v_0, v_1, v_2) f_ST4_0001(T0, v_0, v_1, v_2) f_ST4_0002(T0, v_0, v_1, v_2) f_ST4_0001(T0, v_0, v_1, v_2) f_ST4_0011(T0, v_0, v_1, v_2) f_ST4_0012(T0, v_0, v_1, v_2) f_ST4_0002(T0, v_0, v_1, v_2) f_ST4_0012(T0, v_0, v_1, v_2) f_ST4_0022(T0, v_0, v_1, v_2) f_ST4_0100(T0, v_0, v_1, v_2) f_ST4_0101(T0, v_0, v_1, v_2) f_ST4_0102(T0, v_0, v_1, v_2) f_ST4_0101(T0, v_0, v_1, v_2) f_ST4_0111(T0, v_0, v_1, v_2) f_ST4_0112(T0, v_0, v_1, v_2) f_ST4_0102(T0, v_0, v_1, v_2) f_ST4_0112(T0, v_0, v_1, v_2) f_ST4_0122(T0, v_0, v_1, v_2) f_ST4_0200(T0, v_0, v_1, v_2) f_ST4_0201(T0, v_0, v_1, v_2) f_ST4_0202(T0, v_0, v_1, v_2) f_ST4_0201(T0, v_0, v_1, v_2) f_ST4_0211(T0, v_0, v_1, v_2) f_ST4_0212(T0, v_0, v_1, v_2) f_ST4_0202(T0, v_0, v_1, v_2) f_ST4_0212(T0, v_0, v_1, v_2) f_ST4_0222(T0, v_0, v_1, v_2) f_ST4_0100(T0, v_0, v_1, v_2) f_ST4_0101(T0, v_0, v_1, v_2) f_ST4_0102(T0, v_0, v_1, v_2) f_ST4_0101(T0, v_0, v_1, v_2) f_ST4_0111(T0, v_0, v_1, v_2) f_ST4_0112(T0, v_0, v_1, v_2) f_ST4_0102(T0, v_0, v_1, v_2) f_ST4_0112(T0, v_0, v_1, v_2) f_ST4_0122(T0, v_0, v_1, v_2) f_ST4_1100(T0, v_0, v_1, v_2) f_ST4_1101(T0, v_0, v_1, v_2) f_ST4_1102(T0, v_0, v_1, v_2) f_ST4_1101(T0, v_0, v_1, v_2) f_ST4_1111(T0, v_0, v_1, v_2) f_ST4_1112(T0, v_0, v_1, v_2) f_ST4_1102(T0, v_0, v_1, v_2) f_ST4_1112(T0, v_0, v_1, v_2) f_ST4_1122(T0, v_0, v_1, v_2) f_ST4_1200(T0, v_0, v_1, v_2) f_ST4_1201(T0, v_0, v_1, v_2) f_ST4_1202(T0, v_0, v_1, v_2) f_ST4_1201(T0, v_0, v_1, v_2) f_ST4_1211(T0, v_0, v_1, v_2) f_ST4_1212(T0, v_0, v_1, v_2) f_ST4_1202(T0, v_0, v_1, v_2) f_ST4_1212(T0, v_0, v_1, v_2) f_ST4_1222(T0, v_0, v_1, v_2) f_ST4_0200(T0, v_0, v_1, v_2) f_ST4_0201(T0, v_0, v_1, v_2) f_ST4_0202(T0, v_0, v_1, v_2) f_ST4_0201(T0, v_0, v_1, v_2) f_ST4_0211(T0, v_0, v_1, v_2) f_ST4_0212(T0, v_0, v_1, v_2) f_ST4_0202(T0, v_0, v_1, v_2) f_ST4_0212(T0, v_0, v_1, v_2) f_ST4_0222(T0, v_0, v_1, v_2) f_ST4_1200(T0, v_0, v_1, v_2) f_ST4_1201(T0, v_0, v_1, v_2) f_ST4_1202(T0, v_0, v_1, v_2) f_ST4_1201(T0, v_0, v_1, v_2) f_ST4_1211(T0, v_0, v_1, v_2) f_ST4_1212(T0, v_0, v_1, v_2) f_ST4_1202(T0, v_0, v_1, v_2) f_ST4_1212(T0, v_0, v_1, v_2) f_ST4_1222(T0, v_0, v_1, v_2) f_ST4_2200(T0, v_0, v_1, v_2) f_ST4_2201(T0, v_0, v_1, v_2) f_ST4_2202(T0, v_0, v_1, v_2) f_ST4_2201(T0, v_0, v_1, v_2) f_ST4_2211(T0, v_0, v_1, v_2) f_ST4_2212(T0, v_0, v_1, v_2) f_ST4_2202(T0, v_0, v_1, v_2) f_ST4_2212(T0, v_0, v_1, v_2) f_ST4_2222(T0, v_0, v_1, v_2)
+DEAL:dim 3::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.