多项式乘法(代码)

多项式乘法

模板题, 下面给出一个很方便使用的模板,可当做黑盒代码使用。

#include <set>
#include <cmath>
#include <cctype>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <assert.h>
#include <algorithm>

using namespace std;

#define fir first
#define sec second
#define pb push_back
#define mp make_pair
#define LL long long
#define INF (0x3f3f3f3f)
#define mem(a, b) memset(a, b, sizeof (a))
#define debug(...) fprintf(stderr, __VA_ARGS__)
#define Debug(x) cout << #x << " = " << x << endl
#define travle(i, x) for (register int i = head[x]; i; i = nxt[i])
#define For(i, a, b) for (register int (i) = (a); (i) <= (b); ++ (i))
#define Forr(i, a, b) for (register int (i) = (a); (i) >= (b); -- (i))
#define file(s) freopen(s".in", "r", stdin), freopen(s".out", "w", stdout)
#define ____ debug("go\n")

namespace io {
    static char buf[1<<21], *pos = buf, *end = buf;
    inline char getc()
    { return pos == end && (end = (pos = buf) + fread(buf, 1, 1<<21, stdin), pos == end) ? EOF : *pos ++; }
    inline int rint() {
        register int x = 0, f = 1;register char c;
        while (!isdigit(c = getc())) if (c == '-') f = -1;
        while (x = (x << 1) + (x << 3) + (c ^ 48), isdigit(c = getc()));
        return x * f;
    }
    inline LL rLL() {
        register LL x = 0, f = 1; register char c;
        while (!isdigit(c = getc())) if (c == '-') f = -1;
        while (x = (x << 1ll) + (x << 3ll) + (c ^ 48), isdigit(c = getc()));
        return x * f;
    }
    inline void rstr(char *str) {
        while (isspace(*str = getc()));
        while (!isspace(*++str = getc()))
            if (*str == EOF) break;
        *str = '\0';
    }
    template<typename T> 
        inline bool chkmin(T &x, T y) { return x > y ? (x = y, 1) : 0; }
    template<typename T>
        inline bool chkmax(T &x, T y) { return x < y ? (x = y, 1) : 0; }    
}
using namespace io;


const int N = 4e6 + 10;
const double PI = acos(-1);
struct vec 
{
    double x, y;
    vec(double xx = 0, double yy = 0) { x = xx, y = yy; }
    friend vec operator+ (const vec &a, const vec &b)
    { return vec(a.x + b.x, a.y + b.y); }
    friend vec operator- (const vec &a, const vec &b)
    { return vec(a.x - b.x, a.y - b.y); }
    friend vec operator* (const vec &a, const vec &b)
    { return vec(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x); }
};
vec omega(int n, int k)
{ return vec(cos(2 * PI / n * k), sin(2 * PI / n * k)); }

namespace FFT 
{
    int r[N];
    inline void DFT(vec *f, int type, int lim)
    {
        for (int i = 0; i < lim; ++ i)
            if (i < r[i])
                swap(f[i], f[r[i]]);
        for (int l = 2; l <= lim; l <<= 1) 
        {
            int k = l >> 1;
            for (vec *p = f; p != f + lim; p += l)
                for (int i = 0; i < k; ++ i)
                {
                    vec t = omega(l, type * i) * p[i + k];
                    p[i + k] = p[i] - t;
                    p[i] = p[i] + t;
                }
        }
    }
    inline void Mul(double *a, int n1, double *b, int n2, double *ans) 
    {
        int lim = 1, l = 0;
        while (lim <= n1 + n2) lim <<= 1, l++;
        for (int i = 0; i < lim; ++ i)
            r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
        static vec a1[N], b1[N];
        for (int i = 0; i < lim; ++ i) 
        {
            a1[i] = vec(i <= n1 ? a[i] : 0, 0);
            b1[i] = vec(i <= n2 ? b[i] : 0, 0);
        }
        DFT(a1, 1, lim); DFT(b1, 1, lim);
        for (int i = 0; i <= lim; ++ i)
            a1[i] = a1[i] * b1[i];
        DFT(a1, -1, lim);
        for (int i = 0; i <= n1 + n2; ++ i)
            ans[i] = a1[i].x / lim;
    }
}
int main() {
#ifndef ONLINE_JUDGE
    file("FFT_Template");
#endif
    int n, m;
    n = rint(), m = rint();
    static double a[N], b[N], c[N];
    for (int i = 0; i <= n; ++ i) a[i] = rint();
    for (int i = 0; i <= m; ++ i) b[i] = rint();
    FFT::Mul(a, n, b, m, c);
    for (int i = 0; i <= n + m; ++ i)
        printf("%d ", (int)(c[i] + 0.5));
}
posted @ 2019-03-15 22:30  茶Tea  阅读(325)  评论(0编辑  收藏  举报