luogu P4512 多項式除法 (模板題、FFT、多項式求逆)

手動博客搬家: 本文發表于20181206 14:42:53, 原地址https://blog.csdn.net/suncongbo/article/details/84853342

題目鏈接: https://www.luogu.org/problemnew/show/P4512

沒想到這算法這么蠢。。一點都不難啊。。我連這都推不出來我是不是沒救了

這個多項式滿足\(A(x)=B(x)Q(x)+R(x)\), 如果已知\(R(x)\)是\(0\), 那顯然很好處理，求個逆就行了。
那如果有余數呢？很簡單，如果我們把這個多項式的系數翻轉(reverse)，那么\(R(x)\)就從低次項變成了高次項，低次項就不再受\(R(x)\)的而影響了。
這是我們的基本思路。下面我們來形式化這個過程:
\(A(x)=B(x)Q(x)+R(x)\)
對于\(n\)次多項式\(F(x)\)令\(F_R(x)=x^nF(\frac{1}{x})\) (這就是前面所說的reverse操作)
則有\(x^nA(\frac{1}{x})=x^mB(\frac{1}{x})x^{n-m}Q(\frac{1}{x})+x^{m-1}R(\frac{1}{x})x^{n-m+1}\)
\(A_R(x)=B_R(x)Q_R(x)+x^{n-m+1}R_R(x)\)
\(A_R(x)\equiv B_R(x)Q_R(x) (\mod x^{n-m+1})\)
于是我們求\(B_R(x)\)在\(\mod x^{n-m+1}\)意義下的逆，然后乘以\(A_R(x)\)即可求出\(Q_R(x)\), 從而得到\(Q(x)\).
然后用\(R(x)=A(x)-B(x)Q(x)\)即可求出\(R(x)\).
(雖然算法簡單但是要注意的地方還挺多……容易錯。)
時間復雜度\(O(n\log n)\), 我寫的進行了\(24\)倍常數的ntt.
空間復雜度\(O(n)\), 我的實現好像需要開\(8\)倍。

代碼實現

#include<cstdio> #include<cstdlib> #include<cstring> #include<algorithm> #define llong long long #define modinc(x) {if(x>=P) x-=P;} using namespace std;const int N = 1<<19; const int LGN = 19; const llong G = 3ll; const int P = 998244353; llong tmp1[N+3],tmp2[N+3],tmp3[N+3],tmp4[N+3]; llong tmp5[N+3],tmp6[N+3],tmp7[N+3],tmp8[N+3],tmp9[N+3]; llong a[N+3],b[N+3],q[N+3],r[N+3]; int id[N+3]; int n,m;llong quickpow(llong x,llong y) {llong cur = x,ret = 1ll;for(int i=0; y; i++){if(y&(1ll<<i)) {y-=(1ll<<i); ret = ret*cur%P;}cur = cur*cur%P;}return ret; } llong mulinv(llong x) {return quickpow(x,P-2)%P;}void initid(int dgr) {int len = 0; for(int i=0; i<=LGN; i++) if((1<<i)==dgr) {len = i; break;}id[0] = 0;for(int i=1; i<dgr; i++) id[i] = (id[i>>1]>>1)|(i&1)<<(len-1); }int getdgr(int x) {int ret = 1;while(ret<=x) ret<<=1;return ret; }void ntt(int dgr,int coe,llong poly[],llong ret[]) {initid(dgr);for(int i=0; i<dgr; i++) ret[i] = poly[i];for(int i=0; i<dgr; i++) if(i<id[i]) swap(ret[i],ret[id[i]]);for(int i=1; i<=(dgr>>1); i<<=1){llong tmp = quickpow(G,(P-1)/(i<<1));if(coe==-1) tmp = mulinv(tmp);for(int j=0; j<dgr; j+=(i<<1)){llong expn = 1ll;for(int k=0; k<i; k++){llong x = ret[j+k],y = ret[j+i+k]*expn%P;ret[j+k] = x+y; modinc(ret[j+k]);ret[j+i+k] = x-y+P; modinc(ret[j+i+k]);expn = (expn*tmp)%P;}}}if(coe==-1){llong tmp = mulinv(dgr);for(int j=0; j<dgr; j++) ret[j] = ret[j]*tmp%P;} }void polyinv(int dgr,llong poly[],llong ret[]) {for(int i=0; i<dgr; i++) ret[i] = 0ll;ret[0] = mulinv(poly[0]);for(int i=1; i<=(dgr>>1); i<<=1){for(int j=0; j<(i<<2); j++) tmp1[j] = j<i ? ret[j] : 0ll;for(int j=0; j<(i<<2); j++) tmp2[j] = j<(i<<1) ? poly[j] : 0ll;ntt((i<<2),1,tmp1,tmp3); ntt((i<<2),1,tmp2,tmp4);for(int j=0; j<(i<<2); j++) tmp3[j] = tmp3[j]*tmp3[j]%P*tmp4[j]%P;ntt((i<<2),-1,tmp3,tmp4);for(int j=0; j<(i<<1); j++) ret[j] = (tmp1[j]+tmp1[j]-tmp4[j]+P)%P; }for(int j=dgr; j<(dgr<<1); j++) ret[j] = 0ll; }void polyrev(int dgr,llong poly[],llong ret[]) {for(int i=0; i<dgr; i++) ret[i] = poly[dgr-1-i]; }void polydiv(int dgr1,int dgr2,llong poly1[],llong poly2[],llong ret1[],llong ret2[]) {int _dgr1 = getdgr(dgr1),_dgr2 = getdgr(dgr2);polyrev(dgr2,poly2,tmp5); polyrev(dgr1,poly1,tmp9);polyinv(_dgr1,tmp5,tmp6);for(int i=dgr1-dgr2+1; i<(_dgr1<<1); i++) tmp6[i] = 0ll;ntt(_dgr1<<1,1,tmp9,tmp7); ntt(_dgr1<<1,1,tmp6,tmp8);for(int i=0; i<(_dgr1<<1); i++) tmp7[i] = tmp7[i]*tmp8[i]%P;ntt(_dgr1<<1,-1,tmp7,tmp8);for(int i=dgr1-dgr2+1; i<(_dgr1<<1); i++) tmp8[i] = 0ll;polyrev(dgr1-dgr2+1,tmp8,ret1);ntt(_dgr1<<1,1,poly2,tmp7); ntt(_dgr1<<1,1,ret1,tmp8);for(int i=0; i<(_dgr1<<1); i++) tmp7[i] = tmp7[i]*tmp8[i]%P;ntt(_dgr1<<1,-1,tmp7,ret2);for(int i=dgr2; i<(_dgr1<<1); i++) ret2[i] = 0ll;for(int i=0; i<dgr2-1; i++) ret2[i] = (poly1[i]-ret2[i]+P)%P; }int main() {scanf("%d%d",&n,&m); n++; m++;for(int i=0; i<n; i++) scanf("%lld",&a[i]);for(int i=0; i<m; i++) scanf("%lld",&b[i]);int dgr1 = getdgr(n),dgr2 = getdgr(m);polydiv(n,m,a,b,q,r);for(int i=0; i<=n-m; i++) printf("%lld ",q[i]); puts("");for(int i=0; i<m-1; i++) printf("%lld ",r[i]);return 0; } 發表于 2019-01-23 20:23 suncongbo 閱讀(...) 評論(...) 編輯 收藏 刷新評論刷新頁面返回頂部

總結

以上是生活随笔為你收集整理的luogu P4512 多项式除法 (模板题、FFT、多项式求逆)的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。