Line data Source code
1 0 : // Distributed under the MIT License. 2 : // See LICENSE.txt for details. 3 : 4 : #pragma once 5 : 6 : #include <string> 7 : 8 : #include "DataStructures/DataBox/Tag.hpp" 9 : #include "DataStructures/DataBox/TagName.hpp" 10 : #include "DataStructures/DataVector.hpp" 11 : #include "DataStructures/Tensor/Tensor.hpp" 12 : #include "Utilities/Functional.hpp" 13 : #include "Utilities/Gsl.hpp" 14 : #include "Utilities/Numeric.hpp" 15 : #include "Utilities/TMPL.hpp" 16 : 17 : namespace L2Norm_detail { 18 : template <typename DataType, typename Symm, typename IndexList> 19 : Scalar<DataType> pointwise_l2_norm_square( 20 : const Tensor<DataType, Symm, IndexList>& tensor) { 21 : auto pointwise_l2_normsq = make_with_value<Scalar<DataType>>(tensor, 0.); 22 : for (auto tensor_element = tensor.begin(); tensor_element != tensor.end(); 23 : ++tensor_element) { 24 : // In order to handle tensory symmetries, we multiply the square of each 25 : // stored component with its multiplicity 26 : get(pointwise_l2_normsq) += 27 : tensor.multiplicity(tensor_element) * square(*tensor_element); 28 : } 29 : return pointwise_l2_normsq; 30 : } 31 : } // namespace L2Norm_detail 32 : 33 : /// @{ 34 : /*! 35 : * \ingroup TensorGroup 36 : * \brief Compute point-wise Euclidean \f$L^2\f$-norm of arbitrary Tensors. 37 : * 38 : * \details 39 : * At each grid point \f$p\f$ in the element, this function computes the 40 : * point-wise Frobenius norm of a given Tensor with arbitrary rank. If the 41 : * Tensor \f$A\f$ has rank \f$n\f$ and dimensionality \f$D\f$, then its 42 : * Frobenius norm at point \f$p\f$ is computed as: 43 : * 44 : * \f{equation} 45 : * ||A||_2(p) = 46 : * \left(\sum^{D-1}_{i_1=0}\sum^{D-1}_{i_2=0}\cdots \sum^{D-1}_{i_n=0} 47 : * |A_{i_1 i_2 \cdots i_n}(p)|^2 \right)^{1/2}, 48 : * \f} 49 : * 50 : * where both contra-variant and co-variant indices are shown as lower indices. 51 : */ 52 : template <typename DataType, typename Symm, typename IndexList> 53 1 : Scalar<DataType> pointwise_l2_norm( 54 : const Tensor<DataType, Symm, IndexList>& tensor) { 55 : return Scalar<DataType>{ 56 : sqrt(get(L2Norm_detail::pointwise_l2_norm_square(tensor)))}; 57 : } 58 : template <typename DataType, typename Symm, typename IndexList> 59 1 : void pointwise_l2_norm(const gsl::not_null<Scalar<DataType>*> norm, 60 : const Tensor<DataType, Symm, IndexList>& tensor) { 61 : get(*norm) = sqrt(get(L2Norm_detail::pointwise_l2_norm_square(tensor))); 62 : } 63 : /// @} 64 : 65 : /*! 66 : * \ingroup TensorGroup 67 : * \brief Compute Euclidean \f$L^2\f$-norm of arbitrary Tensors reduced over an 68 : * element. 69 : * 70 : * \details 71 : * Computes the RMS value of the point-wise Frobenius norm of a given Tensor 72 : * with arbitrary rank over all grid points in an element. If the Tensor \f$A\f$ 73 : * has rank \f$n\f$ and dimensionality \f$D\f$, and the element (of order 74 : * \f$N\f$) has \f$N+1\f$ points, then its element-reduced Frobenius norm is 75 : * computed as: 76 : * 77 : * \f{equation} 78 : * ||A||_2 = 79 : * \left(\frac{1}{N+1}\sum^{N}_{p=0} 80 : * \left(\sum^{D-1}_{i_1=0}\sum^{D-1}_{i_2=0}\cdots 81 : * \sum^{D-1}_{i_n=0} |A^p_{i_1 i_2 \cdots i_n}|^2 \right) 82 : * \right)^{1/2}, 83 : * \f} 84 : * 85 : * where both contra-variant and co-variant indices are shown as lower indices, 86 : * and \f$p\f$ indexes grid points in the element. 87 : * 88 : * \warning This function reduces the Frobenius norm over the element, not the 89 : * whole domain. 90 : */ 91 : template <typename DataType, typename Symm, typename IndexList> 92 1 : double l2_norm(const Tensor<DataType, Symm, IndexList>& tensor) { 93 : const auto pointwise_l2_normsq = 94 : L2Norm_detail::pointwise_l2_norm_square(tensor); 95 : using Plus = funcl::Plus<funcl::Identity>; 96 : return sqrt(alg::accumulate(get(pointwise_l2_normsq), 0., Plus{}) / 97 : tensor.begin()->size()); 98 : } 99 : 100 : namespace Tags { 101 : /*! 102 : * \ingroup DataBoxTagsGroup 103 : * \ingroup DataStructuresGroup 104 : * Point-wise Euclidean \f$L^2\f$-norm of a Tensor. 105 : * \see `pointwise_l2_norm()` for details. 106 : */ 107 : template <typename Tag> 108 1 : struct PointwiseL2Norm : db::SimpleTag { 109 0 : using type = Scalar<typename Tag::type::type>; 110 0 : static std::string name() { 111 : return "PointwiseL2Norm(" + db::tag_name<Tag>() + ")"; 112 : } 113 : }; 114 : 115 : /*! 116 : * \ingroup DataBoxTagsGroup 117 : * \ingroup DataStructuresGroup 118 : * Computes the point-wise Euclidean \f$L^2\f$-norm of a Tensor. 119 : * \see `pointwise_l2_norm()` for details. 120 : */ 121 : template <typename Tag> 122 1 : struct PointwiseL2NormCompute : PointwiseL2Norm<Tag>, db::ComputeTag { 123 0 : using base = PointwiseL2Norm<Tag>; 124 0 : using return_type = typename base::type; 125 0 : static constexpr auto function = 126 : static_cast<void (*)(const gsl::not_null<return_type*>, 127 : const typename Tag::type&)>(&pointwise_l2_norm); 128 0 : using argument_tags = tmpl::list<Tag>; 129 : }; 130 : 131 : /*! 132 : * \ingroup DataBoxTagsGroup 133 : * \ingroup DataStructuresGroup 134 : * Euclidean \f$L^2\f$-norm of a Tensor, RMS over grid points in element. 135 : * \see `l2_norm()` for details. 136 : */ 137 : template <typename Tag> 138 1 : struct L2Norm : db::SimpleTag { 139 0 : using type = double; 140 0 : static std::string name() { return "L2Norm(" + db::tag_name<Tag>() + ")"; } 141 : }; 142 : 143 : /*! 144 : * \ingroup DataBoxTagsGroup 145 : * \ingroup DataStructuresGroup 146 : * Computes the Euclidean \f$L^2\f$-norm of a Tensor, RMS over grid points in 147 : * element. \see `l2_norm()` for details. 148 : * 149 : * \warning This compute tag reduces the Frobenius norm over the element, not 150 : * the whole domain. 151 : */ 152 : template <typename Tag> 153 1 : struct L2NormCompute : L2Norm<Tag>, db::ComputeTag { 154 0 : using base = L2Norm<Tag>; 155 0 : using return_type = typename base::type; 156 0 : static void function(const gsl::not_null<return_type*> norm, 157 : const typename Tag::type& tensor) { 158 : *norm = l2_norm(tensor); 159 : } 160 0 : using argument_tags = tmpl::list<Tag>; 161 : }; 162 : } // namespace Tags