luogu P4238 多項式求逆 (模板題、FFT)

手動博客搬家: 本文發表于20181125 13:21:46, 原地址https://blog.csdn.net/suncongbo/article/details/84485718

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

題意: 給定\(n\)次多項式\(A(x)\), 求\(n\)次多項式\(B(x)\)滿足\(B(x)A(x)\equiv 1(\mod x^n)\)

題解:
DFT，每個數對\(998244353\)求逆元。IDFT回來。
發現，錯了。
為什么呢？
因為要對\(x^n\)取模。
例如，\(1-x\)在模\(x^2\)意義下的逆元是\(1+x\), 但是在實際上逆元是\(1+x+x^2+x^3+...\), 是無窮和式。
所以此路不通。

考慮求解多項式問題的常用方法——分治法。
設已求\(B_0(x)\)滿足\(B_0(x)A(x)\equiv 1(\mod x^n)\), 現要求\(B(x)\)滿足\(B(x)A(x)\equiv 1(\mod x^{2n})\)
顯然有\(B(x)-B_0(x)\equiv 0(\mod x^n)\)
為了湊出\(x^{2n}\)兩邊平方得
\(B^2(x)-2B_0(x)B(x)+B_0^2(x)\equiv 0(\mod x^{2n})\)
如何求\(B\)呢？因為\(A(x)B(x)\equiv 0(\mod x^{2n})\), 我們將式子兩邊同乘\(A(x)\)
\(A(x)B(x)B(x)-2B_0(x)A(x)B(x)+A(x)B_0^2(x)\equiv 0(\mod x^{2n})\)
\(B(x)\equiv 2B_0(x)-A(x)B_0^2(x) (\mod x^{2n})\)
右邊的式子FFT計算即可。
時間復雜度？
\(T(n)=T(\frac{n}{2})+O(n\log n)\)
\(T(n)=O(n\log n)\).
常數？首先隱藏在遞歸復雜度背后有一個\(2\)倍常數。
然后我們把兩個多項式相乘需要\(3\)次FFT, 三個就要\(6\)次。
因此常數為\(12\)倍。
如何優化？
\(IDFT(DFT(IDFT(DFT(A)\times DFT(B_0)))\times DFT(B_0))\)
變成\(IDFT(DFT(A)\times DFT^2(B_0)\)
\(3\)次即可！
常數變為\(6\)倍。
UPD: 關于這里的常數問題: 因為我遞歸里DFT的范圍是\(2n\)，最終的復雜度是\(T(2n) = 6(2n)\log (2n)\), 因此我認為應為\(12\)倍常數。
空間？空間復雜度\(O(n)\), 但是要開\(4d\)的數組，其中\(d\)是\(>n\)的最小的\(2\)的冪。

代碼
因為FFT數組清零等原因代碼(對我來說)很難寫。
貼一下我剛剛寫的版本吧，還算是比較簡單。
(話說怎么CSDN突然傲嬌了啊。。縮進變成1格？？)

#include<cstdio> #include<cstdlib> #include<cstring> #include<algorithm> #define llong long long #define ldouble long double #define uint unsigned int #define ullong unsigned long long #define udouble unsigned double #define uldouble unsigned long double #define modinc(x) {if(x>=P) x-=P;} #define pii pair<int,int> #define piii pair<pair<int,int>,int> #define piiii pair<pair<int,int>,pair<int,int> > #define pli pair<llong,int> #define pll pair<llong,llong> #define Memset(a,x) {memset(a,x,sizeof(a));} using namespace std;const int N = 1<<19; const int P = 998244353; const int LGN = 19; const int G = 3; llong a[N+3]; llong b[N+3]; llong tmp1[N+3],tmp2[N+3],tmp3[N+3],tmp4[N+3],tmp5[N+3],tmp6[N+3]; int id[N+2]; int n;void initid(int _len) {id[0] = 0;for(int i=1; i<(1<<_len); i++) id[i] = (id[i>>1]>>1)|((i&1)<<(_len-1)); }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);}void ntt(int dgr,int coe,llong poly[],llong ret[]) {int len = 0; for(int i=0; i<=LGN; i++) if((1<<i)==dgr) {len = i; break;}initid(len); for(int i=0; i<dgr; i++) ret[i] = 0ll;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 = (expn*ret[j+i+k])%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;}}} }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++) tmp5[j] = tmp3[j]*tmp3[j]%P*tmp4[j]%P;ntt((i<<2),-1,tmp5,tmp6); llong tmp = mulinv(i<<2);for(int j=0; j<(i<<2); j++) tmp6[j] = tmp6[j]*tmp%P;for(int j=0; j<(i<<1); j++) ret[j] = (tmp1[j]+tmp1[j]-tmp6[j]+P)%P;} }int main() {scanf("%d",&n); int dgr = 1; while(dgr<=n) dgr<<=1;for(int i=0; i<n; i++) scanf("%lld",&a[i]);polyinv(dgr,a,b);for(int i=0; i<n; i++) printf("%lld ",b[i]);return 0; } 發表于 2019-01-23 20:19 suncongbo 閱讀(...) 評論(...) 編輯 收藏 刷新評論刷新頁面返回頂部

總結

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

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