luogu P4725 多項式對數函數 (模板題、FFT、多項式求逆、求導和積分)

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

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

題目大意: 給定一個\(n\)次多項式\(A(x)\), 求一個\(n\)次多項式\(B(x)\)滿足\(B(x)\equiv \ln A(x) (\mod x^n)\)

題解: 神數學模板題……
數學真奇妙！
前驅知識
導數、積分相關
冪函數的求導
\(f(x)=x^n, f'(x)=nx^{n-1}\)
和的導數等于導數的和
\((f+g)'(x)=f'(x)+g'(x)\)
一般多項式的求導
\(f(x)=\sum^{n-1}_{i=0} a_ix^i, f'(x)=\sum^{n-2}_{i=0} (i+1)a_{i+1}x^i\)
對數函數\(\ln\)的求導
\(f(x)=\ln(x), f'(x)=\frac{1}{x}\)
復合函數求導——鏈式法則
\(f(g(x))'=f'(g(x))g'(x)\)
求導的逆運算——積分

本題解法
\(g(x)\equiv \ln f(x) (\mod x^n)\)
兩邊同時求導可得
\(g'(x)\equiv \frac{f'(x)}{f(x)} (\mod x^n)\)
結束！
多項式求逆算\(\frac{1}{f(x)}\)，再和\(f'(x)\)相乘即可得到\(g'(x)\)。
(多項式求逆見蒟蒻一篇博客 https://blog.csdn.net/suncongbo/article/details/84485718)
\(g'(x)\)積個分得到\(g(x)\). 常數項，直接為\(0\).
時間復雜度\(O(n\log n)\)
常數，我寫的大概\(9\)倍吧，求逆是\(6\)倍，再做個乘法就是\(3\)倍。
UPD: 仔細想了一下我這個常數好像是\(18\)倍（見我多項式求逆那篇博客）

代碼

#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]; //inv llong tmp7[N+3],tmp8[N+3],tmp9[N+3],tmp10[N+3]; //ln 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;}}}if(coe==-1){llong tmp = mulinv(dgr);for(int i=0; i<dgr; i++) ret[i] = ret[i]*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 i=dgr; i<(dgr<<1); i++) ret[i] = 0ll; }void polyder(int dgr,llong poly[],llong ret[]) {for(int i=0; i<dgr-1; i++) ret[i] = poly[i+1]*(i+1)%P; }void polyint(int dgr,llong poly[],llong ret[]) {for(int i=1; i<=dgr; i++) ret[i] = poly[i-1]*mulinv(i)%P; }void polyln(int dgr,llong poly[],llong ret[]) {polyder(dgr,poly,tmp7);polyinv(dgr,poly,tmp8);ntt((dgr<<1),1,tmp8,tmp9); ntt((dgr<<1),1,tmp7,tmp10);for(int i=0; i<(dgr<<1); i++) tmp9[i] = tmp9[i]*tmp10[i]%P;ntt((dgr<<1),-1,tmp9,tmp10);polyint(dgr,tmp10,ret); }int main() {scanf("%d",&n); int dgr = 1; while(dgr<=n) dgr<<=1;for(int i=0; i<n; i++) scanf("%lld",&a[i]);polyln(dgr,a,b);for(int i=0; i<n; i++) printf("%lld ",b[i]);return 0; } 發表于 2019-01-23 20:20 suncongbo 閱讀(...) 評論(...) 編輯 收藏 刷新評論刷新頁面返回頂部

總結

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

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