#pragma once
#include <cassert>
#include <iomanip>             // setw()
#include <iostream>
#include <omp.h>
#include <vector>

/**     Inner product
    @param[in] x    vector
    @param[in] y    vector
    @return         resulting Euclidian inner product <x,y>
*/
double scalar(std::vector<double> const &x, std::vector<double> const &y);
double scalar_trans(std::vector<double> const &x, std::vector<double> const &y);

/**  Second inner product: use #pragma omp parallel (no "for" in pragma)
    The work is split manually inside the parallel region (stride or chunk).
*/
double scalar_manual(std::vector<double> const &x, std::vector<double> const &y);

/**     l2-norm
    @param[in] x    vector
    @return         resulting Euclidian norm
*/
double norm(std::vector<double> const &x);

/**      Vector @p b adds its elements to vector @p a .
    @param[in] a    vector
    @param[in] b    vector
    @return          a+=b componentwise
*/
template<class T>
std::vector<T> &operator+=(std::vector<T> &a, std::vector<T> const &b)
{
    assert(a.size()==b.size());
    for (size_t k = 0; k < a.size(); ++k) {
        a[k] += b[k];
    }
    return a;
}

// Declare the reduction operation in OpenMP for an STL-vector (existing)
#pragma omp declare reduction(VecAdd : std::vector<int>  : omp_out += omp_in) \
 initializer (omp_priv=omp_orig)

/**     Test for vector reduction.
 *  existing: converts thread-private vectors into a componentwise sum
 *  @param[in] n  size of global/private vector
 *  @return  resulting global vector.
*/
std::vector<int> reduction_vec(int n);

/**     New: append per-thread vectors into a single big vector.
 *  The result will have size n * numThreads, where each thread contributes a contiguous block.
 */
std::vector<int> reduction_vec_append(int n);

/**     Output of a vector.
    @param[in,out] s    output stream
    @param[in]     x    vector
    @return        modified output stream
*/
template <class T>
std::ostream &operator<<(std::ostream &s, std::vector<T> const &x)
{
    for (auto const &v : x)  s << std::setw(4) << v << "  ";
    return s;
}