Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines the class templates ProductOf2Maps and ProductOf3Maps.
6 :
7 : #pragma once
8 :
9 : #include <array>
10 : #include <cstddef>
11 : #include <optional>
12 : #include <utility>
13 :
14 : #include "DataStructures/Tensor/Tensor.hpp"
15 : #include "Utilities/MakeWithValue.hpp"
16 : #include "Utilities/TMPL.hpp"
17 :
18 : /// \cond
19 : namespace PUP {
20 : class er;
21 : } // namespace PUP
22 : /// \endcond
23 :
24 : namespace domain {
25 : namespace CoordinateMaps {
26 : /// \ingroup CoordinateMapsGroup
27 : /// \brief Product of two codimension=0 CoordinateMaps.
28 : ///
29 : /// \tparam Map1 the map for the first coordinate(s)
30 : /// \tparam Map2 the map for the second coordinate(s)
31 : template <typename Map1, typename Map2>
32 1 : class ProductOf2Maps {
33 : public:
34 0 : static constexpr size_t dim = Map1::dim + Map2::dim;
35 0 : using map_list = tmpl::list<Map1, Map2>;
36 : static_assert(dim == 2 or dim == 3,
37 : "Only 2D and 3D maps are supported by ProductOf2Maps");
38 :
39 : // Needed for Charm++ serialization
40 0 : ProductOf2Maps() = default;
41 :
42 0 : ProductOf2Maps(Map1 map1, Map2 map2);
43 :
44 : template <typename T>
45 0 : std::array<T, dim> operator()(const std::array<T, dim>& source_coords) const;
46 :
47 : /// The inverse function is only callable with doubles because the inverse
48 : /// might fail if called for a point out of range, and it is unclear
49 : /// what should happen if the inverse were to succeed for some points in a
50 : /// DataVector but fail for other points.
51 1 : std::optional<std::array<double, dim>> inverse(
52 : const std::array<double, dim>& target_coords) const;
53 :
54 : template <typename T>
55 0 : tnsr::Ij<T, dim, Frame::NoFrame> inv_jacobian(
56 : const std::array<T, dim>& source_coords) const;
57 :
58 : template <typename T>
59 0 : tnsr::Ij<T, dim, Frame::NoFrame> jacobian(
60 : const std::array<T, dim>& source_coords) const;
61 :
62 : // NOLINTNEXTLINE(google-runtime-references)
63 0 : void pup(PUP::er& p);
64 :
65 0 : bool is_identity() const { return is_identity_; }
66 :
67 0 : static constexpr bool supports_hessian{Map1::supports_hessian and
68 : Map2::supports_hessian};
69 :
70 : private:
71 0 : friend bool operator==(const ProductOf2Maps& lhs, const ProductOf2Maps& rhs) {
72 : return lhs.map1_ == rhs.map1_ and lhs.map2_ == rhs.map2_ and
73 : lhs.is_identity_ == rhs.is_identity_;
74 : }
75 :
76 0 : Map1 map1_;
77 0 : Map2 map2_;
78 0 : bool is_identity_ = false;
79 : };
80 :
81 : template <typename Map1, typename Map2>
82 0 : bool operator!=(const ProductOf2Maps<Map1, Map2>& lhs,
83 : const ProductOf2Maps<Map1, Map2>& rhs);
84 :
85 : /// \ingroup CoordinateMapsGroup
86 : /// \brief Product of three one-dimensional CoordinateMaps.
87 : template <typename Map1, typename Map2, typename Map3>
88 1 : class ProductOf3Maps {
89 : public:
90 0 : static constexpr size_t dim = Map1::dim + Map2::dim + Map3::dim;
91 0 : using map_list = tmpl::list<Map1, Map2, Map3>;
92 : static_assert(dim == 3, "Only 3D maps are implemented for ProductOf3Maps");
93 :
94 : // Needed for Charm++ serialization
95 0 : ProductOf3Maps() = default;
96 :
97 0 : ProductOf3Maps(Map1 map1, Map2 map2, Map3 map3);
98 :
99 : template <typename T>
100 0 : std::array<T, dim> operator()(const std::array<T, dim>& source_coords) const;
101 :
102 0 : std::optional<std::array<double, dim>> inverse(
103 : const std::array<double, dim>& target_coords) const;
104 :
105 : template <typename T>
106 0 : tnsr::Ij<T, dim, Frame::NoFrame> inv_jacobian(
107 : const std::array<T, dim>& source_coords) const;
108 :
109 : template <typename T>
110 0 : tnsr::Ij<T, dim, Frame::NoFrame> jacobian(
111 : const std::array<T, dim>& source_coords) const;
112 :
113 : // NOLINTNEXTLINE(google-runtime-references)
114 0 : void pup(PUP::er& p);
115 :
116 0 : bool is_identity() const { return is_identity_; }
117 :
118 0 : static constexpr bool supports_hessian{Map1::supports_hessian and
119 : Map2::supports_hessian and
120 : Map3::supports_hessian};
121 :
122 : private:
123 0 : friend bool operator==(const ProductOf3Maps& lhs, const ProductOf3Maps& rhs) {
124 : return lhs.map1_ == rhs.map1_ and lhs.map2_ == rhs.map2_ and
125 : lhs.map3_ == rhs.map3_ and lhs.is_identity_ == rhs.is_identity_;
126 : }
127 :
128 0 : Map1 map1_;
129 0 : Map2 map2_;
130 0 : Map3 map3_;
131 0 : bool is_identity_ = false;
132 : };
133 :
134 : template <typename Map1, typename Map2, typename Map3>
135 0 : bool operator!=(const ProductOf3Maps<Map1, Map2, Map3>& lhs,
136 : const ProductOf3Maps<Map1, Map2, Map3>& rhs);
137 : } // namespace CoordinateMaps
138 : } // namespace domain
|
