#include <bits/stdc++.h>
#include <tbb/tick_count.h>

using namespace std;

const int N = 4e6;
double x[2*N], y[2*N];

void saxpy(double a) {
        for (int i = 0 ; i < 2*N-1 ; i+=2)
                y[i] = a * y[i+1] + y[i];
}

void saxpy_parallel(double a, int kParallelism) {
  auto fun = [=](double a, int st, int en) {
    for (int i = st; i < en; i+=2)
      y[i] = a * y[i+1] * y[i];
  };
  {
    auto parallelizer = [=](decltype(fun) func, int st, int en) {
      std::thread bmvm_threads[kParallelism - 1];
      int tasksize = en - st;
      int batchsize = tasksize / kParallelism;
      for (int i = 0; i < kParallelism - 1; i++)
        bmvm_threads[i] =
            std::thread(func, a, st + i * batchsize, st + (i + 1) * batchsize);
      func(a, st + (kParallelism - 1) * batchsize, en);
      for (int i = 0; i < kParallelism - 1; i++)
        bmvm_threads[i].join();
    };
    parallelizer(fun, 0, 2*N-1);
  }

}

int main() {
        std::uniform_real_distribution<double> unif(0, N);
        std::default_random_engine re;
        for (int i = 0 ; i < N ; i++)
                y[i] = unif(re), x[i] = unif(re);
        std::chrono::time_point<std::chrono::high_resolution_clock> start;
        start = std::chrono::high_resolution_clock::now();
        saxpy(4.0);
        std::chrono::duration<double> time_taken =
                std::chrono::high_resolution_clock::now() - start;
        cout << "Serial Time: " << time_taken.count() << endl;
        start = std::chrono::high_resolution_clock::now();
        saxpy_parallel(4.0, 2);
        time_taken = std::chrono::high_resolution_clock::now() - start;
        cout << "Parallel Time: " << time_taken.count() << endl;
        return 0;
}