norm.hpp

Go to the documentation of this file.

// SPDX-License-Identifier: BSD-2-Clause
 
#ifndef pRC_TENSOR_TRAIN_COMMON_FUNCTIONS_NORM_H
#define pRC_TENSOR_TRAIN_COMMON_FUNCTIONS_NORM_H
 
#include <prc/core/basic/limits.hpp>
#include <prc/core/functors/view.hpp>
#include <prc/core/tensor/functions/reshape.hpp>
#include <prc/tensor_train/algorithms/orthogonalize.hpp>
#include <prc/tensor_train/operator/type_traits.hpp>
#include <prc/tensor_train/tensor/type_traits.hpp>
 
namespace pRC
{
    template<Index P = 2, Index Q = P, class X, class R = RemoveReference<X>,
        If<Any<TensorTrain::IsTensorish<R>, TensorTrain::IsOperatorish<R>>> = 0,
        If<IsSatisfied<(P == 2 && Q == 2)>> = 0>
    static inline constexpr auto norm(X &&a)
    {
        using T = Tensor<typename R::Value>;
 
        return exp(expand(
            makeSeries<Index, typename RemoveReference<X>::Dimension{}>(),
            [&a](auto const... seq)
            {
                return RecursiveLambda(
                    []<class XA, class XB>(auto const &self, XA &&a, XB &&b,
                        auto &&...cores) -> T
                    {
                        auto const [nrm, q, r] =
                            orthogonalize<Position::Left>(forward<XA>(a));
 
                        if(nrm <= identity<T>(
                                      NumericLimits<typename T::Value>::min()))
                        {
                            return zero();
                        }
 
                        Tensor const next = contract<1, 0>(r, forward<XB>(b));
 
                        if constexpr(sizeof...(cores) == 0)
                        {
                            return log(nrm) + log(norm(next));
                        }
                        else
                        {
                            return log(nrm) +
                                self(next, forward<decltype(cores)>(cores)...);
                        }
                    })(forward<X>(a).template core<seq>()...);
            }));
    }
}
#endif // pRC_TENSOR_TRAIN_COMMON_FUNCTIONS_NORM_H