#include <iostream>
#include <vector>
#include <complex>
#include <cmath>

using namespace std;

vector<complex<double>> fft(const vector<complex<double>> &a, bool inverse)
{
    const int n = a.size();

    if (n == 1)
        return a;

    complex<double> w = 1.0;
    complex<double> w_n;
    if (inverse) {
        w_n = exp(1i * (-2.0) * M_PI / (double) n);
    } else {
        w_n = exp(1i * 2.0 * M_PI / (double) n);
    }

    // Delimo
    vector<complex<double>> a_0(n / 2); 
    vector<complex<double>> a_1(n / 2);
    for (int i = 0; i < n / 2; i++) {
        a_0[i] = a[2*i];
        a_1[i] = a[2*i + 1];
    }

    // Rekurzija
    vector<complex<double>> y_0 = fft(a_0, inverse);
    vector<complex<double>> y_1 = fft(a_1, inverse);

    // Objedinjujemo
    vector<complex<double>> y (n);
    for (int k = 0; k < n / 2; k++) {
        y[k] = y_0[k] + w * y_1[k];
        y[k + n / 2] = y_0[k] - w * y_1[k];
        w *= w_n;
    }

    return y;
}

vector<complex<double>> convolution(const vector<complex<double>> &a,
                                    const vector<complex<double>> &b)
{
    const int n = a.size();

    vector<complex<double>> x = fft(a, false);
    vector<complex<double>> y = fft(b, false);

    for (int i = 0; i < n; i++) {
        x[i] *= y[i];
    }

    vector<complex<double>> c = fft(x, true);

    for (int i = 0; i < n; i++) {
        c[i] /= (double) n;
    }

    return c;
}

void print_poly(const vector<complex<double>> &a) 
{
    const int n = a.size();

    for (int i = 0; i < n - 1; i++) {
        cout << real(a[i]) << "*x^" << i << " + ";
    }
    cout << real(a[n - 1]) << "*x^" << n - 1 << endl;
}

int main(void)
{
    vector<complex<double>> a = {1., 2., 1., 1., 0., 0., 0., 0.};
    vector<complex<double>> b = {3., 2., 1., 2., 0., 0., 0., 0.};

    print_poly(a);
    print_poly(b);

    vector<complex<double>> c = convolution(a, b);

    print_poly(c);

    return 0;
}