From f6167fac2b8359415abf70b40c4b67a091b6f65c Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Sat, 27 Apr 2019 10:46:39 +0200 Subject: [PATCH] Add functions to add tensor types to SD subsitution maps --- .../sd/symengine_tensor_operations.h | 143 ++++++++++++++++++ .../sd/symengine_tensor_operations.cc | 31 ++-- .../symengine/substitution_maps_tensor_01.cc | 62 ++++++++ .../substitution_maps_tensor_01.output | 17 +++ 4 files changed, 237 insertions(+), 16 deletions(-) create mode 100644 tests/symengine/substitution_maps_tensor_01.cc create mode 100644 tests/symengine/substitution_maps_tensor_01.output diff --git a/include/deal.II/differentiation/sd/symengine_tensor_operations.h b/include/deal.II/differentiation/sd/symengine_tensor_operations.h index 22ab4c7dde..c325bf9cc8 100644 --- a/include/deal.II/differentiation/sd/symengine_tensor_operations.h +++ b/include/deal.II/differentiation/sd/symengine_tensor_operations.h @@ -366,6 +366,61 @@ namespace Differentiation //@} + /** + * @name Symbol substitution map enlargement + */ + //@{ + + /** + * A convenience function for adding an entry to the @p substitution_map. + * The new entries will have the keys given in the @p symbol_tensor with + * their paired values extracted from the corresponding elements of the + * @p value_tensor. + * + * For more context which this function is used, see the other + * \ref add_to_substitution_map(types::substitution_map &,const + * Expression &, const Expression &) function. + * + * @tparam ignore_invalid_symbols See the other + * \ref add_to_substitution_map(types::substitution_map &,const + * Expression &,const Expression &) function for a detailed + * discussion on the role of this template argument. + */ + template + void + add_to_substitution_map(types::substitution_map &substitution_map, + const Tensor &symbol_tensor, + const Tensor & value_tensor); + + /** + * A convenience function for adding an entry to the @p substitution_map. + * The new entries will have the keys given in the @p symbol_tensor with + * their paired values extracted from the corresponding elements of the @p value_tensor. + * + * For more context which this function is used, see the other + * \ref add_to_substitution_map(types::substitution_map &,const + * Expression &, const Expression &) function. + * + * @tparam ignore_invalid_symbols See the other + * \ref add_to_substitution_map(types::substitution_map &,const + * Expression &,const Expression &) function for a detailed + * discussion on the role of this template argument. + */ + template + void + add_to_substitution_map( + types::substitution_map & substitution_map, + const SymmetricTensor &symbol_tensor, + const SymmetricTensor & value_tensor); + + //@} + } // namespace SD } // namespace Differentiation @@ -384,6 +439,12 @@ namespace Differentiation namespace internal { + template + TableIndices<4> + make_rank_4_tensor_indices(const unsigned int &idx_i, + const unsigned int &idx_j); + + template TableIndices concatenate_indices(const TableIndices &indices_1, @@ -749,6 +810,88 @@ namespace Differentiation return internal::tensor_diff_tensor(symbol_tensor, op); } + + /* ---------------- Symbolic substitution map enlargement --------------*/ + + + namespace internal + { + template class TensorType> + std::vector> + tensor_substitution_map( + const TensorType &symbol_tensor, + const TensorType & value_tensor) + { + std::vector> symbol_values; + for (unsigned int i = 0; i < symbol_tensor.n_independent_components; + ++i) + { + const TableIndices indices( + symbol_tensor.unrolled_to_component_indices(i)); + symbol_values.push_back( + std::make_pair(symbol_tensor[indices], value_tensor[indices])); + } + return symbol_values; + } + + + template + std::vector> + tensor_substitution_map( + const SymmetricTensor<4, dim, Expression> &symbol_tensor, + const SymmetricTensor<4, dim, ValueType> & value_tensor) + { + std::vector> symbol_values; + 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(i, j); + symbol_values.push_back( + std::make_pair(symbol_tensor[indices], value_tensor[indices])); + } + return symbol_values; + } + } // namespace internal + + + template + void + add_to_substitution_map(types::substitution_map &substitution_map, + const Tensor &symbol_tensor, + const Tensor & value_tensor) + { + add_to_substitution_map( + substitution_map, + internal::tensor_substitution_map(symbol_tensor, value_tensor)); + } + + + template + void + add_to_substitution_map( + types::substitution_map & substitution_map, + const SymmetricTensor &symbol_tensor, + const SymmetricTensor & value_tensor) + { + add_to_substitution_map( + substitution_map, + internal::tensor_substitution_map(symbol_tensor, value_tensor)); + } + } // namespace SD } // namespace Differentiation diff --git a/source/differentiation/sd/symengine_tensor_operations.cc b/source/differentiation/sd/symengine_tensor_operations.cc index bd2e743623..e07a3b4911 100644 --- a/source/differentiation/sd/symengine_tensor_operations.cc +++ b/source/differentiation/sd/symengine_tensor_operations.cc @@ -38,24 +38,23 @@ namespace Differentiation namespace internal { - namespace + template + TableIndices<4> + make_rank_4_tensor_indices(const unsigned int &idx_i, + const unsigned int &idx_j) { - template - 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]); - } - + 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]); + } + namespace + { template std::string make_index_string(const TableIndices &indices) diff --git a/tests/symengine/substitution_maps_tensor_01.cc b/tests/symengine/substitution_maps_tensor_01.cc new file mode 100644 index 0000000000..c043a6bf35 --- /dev/null +++ b/tests/symengine/substitution_maps_tensor_01.cc @@ -0,0 +1,62 @@ +// --------------------------------------------------------------------- +// +// 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 functions to add tensor variables to substitution maps work +// correctly. + +#include + +#include + +#include "../tests.h" + +using namespace dealii; +namespace SD = Differentiation::SD; +namespace SE = ::SymEngine; + +int +main() +{ + initlog(); + const unsigned int dim = 2; + + Tensor<2, dim> t; + SymmetricTensor<2, dim> st; + for (unsigned int i = 0; i < t.n_independent_components; ++i) + t[t.unrolled_to_component_indices(i)] = i; + for (unsigned int i = 0; i < st.n_independent_components; ++i) + st[st.unrolled_to_component_indices(i)] = i; + + + SD::types::substitution_map substitution_map; + SD::add_to_substitution_map(substitution_map, + SD::Expression("x1"), + SD::Expression(1)); + SD::add_to_substitution_map(substitution_map, + SD::make_tensor_of_symbols<2, dim>("t1"), + t); + SD::add_to_substitution_map( + substitution_map, SD::make_symmetric_tensor_of_symbols<2, dim>("st1"), st); + + SD::add_to_substitution_map( + substitution_map, + std::make_pair(SD::make_tensor_of_symbols<2, dim>("t2"), t), + std::make_pair(SD::make_symmetric_tensor_of_symbols<2, dim>("st2"), st)); + + SD::Utilities::print_substitution_map(deallog, substitution_map); + + deallog << "OK" << std::endl; +} diff --git a/tests/symengine/substitution_maps_tensor_01.output b/tests/symengine/substitution_maps_tensor_01.output new file mode 100644 index 0000000000..c562f76802 --- /dev/null +++ b/tests/symengine/substitution_maps_tensor_01.output @@ -0,0 +1,17 @@ + +DEAL::x1 = 1 +t2_01 = 1.0 +t2_00 = 0.0 +t2_10 = 2.0 +t2_11 = 3.0 +t1_00 = 0.0 +t1_01 = 1.0 +t1_11 = 3.0 +t1_10 = 2.0 +st2_11 = 1.0 +st2_00 = 0.0 +st2_01 = 2.0 +st1_11 = 1.0 +st1_00 = 0.0 +st1_01 = 2.0 +OK -- 2.39.5