寫了一發多項式exp，picks太神辣

  1 #include <bits/stdc++.h>
  2 
  3 using namespace std;
  4 
  5 typedef long long ll;
  6 const int N = 262144,mod = 998244353,G = 3;
  7 
  8 void read(int&x){
  9     x = 0;char ch = getchar();
 10     while (ch < ‘0‘ || ‘9‘ < ch) ch = getchar();
 11     while (‘
 
0‘ <= ch&& ch <= ‘9‘) x = x * 10 + ch - ‘0‘,ch = getchar(); 
 12 }
 13 int power(int x,int y,int z){
 14     int ans = 1;
 15     while (y){
 16         if (y&1) ans = (ll)ans*x%z;
 17         x = (ll)x*x%z;
 18         y >>= 1;
 19     }
 20     return ans;
 21 }
 22 int add(int
 
 x,int y){
 23     x += y;
 24     while (x >= mod) x -= mod;
 25     return x;
 26 }
 27 
 28 
 29 int a[N+10],b[N+10],c[N+10],d[N+10],p[N];
 30 int pg[N+10],pg_inv[N+10],inv[N+10];
 31 int n,m;
 32 void FFT(int*,int,int);
 33 void poly_exp(int*,int*,int*,int*,int);
 34 void poly_inv(int*,int*,int*,int);
 
 35 void poly_ln(int*,int*,int*,int);
 36 void poly_der(int*,int);
 37 void poly_mot(int*,int); 
 38 int main(){
 39     read(n);
 40     for (int i = 0;i < n;i++) read(a[n-i-1]);
 41     for (m = 1;m < n;m<<=1);
 42     for(int i = 1;i <= (m<<1);i <<= 1){
 43         pg[i] = power(G,(mod-1)/i,mod);
 44         pg_inv[i] = power(pg[i],mod-2,mod);
 45     }
 46     inv[1] = 1;
 47     for (int i = 2;i <= N;i++)inv[i] = mod-(ll)mod/i*inv[mod%i]%mod;
 48     /* 
 49     poly_inv(a,b,c,m); 
 50 //    for (int i = 0;i < m<<1;i++)cout << a[i]<<‘ ‘;cout<<endl;
 51 //    for (int i = 0;i < m<<1;i++)cout << b[i]<<‘ ‘;cout<<endl;
 52 //    cout <<(ll)a[0]*b[0]%mod<<endl;
 53     FFT(a,1,m<<1);
 54 //    for (int i = 0;i < m<<1;i++)cout <<a[i]<<‘*‘;cout<<endl;
 55     FFT(b,1,m<<1);
 56     for(int i = 0;i <m<<1;i++) a[i] =(ll)a[i]*b[i]%mod;
 57     FFT(a,-1,m<<1);
 58 //    for (int i = 0;i < m;i++)cout << a[i]<<‘ ‘;cout<<endl;
 59     bool flag = 1;
 60     for(int i = 1;i < m;i++)if (a[i] != 0) flag = 0;
 61     if(a[0] != 1) flag = 0;
 62     puts(flag ? "YES" : "NO");
 63     */ 
 64 //    /*
 65     poly_exp(a,b,c,d,m);
 66 //    for (int i = 0;i < m<<1;i++) cout << a[i]<<‘ ‘;cout << endl;
 67 //    for (int i = 0;i < m<<1;i++)cout << b[i] <<‘ ‘;cout<<endl;
 68     /*
 69 //    b[0] = 1;b[1] = 1;b[2] = inv[2];b[3] = inv[6];
 70     poly_ln(b,c,d,m);bool flag = 1;
 71 //    for (int i = 0;i < m;i++) cout << c[i]<<‘ ‘;cout<<endl;
 72     for (int i = 0;i <m;i++)if(c[i] != a[i]) flag = 0;
 73     puts(flag ? "YES" : "NO");
 74     */
 75 //    */
 76     return 0;
 77 }
 78 void FFT(int *a,int d,int deg){
 79     for (int i = 1;i < deg;i <<= 1)
 80         for (int j = 0;j < i;j++)
 81             p[i+j] = p[j]+deg/i/2;
 82     for (int i = 0;i <deg;i++) if (p[i] > i) swap(a[i],a[p[i]]);
 83     for (int n = 2;n <= deg;n <<= 1){
 84         int wn = (d == 1 ? pg[n] : pg_inv[n]);
 85         for (int i = 0;i <deg;i += n){
 86             int w = 1;
 87             for (int j = i;j < i + n/2;j++){
 88                 int x = a[j],y = (ll)a[j+n/2]*w%mod;
 89                 a[j] = add(x,y);
 90                 a[j+n/2] = add(x,mod-y);
 91                 w = (ll)w*wn%mod;
 92             }
 93         }
 94     }
 95     if (d == -1){
 96         for (int i = 0;i < deg;i++)
 97             a[i] = (ll)a[i]*inv[deg]%mod;
 98     }
 99 }
100 void poly_der(int*a,int deg){
101     for (int i = 0;i < deg;i++)
102         a[i] = (ll)a[i+1]*(i+1)%mod;
103 }
104 void poly_mot(int*a,int deg){
105     for (int i = deg;i > 0;i--)
106         a[i] = (ll)a[i-1]*inv[i]%mod;
107     a[0] = 0;
108 }
109 void poly_inv(int*a,int*b,int*c,int deg){
110     if (deg == 1){
111         if (a[0] <= N) b[0] = inv[a[0]];else 
112         b[0] = power(a[0],mod-2,mod);
113         return;
114     }
115     poly_inv(a,b,c,deg/2);
116     int deg2 = deg<<1;
117     //B = B(2-BA)
118     for (int i = 0;i < deg;i++) c[i] = a[i];
119     for (int i = deg;i < deg2;i++) c[i] = 0;
120     FFT(c,1,deg2);
121     FFT(b,1,deg2);
122     for (int i = 0;i < deg2;i++) b[i] = (ll)b[i]*add(2,mod-(ll)b[i]*c[i]%mod)%mod;
123     FFT(b,-1,deg2);
124     for (int i = deg;i < deg2;i++) b[i] = 0;
125 }
126 void poly_ln(int*a,int*b,int*c,int deg){
127     /*
128         b = lna
129         b‘ = a‘/a
130     */
131     int deg2 = deg<<1;
132     for (int i = 0;i < deg2;i++)c[i] = 0;
133     poly_inv(a,c,b,deg);
134     for(int i = 0;i <= deg;i++) b[i] = a[i];
135     for(int i = deg+1;i < deg2;i++) b[i] = 0;
136     poly_der(b,deg);
137     FFT(b,1,deg2);
138     FFT(c,1,deg2);
139     for (int i = 0;i < deg2;i++) b[i] = (ll)b[i]*c[i]%mod;
140     FFT(b,-1,deg2);
141     for (int i = deg;i < deg2;i++) b[i] = 0;
142     poly_mot(b,deg);
143     for(int i = deg;i < deg2;i++)b[i] = 0;
144 }
145 void poly_exp(int *a,int *b,int *c,int*d,int deg){
146     //B = B(1-lnB+A)
147     if(deg == 1){
148         //if a[0] = 0
149         b[0] = 1;
150         return;
151     }
152     int deg2=deg<<1;    
153     poly_exp(a,b,c,d,deg/2);
154     for (int i = deg;i < deg2;i++) b[i] = 0;
155     poly_ln(b,c,d,deg);
156 
157     for (int i = 0;i < deg;i++) c[i] = add(a[i]+(i==0),mod-c[i]);
158 //    cout <<"1-lnB+A ";for (int i = 0;i <deg2;i++) cout << c[i]<<‘ ‘;cout<<endl;    
159     for (int i = deg;i < deg2;i++) c[i] = 0;
160     FFT(c,1,deg2);
161     FFT(b,1,deg2);
162     for(int i = 0;i < deg2;i++) b[i] = (ll)b[i]*c[i]%mod;
163     FFT(b,-1,deg2);
164     for (int i = deg;i < deg2;i++) b[i] = 0;
165 //    for (int i = 0;i < deg2;i++)cout << b[i]<<‘ ‘;cout<<"deg:"<<deg<<"\n";
166 }

因為只會求e^0的值，所以只能解決常數項為0的多項式的exp。。。。

多項式exp其實求的是e^A的泰勒展開的前若幹項，

用牛頓叠代法，求出解的遞推式

最後一行在多項式的意義下在這個情況中證明了這個方法至少是平方收斂（多項式下的牛頓叠代應該都是可以類似證精度的）（數意義下的不會證orz）

然後。。我試圖用原始的初等的推多項式求逆的那一套方法推這個式子，成功GG。。。。

然後寫的時候發現一些問題，就是這個求原函數啊，求模x^n意義下的原函數需要模x^(n+1)意義下的函數，這題因為需要求原函數的函數都是沒取模的所以不影響，不然還有註意一下細節。。。

多項式exp