虾米-K的其它豆列
&&&&&&&&&&&&
&(1人收藏)雪舞喵和FFT=。= | 雪舞喵窝(^=。=^)
当前位置&:& / /雪舞喵和FFT=。=
于是今天雪舞喵听了一上午的fft，终于把迭代版预处理搞会了=。=然后就来更新一下
那个迭代版预处理其实就是观察一下预处理前后序列的性质=。=比如说原来一个8长度的序列， 下标分别是
0 1 2 3 4 5 6 7
递归预处理之后呢 重新排列原序列的下标 变成了这样
0 4 2 6 1 5 3 7
考虑一下把它们写成二进制形式
000 001 010 011 100 101 110 111 ==&
000 100 010 110 001 101 011 111
不就是把二进制形式下的高位变成低位，低位变成高位么=。=
换句话说 在二进制下倒过来读
这个很简单啊 只要先预处理出0~n-1的数字二进制反转后是多少，然后再对原序列应用这种变换就好了喵。
而且0~n-1的数字二进制反转可以线性预处理。首先0反转之后是0，1反转之后是n&&1
接下来 对于一个二进制数x(比如，)，设它的反转为rev(x)，我们考虑rev(x&&1)(是)，只要把它&&1在把x的最低位填到它的最高位就可以啦
最后附上新(O(n))的狗带函数代码
void godie(int *v, int l)//预处理
int au[l], ag[l];
au[0] = 0;
au[1] = l&&1;
memcpy(ag, v, l&&2);
yuki(2, l)
au[i] = (au[i&&1]&&1) | ((i&1)*(l&&1));
yuki(0, l)
v[i] = ag[au[i]];
1234567891011
void godie(int *v, int l)//预处理{&&int au[l], ag[l];&&au[0] = 0;&&au[1] = l&&1;&&memcpy(ag, v, l&&2);&&yuki(2, l)&&&&au[i] = (au[i&&1]&&1) | ((i&1)*(l&&1));&&yuki(0, l)&&&&v[i] = ag[au[i]];}
顺便一提，原来那个递归狗带在codevs上交超大整数乘法是合计543ms，最慢175ms。改成这个之后变成合计337ms，最慢106ms了。
喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵(分割线)喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵喵
于是经过昨天一天的艰苦奋斗，今天总算是把fft这个噩梦解决了=。=
以下是雪舞喵的个人理解，如果有不清楚的地方，原因一定是你是神犇，我是苣蒻，我们思维不是一个级别的。
所以呢，为什么要学FFT呢，因为我们要用它计算DFT，为什么又要计算DFT呢，因为我们要算两个序列的卷积。然而暴力算卷积是 \(O(n^2)\) 的……喵？卷积是什么？
嗯……大概就是这样，给定两个长度为n的序列A和B， 则它们的卷积为一个长度为n的序列C，并满足以下关系：
\[C_{k}=\sum_{i=0}^{k}A_{i}B_{k-i}\]
那么假如现在有两个n项的多项式 \(A(x)=\sum_{i=0}^{n}a_ix^i, B(x)=\sum_{i=0}^{n}b_ix^i\)
我们得到它们的系数序列 \(A=[a_1, a_2, ..., a_n], B=[b_1, b_2, ..., b_n]\) ，然后手动推一下(或者脑补一下)两个多项式相乘得到的系数序列C，可以发现C序列正好是A，B序列的卷积(只要在A,B序列后各补上n个0)。
因此，多项式乘法就是多项式系数序列的卷积。
然而暴力算卷积是 \(O(n^2)\) 的……为了优化，我们需要知道一个定理，叫做卷积定理。
这个定理的内容是这样的，设F(A)为A序列进行DFT后的结果，C序列是A和B序列的卷积，那么通过定理我们可以推出F(C)和F(A),F(B)的关系 \(F(C)_k = F(A)_kF(B)_k\) ，这个关系就简单多了啊，这样，我们如果知道了F(A)和F(B)，就能在O(n)的时间内就能算出来F(C)。
至于证明雪舞喵一向非常头疼，如果需要请自行百度=。=
于是下一步来了，怎么计算A和B的DFT呢？
根据定义， \(F(A)_{k} = \sum_{i=0}^{n-1}\omega _{n}^{ik}x_i\)
这个 \(\omega _n\) 是一个叫做n次单位根的东西，他等于 \(sin(2\pi /n)+cos(2\pi /n)i\) ，定义就是它的n次幂等于1。
这是个复数啊，因此FFT中的所有计算(+,-,*)都是复数计算。
**对于不知道什么是复数的，雪舞喵也科普一下吧
**首先呢，在实数域内，对一个负数开平方根是无意义的，这个就比较麻烦，总是需要乱七八糟讨论。
**数学家们和雪舞喵一样非常懒，觉的这讨论太麻烦了，就让他有意义吧，然后就发明了复数。
**一个复数由两个实数组成，称为实部和虚部，假设它的实部为a，虚部为b，那么它形式上表示为a+bi
**这个i是什么东西呢，可以理解为一个常量，定义就是 \(i^2 = -1\)
**因此他就类似一个二项式，那么两个复数(a1, b1)+(a2, b2)(用一种比较方便的表示方法)就等于(a1+a2, b1+b2)，减法同理
**乘法(a1, b1)*(a2, b2) = (a1*a2-b1*b2, a1*b2+a2*b1)
然而，这公式喵了个咪的不是 \(O(n^2)\) 的么！
不要急，接下来才是真正的FFT。
刚才说过FFT是个算法，他是用来计算DFT的，并且可以在O(nlogn)的时间内算出一个序列的DFT
不过他有两个小缺点。
小缺点待会再说，先研究一下他是怎么回事=。=
对于一个序列A，我们把他拆偶数下标元素组成的偶子序列E和奇数下标元素组成的奇子序列O，那么如果我们算出F(E)和F(O)，就可以在线性时间内将两个DFT合并成为F(A)。这个证明起来比较容易
假定n为偶数，则 \(F(A)_{k} = \sum_{i = 0}^{n-1} \omega _{n}^{ik}x_i = \sum_{i = 2t}^{n-2}\omega _{n}^{ik}x_i+\sum_{i = 2t+1}^{n-1}\omega _n^{ik}x_i\)
\[=\sum_{t = 0}^{n/2}\omega _{n}^{2tk}x_{2t}+\sum_{t=0}^{n/2}\omega _n^{2tk+k}x_{2t+1}=\sum_{t=0}^{n-1}\omega _{n/2}^{tk}x_{2t}+\omega _n^k\sum_{t=0}^{n-1}\omega_{n/2}^{tk}x_{2t+1}\]
\[=F(E)_k+\omega _n^kF(O)_k\]
以上证明用到了单位根的一些奇怪的性质，比如说 \(\omega _n^{2t} = \omega _{n/2}^t\) 之类的，它还有一条性质 \(\omega _n^k = -\omega _n^{n/2+k}\) 。
这性质有啥用呢，用来计算 \(F(A)_{k+n/2}\) 。求出的两个子序列长度只有原序列的一半，因此如果一次计算序列的一个元素的话需要扫描两个子序列两遍，这样就不能在一个大数组中搞了，本着节省空间的原则，我们同时计算出 \(F(A)_k\) 和 \(F(A)_{k+n/2}\) 就可以覆盖 \(F(E)_k\) 和 \(F(O)_k\) ，就可以不断的在一个数组里来回搞了。
由于要求n为偶数，而且我们需要递归解决问题，因此我们必须先把序列补到2的幂，这是第一个小缺点。
怎么补呢，在低位方向不断补0就可以了。(低位在后的话就往后插，低位在前就往前插)
这样我们就可以使用分治法了，概括一下FFT的流程：
1.划分问题：把给定的长度为2的幂的序列划分为偶子序列和奇子序列两部分。
2.递归求解：递归的计算偶子序列和奇子序列的DFT。
3.合并问题：把两个子序列线性合并。
显然FFT的复杂度是O(nlogn)
边界就是一个元素的序列DFT后还是它本身。
然而在实现上还有一个细节需要注意，就是我们需要对这个序列进行预处理。
为什么呢，因为我们如果正常的计算这个序列，比如说计算一个长度为8的序列，在算第0位时可以同时算出第4位，然而 \(F(E)_0\) 保存在第0位， \(F(O)_0\) 却保存在第1位，不能正好覆盖原来的计算结果。为了避免这个问题，我们操作一下这个序列，只要把偶子序列放在序列的左半边，奇子序列放在序列的右半边，在递归的预处理左半边序列和右半边序列。这样就可以避免这个问题啦。
同时这样处理之后我们还可以把FFT的主过程写成迭代版的。
虽然市面上流行版本的FFT代码这个预处理也是迭代版的，用了一些奇怪的位运算，把空间复杂度彻底压到O(n)了。
然而我并不会奇怪的位运算，就写了一个递归版预处理。(这样的话空间复杂度其实也是O(n)，因为n+n/2+n/4+n/8+......=2n，不过new的时间常数稍大)(见下面代码)
接下来研究一下实现多项式乘法。
首先得到系数序列。如果给定的两个序列不等长，需要在较短的序列高位方向补0补到一样长。
接下来向后补到序列长度的二倍。
再向后补到2的幂。(这个后事实上是看你怎么保存这个序列，如果高位在前就向低位方向，如果低位在前就向高位方向，这两部补位只要保证算卷积的时候如果按顺序计算，必须先计算的是保存数的部分，在计算保存0的部分就可以啦)
然后计算两个序列的DFT。
根据卷积定理，将DFT对应项相乘。
计算结果序列的IDFT(即DFT的逆运算，定义是 \(F^{-1}(A)_k = \frac{1}{n}\sum_{i=0}^{n-1}\omega _{n}^{-ik}x_i\) ，实现上只需要把 \(2\pi\) 改成 \(-2\pi\) 就可以啦。
所以上代码了(我以大数乘法为例)
#define yuki(x, y) for(int i = i & (y); ++i)
#define yukj(x, y) for(int j = j & (y); ++j)
#define yukii(x, y) for(int i = i &= (y); ++i)
#define yukji(x, y) for(int j = j &= (y); ++j)
#define yuk(x, y, z) for(int x = x & (z); ++x)
#define yui(x, y, z) for(int x = x &= (z); --x)
#define sclr(x) memset(x, 0, sizeof(x))
#define sclr1(x) memset(x, -1, sizeof(x))
#define scl(x, y) memset(x, y, sizeof(x))
#define ft first
#define sc second
#define mp make_pair
#include &cstdio&
#include &cstring&
#include &algorithm&
#include &utility&
#include &cmath&
#include &vector&
#include &stack&
char sa[101000], sb[101000];
int la, lb,
stack&char&
const double _2pi = 2*3.;
typedef pair&double, double&
cpx operator+(const cpx &lhs, const cpx &rhs)//复数运算
return mp(lhs.ft+rhs.ft, lhs.sc+rhs.sc);
cpx operator-(const cpx &lhs, const cpx &rhs)
return mp(lhs.ft-rhs.ft, lhs.sc-rhs.sc);
cpx operator*(const cpx &lhs, const cpx &rhs)
return mp(lhs.ft*rhs.ft-lhs.sc*rhs.sc, lhs.ft*rhs.sc+lhs.sc*rhs.ft);
cpx va[1048600], vb[1048600], vc[1048600];
void sread()//输入
scanf("%s%s", sa, sb);
la = strlen(sa);
lb = strlen(sb);
if(la & lb)
swap(la, lb);
swap(sa, sb);
int t = la-
yuki(0, la) va[i] = mp((double)(sa[i]-'0'), 0.0);
yuki(t, la) vb[i] = mp((double)(sb[i-t]-'0'), 0.0);//在较短的序列高位方向补0
for(l = 1; l & l &&= 1);
void godie(cpx *v, int l)//递归预处理序列
if(l &= 2)
int m = l&&1;
cpx *aug[2];
aug[0] = new cpx[m+10];
aug[1] = new cpx[m+10];
yuki(0, l) aug[i&1][i&&1] = v[i];//把偶子序列放在aug[0]，奇子序列放在aug[1]
godie(aug[0], m);
godie(aug[1], m);
memcpy(v, aug[0], m*sizeof(cpx));
memcpy(v+m, aug[1], m*sizeof(cpx));
delete[] aug[0];
delete[] aug[1];
void fft(cpx *v, int l, bool flg)//flg==true表示DFT,flg==false表示IDFT
const double pi2 = flg ? _2pi : -_2
godie(v, l);//预处理
for(int k = 2; k &= k &&= 1)//k表示待合并序列长度
cpx wn = mp(cos(pi2/k), sin(pi2/k));
for(int i = 0; i & i += k)//合并每个待合并序列
cpx w = mp(1.0, 0.0);
yukj(0, k&&1)
cpx t1 = v[i+j], t2 = v[i+j+(k&&1)]*w;
v[i+j] = t1+t2;
v[i+j+(k&&1)] = t1-t2;
cpx inv = mp(1.0/l, 0.0);//如果计算IDFT，需要把序列乘1/n
yuki(0, l) v[i] = v[i]*
void printans()
int cura = 0;
yui(i, la-1, 0)
cura += round(vc[i].ft);
int t = cura%10;
cura /= 10;
if(i != la-1)ans.push(t+'0');
while(cura)
ans.push(cura%10+'0');
cura /= 10;
while(!ans.empty())
if(!ok && ans.top() == '0') ans.pop();
putchar(ans.top());
ans.pop();
int main(int argc, char **argv)
fft(va, l, true);
fft(vb, l, true);//计算A和B的DFT
yuki(0, l) vc[i] = va[i]*vb[i];//对应项相乘
fft(vc, l, false);//计算答案的IDFT
printans();
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133
#define yuki(x, y) for(int i = i & (y); ++i)#define yukj(x, y) for(int j = j & (y); ++j)#define yukii(x, y) for(int i = i &= (y); ++i)#define yukji(x, y) for(int j = j &= (y); ++j)#define yuk(x, y, z) for(int x = x & (z); ++x)#define yui(x, y, z) for(int x = x &= (z); --x)#define sclr(x) memset(x, 0, sizeof(x))#define sclr1(x) memset(x, -1, sizeof(x))#define scl(x, y) memset(x, y, sizeof(x))#define ft first#define sc second#define mp make_pair#include &cstdio&#include &cstring&#include &algorithm&#include &utility&#include &cmath&#include &vector&#include &stack&using namespace std;typedef long long lol;char sa[101000], sb[101000];int la, lb, l;stack&char& ans;const double _2pi = 2*3.;typedef pair&double, double& cpx;cpx operator+(const cpx &lhs, const cpx &rhs)//复数运算{&&return mp(lhs.ft+rhs.ft, lhs.sc+rhs.sc);}cpx operator-(const cpx &lhs, const cpx &rhs){&&return mp(lhs.ft-rhs.ft, lhs.sc-rhs.sc);}cpx operator*(const cpx &lhs, const cpx &rhs){&&return mp(lhs.ft*rhs.ft-lhs.sc*rhs.sc, lhs.ft*rhs.sc+lhs.sc*rhs.ft);}cpx va[1048600], vb[1048600], vc[1048600];void sread()//输入{
scanf("%s%s", sa, sb);
la = strlen(sa);
lb = strlen(sb);
if(la & lb)
swap(la, lb);
swap(sa, sb);
int t = la-lb;
yuki(0, la) va[i] = mp((double)(sa[i]-'0'), 0.0);
yuki(t, la) vb[i] = mp((double)(sb[i-t]-'0'), 0.0);//在较短的序列高位方向补0
for(l = 1; l & la; l &&= 1);}void godie(cpx *v, int l)//递归预处理序列{&&if(l &= 2) return;&&int m = l&&1;&&cpx *aug[2];&&aug[0] = new cpx[m+10];&&aug[1] = new cpx[m+10];&&yuki(0, l) aug[i&1][i&&1] = v[i];//把偶子序列放在aug[0]，奇子序列放在aug[1]&&godie(aug[0], m);&&godie(aug[1], m);&&memcpy(v, aug[0], m*sizeof(cpx));&&memcpy(v+m, aug[1], m*sizeof(cpx));&&delete[] aug[0];&&delete[] aug[1];}void fft(cpx *v, int l, bool flg)//flg==true表示DFT,flg==false表示IDFT{&&const double pi2 = flg ? _2pi : -_2pi;&&godie(v, l);//预处理&&for(int k = 2; k &= l; k &&= 1)//k表示待合并序列长度&&{&&&&cpx wn = mp(cos(pi2/k), sin(pi2/k));&&&&for(int i = 0; i & l; i += k)//合并每个待合并序列&&&&{&&&&&&cpx w = mp(1.0, 0.0);&&&&&&yukj(0, k&&1)&&&&&&{&&&&cpx t1 = v[i+j], t2 = v[i+j+(k&&1)]*w;&&&&v[i+j] = t1+t2;&&&&v[i+j+(k&&1)] = t1-t2;&&&&w = w * wn;&&&&&&}&&&&}&&}&&if(!flg)&&{&&&&cpx inv = mp(1.0/l, 0.0);//如果计算IDFT，需要把序列乘1/n&&&&yuki(0, l) v[i] = v[i]*inv;&&}}void printans(){&&int cura = 0;&&yui(i, la-1, 0)&&{&&&&cura += round(vc[i].ft);&&&&int t = cura%10;&&&&cura /= 10;&&&&if(i != la-1)ans.push(t+'0');&&}&&while(cura)&&{&&&&ans.push(cura%10+'0');&&&&cura /= 10;&&}&&bool ok = false;&&while(!ans.empty())&&{&&&&if(!ok && ans.top() == '0') ans.pop();&&&&else&&&&{&&&&&&putchar(ans.top());&&&&&&ans.pop();&&&&&&ok = true;&&&&}&&}&&puts("");}int main(int argc, char **argv){&&sread();&&fft(va, l, true);&&fft(vb, l, true);//计算A和B的DFT&&yuki(0, l) vc[i] = va[i]*vb[i];//对应项相乘&&fft(vc, l, false);//计算答案的IDFT&&printans();&&return 0;}
从代码中可以看出来，由于使用了double类型，精度误差是不小的，这是第二个小缺点。
然而这第二个小缺点是可以避免的，那就是使用NTT。
NTT是使用整数(而不是会掉精度的复数)在模一个素数的意义下进行的(因此为了保证结果正确，这个素数要尽量大)。
其实实现上和FFT几乎没区别，就是把复数运算换成带模数的整数运算，单位根换成一个整数而已。
这样的话我们需要选择一个数 \(a_n\) 来代替原来使用的 \(\omega _n\) ，选择什么样的数呢，只要满足单位根的3条性质就可以啦。
刚才已经说过两条了，其实第3条也说过，就是它的定义，需要保证 \(a_n^n=1, a_n^p \ne 1 (p \ne kn, k \in \mathbb{Z})\)
然后经过一系列复杂的推理证明 \(a_n=g^{\frac{M-1}{n}}\) 是可行的(其中M是那个素数，g是M的原根。注意NTT中所有运算全是在模M的意义下进行的，以下部分就省略了)。
**关于原根：原根的定义是满足方程 \(g^{x}\equiv 1(mod\,M)\) 的最小整数解为 \(x=\varphi(M)\) 的所有g。对于一个素数就是x=n-1
**由于一个数n他的最小原根是loglogloglogloglogn级别的(6个log)，所以可以直接从2枚举
**判断他是不是原根怎么搞呢，首先我们对M-1进行质因数分解，得到它的所有不同的质因数p1, p2, p3……pk
**只要判断出 \(g^{\frac{M-1}{p_i}}\not\equiv 1(mod\,M)\) 对于 \(i \in [1, k]\) 均成立即可。证明较复杂请自行百度。
**事实上我们可以直接暴力判断1到n-1
**因为这个素数肯定是写代码的时候就知道是多少(要么是给定的，要么是自己选择的)
然而想让这个 \(a_n\) 可行我们对M还有一些要求。
首先M是个素数。
其次M-1必须能除开n，即 \(M=k2^p+1(k, p \in \mathbb{Z})\) 。
然后M还要尽量大(但是M的平方还不能爆long long)
因此如果你看一道题它的模数特别的奇怪(比如之类的)(他等于119*2^23+1)而不是常规的1e9+7,1e8+7啥的，就大概需要FFT(然而知道这个并没有什么卵用，还是不会)
好了这样终于不用担心掉精度啦，只要把刚才的代码稍加修改即可～
其他的只需要把INTT中的1/n换成n的逆元，g换成g的逆元即可
#define yuki(x, y) for(int i = i & (y); ++i)
#define yukj(x, y) for(int j = j & (y); ++j)
#define yukii(x, y) for(int i = i &= (y); ++i)
#define yukji(x, y) for(int j = j &= (y); ++j)
#define yuk(x, y, z) for(int x = x & (z); ++x)
#define yui(x, y, z) for(int x = x &= (z); --x)
#define sclr(x) memset(x, 0, sizeof(x))
#define sclr1(x) memset(x, -1, sizeof(x))
#define scl(x, y) memset(x, y, sizeof(x))
#define ft first
#define sc second
#define mp make_pair
#include &cstdio&
#include &cstring&
#include &algorithm&
#include &utility&
#include &cmath&
#include &vector&
#include &stack&
char sa[101000], sb[101000];
int la, lb,
stack&char&
const int M = , g = 5;//M：模数，g：原根
int splus(int a, int b){return (a+b)%M;}//模数意义下的运算
int &eplus(int &a, int b){return a = splus(a, b);}
int stimes(int a, int b){return (lol)a*b%M;}
int &etimes(int &a, int b){return a = stimes(a, b);}
int sminus(int a, int b)
int t = (a-b)%M;
if(t & 0) t += M;
int &eminus(int &a, int b){return a = sminus(a, b);}
int ksm(int a, int b)//快速幂
if(b == 1)
int t = ksm(a, b&&1);
etimes(t, t);
if(b&1) etimes(t, a);
int va[1048600], vb[1048600], vc[1048600];
void sread()//读入
scanf("%s%s", sa, sb);
la = strlen(sa);
lb = strlen(sb);
if(la & lb)
swap(la, lb);
swap(sa, sb);
int t = la-
yuki(0, la) va[i] = sa[i]-'0';
yuki(t, la) vb[i] = sb[i-t]-'0';
for(l = 1; l & l &&= 1);
void godie(int *v, int l)//预处理
if(l &= 2)
int m = l&&1;
int *aug[2];
aug[0] = new int[m+10];
aug[1] = new int[m+10];
yuki(0, l) aug[i&1][i&&1] = v[i];
godie(aug[0], m);
godie(aug[1], m);
memcpy(v, aug[0], m&&2);
memcpy(v+m, aug[1], m&&2);
delete[] aug[0];
delete[] aug[1];
void ntt(int *v, int l, bool flg)//flg==true表示NTT,flg==false表示INTT
const int g2 = flg ? g : ksm(g, M-2);
godie(v, l);
for(int k = 2; k &= k &&= 1)
int an = ksm(g2, (M-1)/k);
for(int i = 0; i & i += k)
int a = 1;
yukj(0, k&&1)
int t1 = v[i+j], t2 = stimes(v[i+j+(k&&1)],a);
v[i+j] = splus(t1, t2);
v[i+j+(k&&1)] = sminus(t1, t2);
etimes(a, an);
int inv = ksm(l, M-2);
yuki(0, l) etimes(v[i], inv);
void printans()
int cura = 0;
yui(i, la-1, 0)
cura += vc[i];
int t = cura%10;
cura /= 10;
if(i != la-1)ans.push(t+'0');
while(cura)
ans.push(cura%10+'0');
cura /= 10;
while(!ans.empty())
if(!ok && ans.top() == '0') ans.pop();
putchar(ans.top());
ans.pop();
int main(int argc, char **argv)
ntt(va, l, true);
ntt(vb, l, true);//这些注释就省略了吧，见上面FFT代码
yuki(0, l) vc[i] = stimes(va[i], vb[i]);
ntt(vc, l, false);
printans();
123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139
#define yuki(x, y) for(int i = i & (y); ++i)#define yukj(x, y) for(int j = j & (y); ++j)#define yukii(x, y) for(int i = i &= (y); ++i)#define yukji(x, y) for(int j = j &= (y); ++j)#define yuk(x, y, z) for(int x = x & (z); ++x)#define yui(x, y, z) for(int x = x &= (z); --x)#define sclr(x) memset(x, 0, sizeof(x))#define sclr1(x) memset(x, -1, sizeof(x))#define scl(x, y) memset(x, y, sizeof(x))#define ft first#define sc second#define mp make_pair#include &cstdio&#include &cstring&#include &algorithm&#include &utility&#include &cmath&#include &vector&#include &stack&using namespace std;typedef long long lol;char sa[101000], sb[101000];int la, lb, l;stack&char& ans;const int M = , g = 5;//M：模数，g：原根int splus(int a, int b){return (a+b)%M;}//模数意义下的运算int &eplus(int &a, int b){return a = splus(a, b);}int stimes(int a, int b){return (lol)a*b%M;}int &etimes(int &a, int b){return a = stimes(a, b);}int sminus(int a, int b){&&int t = (a-b)%M;&&if(t & 0) t += M;&&return t;}int &eminus(int &a, int b){return a = sminus(a, b);}int ksm(int a, int b)//快速幂{&&if(b == 1) return a;&&int t = ksm(a, b&&1);&&etimes(t, t);&&if(b&1) etimes(t, a);&&return t;}int va[1048600], vb[1048600], vc[1048600];void sread()//读入{&&scanf("%s%s", sa, sb);&&la = strlen(sa);&&lb = strlen(sb);&&if(la & lb)&&{&&&&swap(la, lb);&&&&swap(sa, sb);&&}&&int t = la-lb;&&yuki(0, la) va[i] = sa[i]-'0';&&yuki(t, la) vb[i] = sb[i-t]-'0';&&la &&= 1;&&for(l = 1; l & la; l &&= 1);}void godie(int *v, int l)//预处理{&&if(l &= 2) return;&&int m = l&&1;&&int *aug[2];&&aug[0] = new int[m+10];&&aug[1] = new int[m+10];&&yuki(0, l) aug[i&1][i&&1] = v[i];&&godie(aug[0], m);&&godie(aug[1], m);&&memcpy(v, aug[0], m&&2);&&memcpy(v+m, aug[1], m&&2);&&delete[] aug[0];&&delete[] aug[1];}void ntt(int *v, int l, bool flg)//flg==true表示NTT,flg==false表示INTT{&&const int g2 = flg ? g : ksm(g, M-2);&&godie(v, l);&&for(int k = 2; k &= l; k &&= 1)&&{&&&&int an = ksm(g2, (M-1)/k);&&&&for(int i = 0; i & l; i += k)&&&&{&&&&&&int a = 1;&&&&&&yukj(0, k&&1)&&&&&&{&&&&int t1 = v[i+j], t2 = stimes(v[i+j+(k&&1)],a);&&&&v[i+j] = splus(t1, t2);&&&&v[i+j+(k&&1)] = sminus(t1, t2);&&&&etimes(a, an);&&&&&&}&&&&}&&}&&if(!flg)&&{&&&&int inv = ksm(l, M-2);&&&&yuki(0, l) etimes(v[i], inv);&&}}void printans(){&&int cura = 0;&&yui(i, la-1, 0)&&{&&&&cura += vc[i];&&&&int t = cura%10;&&&&cura /= 10;&&&&if(i != la-1)ans.push(t+'0');&&}&&while(cura)&&{&&&&ans.push(cura%10+'0');&&&&cura /= 10;&&}&&bool ok = false;&&while(!ans.empty())&&{&&&&if(!ok && ans.top() == '0') ans.pop();&&&&else&&&&{&&&&&&putchar(ans.top());&&&&&&ans.pop();&&&&&&ok = true;&&&&}&&}&&puts("");}int main(int argc, char **argv){&&sread();&&ntt(va, l, true);&&ntt(vb, l, true);//这些注释就省略了吧，见上面FFT代码&&yuki(0, l) vc[i] = stimes(va[i], vb[i]);&&ntt(vc, l, false);&&printans();&&return 0;}
好了，到这就差不多了吧，然而依然没有什么卵用，原因是FFT的题依然不会喵喵喵～
Snow halationby μ's
Your browser does not support the audio element.
最新(^~。~^)
最热(^=。=^)
- 179 views - 166 views - 107 views - 99 views - 93 views
随便看看喵
按月查看喵
按月查看喵
2017年六月 &(1)
2017年四月 &(1)
2017年二月 &(1)
2016年十二月 &(4)
2016年十月 &(1)
2016年九月 &(1)
2016年七月 &(1)
2016年五月 &(8)
2016年四月 &(4)
2016年三月 &(7)
分类查看喵分类查看喵
选择分类目录
代码区&&(5)
公告栏&&(2)
日志区&&(18)
水区&&(18)
我是日历喵
2017年十月
891011121314
15161718192021
22232425262728更多频道内容在这里查看
爱奇艺用户将能永久保存播放记录
过滤短视频
暂无长视频（电视剧、纪录片、动漫、综艺、电影）播放记录，
小鱼儿滴喵麻麻(=^_^=)
Ta还未创建任何播单～
正在加载中，请稍候...
没有更多内容了~
Copyright (C) 2017
All Rights Reserved}