#include <algorithm>
#include <iterator>
#include <iostream>
#include <fstream>
#include <numeric>
#include <cstdlib>
#include <vector>
#include <chrono>
#include <iomanip>
#include "assert.h"

using namespace std;

struct Complex
{
    double re;
    double im;
};

Complex operator+(Complex x, Complex y)
{
    Complex ret_val = { x.re + y.re, x.im + y.im };
    return ret_val;
}

Complex operator*(Complex x, Complex y)
{
    Complex ret_val = { x.re*y.re - x.im*y.im, x.re*y.im + x.im*y.re };
    return ret_val;
}

Complex random_complex()
{
    Complex ret_val = { rand() / 1000. - 1000., rand() / 1000. - 1000. };
    return ret_val;
}

template<typename F>
void benchmark(F f)
{
    using namespace chrono;
    auto start = high_resolution_clock::now();
    f();
    auto end = high_resolution_clock::now();
    auto elapsed = end - start;
    cout << "Elapsed = " << duration_cast<duration<double>>(elapsed).count() * 1000 << " ms" << endl;
}

basic_ostream<char> & operator<<(basic_ostream<char> &stm, Complex c)
{
    stm << c.re << " " << c.im << std::endl;
    return stm;
}

basic_istream<char> & operator>>(basic_istream<char> &stm, Complex &c)
{
    stm >> c.re;
    stm >> c.im;
    return stm;
}

void generate_complex_sequence(vector<Complex> &v, vector<Complex> &u)
{
    constexpr auto N = 1u << 24;
    v.reserve(N);
    generate_n(back_inserter(v), N, random_complex);
    u = v;
    random_shuffle(begin(v), end(v));
    random_shuffle(begin(u), end(u));
}

void save_complex_sequence(vector<Complex> &v, vector<Complex> &u)
{
    ofstream stm("complex.dat", ios_base::out);
    assert(stm);
    stm << std::setprecision(14);
    ostream_iterator<Complex> it(stm);
    std::copy(v.begin(), v.end(), it);
    std::copy(u.begin(), u.end(), it);
}

void load_complex_sequence(vector<Complex> &v, vector<Complex> &u)
{
    constexpr auto N = 1u << 24;
    constexpr auto buf_len = 5u << 20;
    char *buffer = new char[buf_len];
    v.reserve(N);
    u.reserve(N);
    ifstream stm;
    stm.rdbuf()->pubsetbuf(buffer, buf_len);
    stm.open("complex.dat", ios_base::in);
    assert(stm);
    stm.precision(16);
    std::string str;
    unsigned i = 0;
    vector<Complex> *cur = &v;
    while(std::getline(stm, str))
    {
        char *ptr = nullptr;
        double re = ::strtod(str.data(), &ptr);
        double im = ::strtod(ptr, &ptr);
        Complex c = {re, im};
        cur->push_back(c);
        i++;
        if(i == N)
            cur = &u;
    }
    // istream_iterator<Complex> it(stm);
    // std::copy_n(it, N, back_inserter(v)); 
    // std::copy_n(it, N, back_inserter(u)); 
    std::cout << "v.size() = " << v.size() << std::endl;
    std::cout << "u.size() = " << u.size() << std::endl;
}

int main()
{
    vector<Complex> v, u;
    std::cout << "Benchmark loading:" << std::endl;
    // load_complex_sequence(v, u);
    generate_complex_sequence(v, u);

    benchmark([&]
    {
        transform(begin(v), end(v), begin(u), begin(v), [](Complex x, Complex y) { return x*y; });
    });

    auto result = accumulate(begin(v), end(v), Complex{});
    std::cout << std::setprecision(16) << " Re: " << result.re << " Im: " << result.im << std::endl;
    return 0;
}

using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.IO;
using System.Linq;
using System.Runtime.InteropServices;
using System.Text;

namespace ConsoleApplication1
{
    struct Complex
    {
        public Complex(double re, double im)
        {
            Re = re;
            Im = im;
        }

        public double Re;
        public double Im;

        public static Complex operator +(Complex x, Complex y)
        {
            return new Complex(x.Re + y.Re, x.Im + y.Im);
        }

        public unsafe static Complex operator *(Complex x, Complex y)
        {
            return new Complex(x.Re * y.Re - x.Im * y.Im, x.Re * y.Im + x.Im * y.Re);
        }

        private const int RandMax = 32767;
        private static Random random = new Random();

        public static Complex RandomComplex()
        {
            return new Complex(random.Next(RandMax) / 1000D, random.Next(RandMax) / 1000D);
        }
    }

    static class Funcs
    {
        public static void Shuffle<T>(IList<T> list)
        {
            Random rng = new Random();
            int n = list.Count;
            while (n > 1)
            {
                n--;
                int k = rng.Next(n + 1);
                T value = list[k];
                list[k] = list[n];
                list[n] = value;
            }
        }

        public static void PopulateWithRandom(Complex[] v, Complex[] u)
        {
            const int n = 1 << 24;
            for (int i = 0; i < n; i++ )
            {
                Complex random = Complex.RandomComplex();
                v[i] = random;
                u[i] = random;
            }

            Funcs.Shuffle(u);
            Funcs.Shuffle(v);
        }

        public static void LoadComplexSequence(Complex[] v, Complex[] u)
        {
            const int n = 1 << 24;

            using (var stm = File.OpenRead(@".\NativeBenchmark\complex.dat"))
            {
                using (var tr = new StreamReader(stm))
                {
                    foreach(var list in new Complex[][]{v, u})
                    {
                        int i = 0;
                        string line;
                        while (i < n && (line = tr.ReadLine()) != null)
                        {
                            var components = line.Split(' ');

                            if (components.Length != 2)
                                throw new ApplicationException("There should be two components");

                            list[i] = new Complex(Convert.ToDouble(components[0]), Convert.ToDouble(components[1]));
                            i++;
                        }

                        Console.WriteLine("Added {0} elements", i);
                    }
                }
            }
        }
    }

    class Program
    {
        [DllImport("winmm.dll")]
        internal static extern uint timeBeginPeriod(uint period);
        [DllImport("winmm.dll")]
        internal static extern uint timeEndPeriod(uint period);

        static void Main(string[] args)
        {

            const int n = 1 << 24;
            Console.WriteLine("n={0}", n);
            var v = new Complex[n];
            var u = new Complex[n];

            var watch = Stopwatch.StartNew();
            Funcs.PopulateWithRandom(v, u);
            //Funcs.LoadComplexSequence(v, u);
            watch.Stop();

            Console.WriteLine("Loading time: {0:F4}", watch.ElapsedMilliseconds);

            timeBeginPeriod(5);
            unsafe
            {
                watch = Stopwatch.StartNew();
                for (int i = 0; i < n; i++)
                {
                    v[i] = v[i] * u[i];
                }
            }
            watch.Stop();
            timeEndPeriod(5);

            var result = new Complex(0D, 0D);
            v.ToList().ForEach(c=>{result = result + c;});
            Console.WriteLine(".NET Elapsed: {0:F4} ms", watch.ElapsedMilliseconds);
            Console.WriteLine("Result: {0}, {1}", result.Re, result.Im);

        }
    }
}

.L94:
    vmovupd    (%rcx), %xmm3
    addq    $1, %r10
    addq    $64, %rcx
    addq    $64, %rdi
    vinsertf128    $0x1, -48(%rcx), %ymm3, %ymm1
    vmovupd    -32(%rcx), %xmm3
    vinsertf128    $0x1, -16(%rcx), %ymm3, %ymm3
    vmovupd    -64(%rdi), %xmm4
    vinsertf128    $1, %xmm3, %ymm1, %ymm0
    vperm2f128    $49, %ymm3, %ymm1, %ymm3
    vunpcklpd    %ymm3, %ymm0, %ymm2
    vunpckhpd    %ymm3, %ymm0, %ymm3
    vinsertf128    $0x1, -48(%rdi), %ymm4, %ymm0
    vmovupd    -32(%rdi), %xmm4
    vinsertf128    $0x1, -16(%rdi), %ymm4, %ymm4
    vinsertf128    $1, %xmm4, %ymm0, %ymm5
    vperm2f128    $49, %ymm4, %ymm0, %ymm4
    vunpcklpd    %ymm4, %ymm5, %ymm1
    vunpckhpd    %ymm4, %ymm5, %ymm4
    vmulpd    %ymm1, %ymm2, %ymm0
    vmulpd    %ymm4, %ymm3, %ymm5
    vmulpd    %ymm4, %ymm2, %ymm2
    vmulpd    %ymm1, %ymm3, %ymm1
    vsubpd    %ymm5, %ymm0, %ymm0
    vaddpd    %ymm1, %ymm2, %ymm1
    vinsertf128    $1, %xmm0, %ymm0, %ymm2
    vperm2f128    $49, %ymm0, %ymm0, %ymm0
    vinsertf128    $1, %xmm1, %ymm1, %ymm3
    vperm2f128    $49, %ymm1, %ymm1, %ymm1
    vshufpd    $12, %ymm3, %ymm2, %ymm2
    vshufpd    $12, %ymm1, %ymm0, %ymm0
    vmovups    %xmm2, -64(%rcx)
    vextractf128    $0x1, %ymm2, -48(%rcx)
    vmovups    %xmm0, -32(%rcx)
    vextractf128    $0x1, %ymm0, -16(%rcx)
    cmpq    %rsi, %r10
    jb    .L94