mapped_functions.hpp

Go to the documentation of this file.

// Copyright (c) 2019-2023 Simons Foundation
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
//     http://www.apache.org/licenses/LICENSE-2.0.txt
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
//
// Authors: Olivier Parcollet, Nils Wentzell
 
#pragma once
 
#include "./concepts.hpp"
#include "./map.hpp"
#include "./traits.hpp"
 
#include <cmath>
#include <complex>
#include <type_traits>
#include <utility>
 
namespace nda {
 
  template <typename T>
  auto real(T t)
    requires(nda::is_scalar_v<T>)
  {
    if constexpr (is_complex_v<T>) {
      return std::real(t);
    } else {
      return t;
    }
  }
 
  template <typename T>
  auto conj(T t)
    requires(nda::is_scalar_v<T>)
  {
    if constexpr (is_complex_v<T>) {
      return std::conj(t);
    } else {
      return t;
    }
  }
 
  inline double abs2(double x) { return x * x; }
 
  inline double abs2(std::complex<double> z) { return (conj(z) * z).real(); }
 
  inline bool isnan(std::complex<double> const &z) { return std::isnan(z.real()) or std::isnan(z.imag()); }
 
  template <typename T>
  T pow(T x, int n)
    requires(std::is_integral_v<T>)
  {
    T r = 1;
    for (int i = 0; i < n; ++i) r *= x;
    return r;
  }
 
  template <Array A>
  auto pow(A &&a, double p) {
    return nda::map([p](auto const &x) {
      using std::pow;
      return pow(x, p);
    })(std::forward<A>(a));
  }
 
  struct conj_f {
    auto operator()(auto const &x) const { return conj(x); };
  };
 
  template <Array A>
  decltype(auto) conj(A &&a) {
    if constexpr (is_complex_v<get_value_t<A>>)
      return nda::map(conj_f{})(std::forward<A>(a));
    else
      return std::forward<A>(a);
  }
 
} // namespace nda