$1004535809=479\times 2^{21}+1$ 相加不会超过 $\text{INT_MAX}$。
$2281701377=17\times 2^{27}+1$ 相乘不会超过 $\text{LONG_LONG_MAX}$。
/*
* FFT.cc
*
* @date 2018.2.13
* @author Kvar_ispw17
* @tag Fast Fourier Transform
* @brief This is a Fast Fourier Transform->
* Implement by Kvar_ispw17(enkerewpo@gmail.com)
*
* Copyright 2018 Kvar_ispw17 __Hyperx All rights reserved.
*
*/
#include <algorithm>
#include <stdexcept>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#define S
const int N = 1e5 + 5;
// COMMENT THE DEFINATION YOU DON'T NEED BELOW TO RUN THE PROGRAM
#ifdef S
// NUMBER THEORIETIC TRANSFORM
#define NTT_INTEGER_VERSION
#else
// FAST FOURIER TRANSFORM
#define FFT_COMPLEX_VERSION
#endif
#if defined(FFT_COMPLEX_VERSION) && defined(NTT_INTEGER_VERSION)
#undef FFT_COMPLEX_VERSION
#undef NTT_INTEGER_VERSION
int main() {
try { throw std::runtime_error("<<<PATTERN NOT SELECTED>>") } catch (...) {
fprintf(stderr, "[WARNING] : PLEASE SELECT A PROGRAM TO RUN INTO.\n");
} return 1;
}
#endif
#ifdef FFT_COMPLEX_VERSION
typedef long long ll;
const long double PI = acos(-1.0);
const long double EPS = 1E-8;
class complex {
public:
long double re, im;
complex() {}
complex(long double re, long double im) : re(re), im(im) {}
complex operator + (const complex &x) { return complex(re + x.re, im + x.im); }
complex operator - (const complex &x) { return complex(re - x.re, im - x.im); }
complex operator * (const complex &x) { return complex(re * x.re - im * x.im, im * x.re + re * x.im); }
complex operator / (const complex &x) {
return complex((re * x.re + im * x.im) / (x.re * x.re + x.im * x.im),
(im * x.re - re * x.im) / (x.re * x.re + x.im * x.im));
}
};
void print(complex x, char c) {
printf("\t%.10lf", fabs(x.re) < EPS ? 0.0 : (double)x.re);
printf("\t\t(%.10lf)", fabs(x.im) < EPS ? 0.0 : (double)x.im);
printf("%c", c);
}
int n, rev[N];
complex F[N], w[N];
void FFT(complex * F, int n, int offset) {
for (int i = 0; i < n; i++) if (rev[i] > i) std::swap(F[i], F[rev[i]]);
for (int i = 2; i <= n; i <<= 1) {
complex wi(cos(offset * 2 * PI / i), sin(offset * 2 * PI / i));
for (int j = 0; j < n; j += i) {
complex w(1, 0);
for (int k = j, h = i >> 1; k < j + h; k++) {
complex t = w * F[k + h], u = F[k];
F[k] = u + t;
F[k + h] = u - t;
w = w * wi;
}
}
}
if (offset == -1) for (int i = 0; i < n; i++) F[i].re = F[i].re / n;
}
int main() {
try {
printf("FFT MODE RUNNING...\n");
scanf("%d", &n);
for (int i = 0; i < n; i++) {
double x;
scanf("%lf", &x);
F[i].re = x;
}
n = 1 << (int)ceil(log2(n));
for (int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << ((ll)log2(n) - 1));
printf("\nAFTER INPUT:\n");
for (int i = 0; i < n; i++) print(F[i], '\n');
FFT(F, n, 1);
printf("\nAFTER FFT:\n");
for (int i = 0; i < n; i++) print(F[i], '\n');
FFT(F, n, -1);
printf("\nAFTER IFFT:\n");
for (int i = 0; i < n; i++) print(F[i], '\n');
return 0;
} catch (...) { return 1; }
}
#endif
#ifdef NTT_INTEGER_VERSION
typedef long long ll;
ll n, inv_n, F[N], rev[N], q;
const ll MOD = 1004535809; // = 479 * 2 ^ 21 + 1
const ll g = 3;
ll q_pow(ll a, ll b) {
ll ret = 1;
while (b) {
if (b & 1) ret = ret * a % MOD;
a = a * a % MOD;
b >>= 1;
}
return ret;
}
void NTT(ll F[], ll n, int offset) {
for (int i = 0; i < n; i++) if (rev[i] > i) std::swap(F[i], F[rev[i]]);
for (int i = 2; i <= n; i <<= 1) {
ll wi = q_pow(g, offset == 1 ? (MOD - 1) / i : MOD - 1 - (MOD - 1) / i);
for (int j = 0; j < n; j += i) {
ll w = 1;
for (int k = j, h = i >> 1; k < j + h; k++) {
ll t = w * F[k + h], u = F[k];
F[k] = (u + t) % MOD;
F[k + h] = ((u - t) % MOD + MOD) % MOD;
w = w * wi % MOD;
}
}
}
if (offset == -1) for (int i = 0; i < n; i++) F[i] = F[i] * inv_n % MOD;
}
int main() {
try {
printf("NTT MODE RUNNING...\n");
scanf("%lld", &n);
for (int i = 0; i < n; i++) scanf("%lld", &F[i]);
n = 1 << (int)ceil(log2(n));
printf("%lld\n", n);
inv_n = q_pow(n, MOD - 2);
printf("[DEBUG: inverse of n is %lld]\n", inv_n);
for (int i = 0; i < n; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << ((ll)log2(n) - 1));
printf("\nAFTER INPUT:\n");
for (int i = 0; i < n; i++) printf("\t%lld\n", F[i]);
NTT(F, n, 1);
printf("\nAFTER NTT:\n");
for (int i = 0; i < n; i++) printf("\t%lld\n", F[i]);
NTT(F, n, -1);
printf("\nAFTER INTT:\n");
for (int i = 0; i < n; i++) printf("\t%lld\n", F[i]);
return 0;
} catch (...) { return 1; }
}
#endif
/*
[TEST CASE]
-
8
1 2 3 4 5 6 7 8
-
9
1.35131334 3.56336636 -3.64242251 -8.3575486 7.86623785 -2.4457546 1.35332446 -9.94949959 1.23434524
-
*/Code language: PHP (php)