复根的设置:typedef complex<double> E;
π=acos(-1):
#define pi acos(-1)
const double pi=3.14159265358979;
单位复根:
E wn(cos(pi/i),f*sin(pi/i)); (快)
E wn=exp(E(0,f*pi/i));(慢)exp(x)为e的x次方
//ω单位复根
FFT
#include<bits/stdc++.h>
#define N 2621450
#define pi acos(-1)
using namespace std;
typedef complex<double> E;
int n,m,l,r[N];
E a[N],b[N];
void fft(E *a,int f){
for(int i=0;i<n;i++)if(i<r[i])swap(a[i],a[r[i]]);
for(int i=1;i<n;i<<=1){
E wn(cos(pi/i),f*sin(pi/i));
for(int p=i<<1,j=0;j<n;j+=p){
E w(1,0);
for(int k=0;k<i;k++,w*=wn){
E x=a[j+k],y=w*a[j+k+i];
a[j+k]=x+y;a[j+k+i]=x-y;
}
}
}
}
inline int read(){
int f=1,x=0;char ch;
do{ch=getchar();if(ch=='-')f=-1;}while(ch<'0'||ch>'9');
do{x=x*10+ch-'0';ch=getchar();}while(ch>='0'&&ch<='9');
return f*x;
}
int main(){
n=read();m=read();
for(int i=0;i<=n;i++)a[i]=read();
for(int i=0;i<=m;i++)b[i]=read();
m+=n;for(n=1;n<=m;n<<=1)l++;
for(int i=0;i<n;i++)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
fft(a,1);fft(b,1);
for(int i=0;i<=n;i++)a[i]=a[i]*b[i];
fft(a,-1);
for(int i=0;i<=m;i++)printf("%d ",(int)(a[i].real()/n+0.5));
}
NFT(快速数论变换)
1
#include<cstdio>
#define getchar() (p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 1<<21, stdin), p1 == p2) ? EOF : *p1++)
#define swap(x,y) x ^= y, y ^= x, x ^= y
#define LL long long
const int MAXN = 3 * 1e6 + 10, P = 998244353, G = 3, Gi = 332748118;
char buf[1<<21], *p1 = buf, *p2 = buf;
inline int read() {
char c = getchar(); int x = 0, f = 1;
while(c < '0' || c > '9') {if(c == '-') f = -1; c = getchar();}
while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar();
return x * f;
}
int N, M, limit = 1, L, r[MAXN];
LL a[MAXN], b[MAXN];
inline LL fastpow(LL a, LL k) {
LL base = 1;
while(k) {
if(k & 1) base = (base * a ) % P;
a = (a * a) % P;
k >>= 1;
}
return base % P;
}
inline void NTT(LL *A, int type) {
for(int i = 0; i < limit; i++)
if(i < r[i]) swap(A[i], A[r[i]]);
for(int mid = 1; mid < limit; mid <<= 1) {
LL Wn = fastpow( type == 1 ? G : Gi , (P - 1) / (mid << 1));
for(int j = 0; j < limit; j += (mid << 1)) {
LL w = 1;
for(int k = 0; k < mid; k++, w = (w * Wn) % P) {
int x = A[j + k], y = w * A[j + k + mid] % P;
A[j + k] = (x + y) % P,
A[j + k + mid] = (x - y + P) % P;
}
}
}
}
int main() {
N = read(); M = read();
for(int i = 0; i <= N; i++) a[i] = (read() + P) % P;
for(int i = 0; i <= M; i++) b[i] = (read() + P) % P;
while(limit <= N + M) limit <<= 1, L++;
for(int i = 0; i < limit; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (L - 1));
NTT(a, 1);NTT(b, 1);
for(int i = 0; i < limit; i++) a[i] = (a[i] * b[i]) % P;
NTT(a, -1);
LL inv = fastpow(limit, P - 2);
for(int i = 0; i <= N + M; i++)
printf("%d ", (a[i] * inv) % P);
return 0;
}
2
#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cmath>
//#include<complex>
using namespace std;
//typedef complex<double> cd;
typedef long long LL;
void exgcd(int a,int b,int& x,int& y){
if(b==0){
x=1;
y=0;
return;
}
int x0,y0;
exgcd(b,a%b,x0,y0);
x=y0;y=x0-int(a/b)*y0;
}
int Inv(int a,int p){
int x,y;
exgcd(a,p,x,y);
x%=p;
while(x<0)x+=p;
return x;
}
int qpow(int a,int b,int p){
if(b<0){
b=-b;
a=Inv(a,p);
}
LL ans=1,mul=a%p;
while(b){
if(b&1)ans=ans*mul%p;
mul=mul*mul%p;
b>>=1;
}
return ans;
}
#define maxn (65537*2)
const int MOD=479*(1<<21)+1,G=3;
int rev[maxn];
void get_rev(int bit){
for(int i=0;i<(1<<bit);i++){
rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
}
}
//from internet
//for(int i=0; i<NUM; i++)
//{
// int t = 1 << i;
// wn[i] = quick_mod(G, (P - 1) / t, P);
//}
LL a[maxn],b[maxn];
void ntt(LL* a,int n,int dft){
for(int i=0;i<n;i++){
if(i<rev[i])
swap(a[i],a[rev[i]]);
}
for(int step=1;step<n;step<<=1){
LL wn;
wn=qpow(G,dft*(MOD-1)/(step*2),MOD);
for(int j=0;j<n;j+=step<<1){
LL wnk=1;//这里一定要用long long不然会迷之溢出
for(int k=j;k<j+step;k++){
LL x=a[k]%MOD,y=(wnk*a[k+step])%MOD;//这里也要用long long
a[k]=(x+y)%MOD;
a[k+step]=((x-y)%MOD+MOD)%MOD;
wnk=(wnk*wn)%MOD;
}
}
}
if(dft==-1){
int nI=Inv(n,MOD);
for(int i=0;i<n;i++)
a[i]=a[i]*nI%MOD;
}
}
#include<cstring>
char s1[maxn],s2[maxn];
int main(){
//scanf("%*d");
scanf("%s%s",s1,s2);
int l1=strlen(s1),l2=strlen(s2);
for(int i=0;i<l1;i++)a[i]=s1[l1-i-1]-'0';
for(int i=0;i<l2;i++)b[i]=s2[l2-i-1]-'0';
int bit,s=2;
for(bit=1;(1<<bit)<(l1+l2-1);bit++)s<<=1;
get_rev(bit);ntt(a,s,1);ntt(b,s,1);
for(int i=0;i<s;i++)
a[i]=a[i]*b[i]%MOD;
ntt(a,s,-1);
for(int i=0;i<s;i++){
a[i+1]+=a[i]/10;
a[i]%=10;
}
int cnt=s;
while(cnt>=0 && a[cnt]==0)cnt--;
if(cnt==-1)printf("0");
for(int i=cnt;i>=0;i--){
printf("%d",a[i]);
}
putchar('\n');
return 0;
}