template <typename T> T pow21(T x) { return x * x; }

template <typename T> constexpr T pow21(T x) noexcept { return x * x; }

  reduction();
  reduction(const vect<T, D> &p, const T &V, const T &x);
  reduction(const reduction &x, const reduction &y);
  reduction operator+(const reduction &y) const { return reduction(*this, y); }
  reduction &operator+=(const reduction &y) { return *this = *this + y; }

  T avg() const { return sum / vol; }
  T sdv() const { return sqrt(fmax1(T(0), vol * sum2 - pow21(sum))); }
  T norm0() const { return vol; }
  T norm1() const { return sumabs / vol; }
  T norm2() const { return sqrt(sum2abs / vol); }
  T norm_inf() const { return maxabs; }

  constexpr reduction();
  constexpr reduction(const vect<T, D> &p, const T &V, const T &x);
  constexpr reduction(const reduction &x, const reduction &y);
  constexpr reduction operator+(const reduction &y) const {
    return reduction(*this, y);
  }
  constexpr reduction &operator+=(const reduction &y) noexcept {
    return *this = *this + y;
  }

  constexpr T avg() const noexcept { return sum / vol; }
  constexpr T sdv() const noexcept {
    return sqrt(fmax1(T(0), vol * sum2 - pow21(sum)));
  }
  constexpr T norm0() const noexcept { return vol; }
  constexpr T norm1() const noexcept { return sumabs / vol; }
  constexpr T norm2() const noexcept { return sqrt(sum2abs / vol); }
  constexpr T norm_inf() const noexcept { return maxabs; }

reduction<T, D>::reduction()

constexpr reduction<T, D>::reduction()

reduction<T, D>::reduction(const vect<T, D> &p, const T &V, const T &x)

constexpr reduction<T, D>::reduction(const vect<T, D> &p, const T &V,
                                     const T &x)

reduction<T, D>::reduction(const reduction &x, const reduction &y)

constexpr reduction<T, D>::reduction(const reduction &x, const reduction &y)