• 关于万恶的多项式


    模板:

     1 //Achen
     2 #include<algorithm>
     3 #include<iostream>
     4 #include<cstring>
     5 #include<cstdlib>
     6 #include<vector>
     7 #include<cstdio>
     8 #include<queue>
     9 #include<cmath>
    10 #include<set>
    11 #include<map>
    12 #define pi acos(-1)
    13 #define Formylove return 0
    14 #define For(i,a,b) for(int i=(a);i<=(b);i++)
    15 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    16 const int N=2097155; 
    17 typedef long long LL;
    18 typedef double db;
    19 using namespace std;
    20 int n,m,ans[N];
    21 
    22 template<typename T>void read(T &x)  {
    23     char ch=getchar(); x=0; T f=1;
    24     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
    25     if(ch=='-') f=-1,ch=getchar();
    26     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    27 }
    28 
    29 struct E {
    30     db a,b;
    31     E(db a=0,db b=0):a(a),b(b){}
    32 }a[N],b[N];
    33 E operator +(const E&A,const E&B) { return E(A.a+B.a,A.b+B.b); }
    34 E operator -(const E&A,const E&B) { return E(A.a-B.a,A.b-B.b); }
    35 E operator *(const E&A,const E&B) { return E(A.a*B.a-A.b*B.b,A.a*B.b+A.b*B.a); }
    36 E operator /(const E&A,const db&B) { return E(A.a/B,A.b/B); }
    37 
    38 int rev[N];
    39 void FFT(E a[],int n,int l,int f) {
    40     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    41     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
    42     for(int i=1;i<n;i<<=1) {
    43         E wi(cos(pi/i),f*sin(pi/i));
    44         for(int j=0,pp=(i<<1);j<n;j+=pp) {
    45             E w(1,0);
    46             for(int k=0;k<i;k++,w=w*wi) {
    47                 E x=a[j+k],y=w*a[j+k+i];
    48                 a[j+k]=x+y; a[j+k+i]=x-y;
    49             }
    50         }
    51     }
    52     if(f==-1) For(i,0,n) a[i].a=(int)(a[i].a/n+0.5);
    53 }
    54 
    55 void mul(E a[],E b[],int c[],int n,int m) {
    56     static E A[N],B[N];
    57     int nn,l;
    58     for(nn=1,l=0;nn<=n+m;nn<<=1) l++;
    59     For(i,0,n) A[i]=a[i]; For(i,n+1,nn) A[i]=E(0,0);
    60     For(i,0,m) B[i]=b[i]; For(i,m+1,nn) B[i]=E(0,0);
    61     FFT(A,nn,l,1); FFT(B,nn,l,1);
    62     For(i,0,nn) A[i]=A[i]*B[i];
    63     FFT(A,nn,l,-1);
    64     For(i,0,n+m) c[i]=A[i].a;
    65 }
    66 
    67 int main() {
    68 #ifdef ANS
    69     freopen(".in","r",stdin);
    70     freopen(".out","w",stdout);
    71 #endif
    72     read(n); read(m);
    73     For(i,0,n) read(a[i].a);
    74     For(i,0,m) read(b[i].a);
    75     mul(a,b,ans,n,m);
    76     For(i,0,n+m) printf("%d ",ans[i]);
    77     Formylove;
    78 }
    FFT
     1 //Achen
     2 #include<algorithm>
     3 #include<iostream>
     4 #include<cstring>
     5 #include<cstdlib>
     6 #include<vector>
     7 #include<cstdio>
     8 #include<queue>
     9 #include<cmath>
    10 #include<set>
    11 #include<map>
    12 #define Formylove return 0
    13 #define For(i,a,b) for(int i=(a);i<=(b);i++)
    14 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    15 const int N=2097155,p=998244353,g=3,gi=332748118; 
    16 typedef long long LL;
    17 typedef double db;
    18 using namespace std;
    19 int n,m;
    20 LL a[N],b[N],ans[N];
    21 
    22 template<typename T>void read(T &x)  {
    23     char ch=getchar(); x=0; T f=1;
    24     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
    25     if(ch=='-') f=-1,ch=getchar();
    26     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    27 }
    28 
    29 LL ksm(LL a,LL b) {
    30     LL rs=1,bs=a;
    31     while(b) {
    32         if(b&1) rs=rs*bs%p;
    33         bs=bs*bs%p;
    34         b>>=1;
    35     }
    36     return rs;
    37 }
    38 
    39 LL mo(LL x) { return x<0?x+p:(x>=p?x-p:x); }
    40 
    41 int rev[N];
    42 void FNT(LL a[],int n,int l,int f) {
    43     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    44     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
    45     for(int i=1;i<n;i<<=1) {
    46         LL wi=ksm((f==1)?g:gi,(p-1)/(i<<1));
    47         for(int j=0,pp=(i<<1);j<n;j+=pp) {
    48             LL w=1;
    49             for(int k=0;k<i;k++,w=w*wi%p) {
    50                 LL x=a[j+k],y=w*a[j+k+i]%p;
    51                 a[j+k]=mo(x+y); a[j+k+i]=mo(x-y);
    52             }
    53         }
    54     }
    55     if(f==-1) {
    56         LL inv_n=ksm(n,p-2);
    57         For(i,0,n) a[i]=a[i]*inv_n%p;
    58     }
    59 }
    60 
    61 void mul(LL a[],LL b[],LL c[],int n,int m) {
    62     static LL A[N],B[N];
    63     int nn,l;
    64     for(nn=1,l=0;nn<=n+m;nn<<=1) l++;
    65     For(i,0,n) A[i]=a[i]; For(i,n+1,nn) A[i]=0;
    66     For(i,0,m) B[i]=b[i]; For(i,m+1,nn) B[i]=0;
    67     FNT(A,nn,l,1); FNT(B,nn,l,1);
    68     For(i,0,nn) A[i]=A[i]*B[i]%p;
    69     FNT(A,nn,l,-1);
    70     For(i,0,n+m) c[i]=A[i];
    71 }
    72 
    73 int main() {
    74 #ifdef ANS
    75     freopen(".in","r",stdin);
    76     freopen(".out","w",stdout);
    77 #endif
    78     read(n); read(m);
    79     For(i,0,n) read(a[i]);
    80     For(i,0,m) read(b[i]);
    81     mul(a,b,ans,n,m);
    82     For(i,0,n+m) printf("%lld ",ans[i]);
    83     Formylove;
    84 }
    FNT
     1 //易错:get_inv中传入的b必须是空的,空间开到 get_inv中的n的2倍 
     2 //Achen
     3 #include<algorithm>
     4 #include<iostream>
     5 #include<cstring>
     6 #include<cstdlib>
     7 #include<vector>
     8 #include<cstdio>
     9 #include<queue>
    10 #include<cmath>
    11 #include<set>
    12 #include<map>
    13 #define Formylove return 0
    14 #define For(i,a,b) for(int i=(a);i<=(b);i++)
    15 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    16 const int N=262155,p=998244353,g=3,gi=332748118; 
    17 typedef long long LL;
    18 typedef double db;
    19 using namespace std;
    20 int n;
    21 LL a[N],ans[N];
    22 
    23 template<typename T>void read(T &x)  {
    24     char ch=getchar(); x=0; T f=1;
    25     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
    26     if(ch=='-') f=-1,ch=getchar();
    27     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    28 }
    29  
    30 LL ksm(LL a,LL b) {
    31     LL rs=1,bs=a;
    32     while(b) {
    33         if(b&1) rs=rs*bs%p;
    34         bs=bs*bs%p;
    35         b>>=1;
    36     }
    37     return rs;
    38 }
    39 
    40 LL mo(LL x) { return x<0?x+p:(x>=p?x-p:x); }
    41 
    42 int rev[N];
    43 void FNT(LL a[],int n,int l,int f) {
    44     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    45     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
    46     for(int i=1;i<n;i<<=1) {
    47         LL wi=ksm((f==1)?g:gi,(p-1)/(i<<1));
    48         for(int j=0,pp=(i<<1);j<n;j+=pp) {
    49             LL w=1;
    50             for(int k=0;k<i;k++,w=w*wi%p) {
    51                 LL x=a[j+k],y=w*a[j+k+i]%p;
    52                 a[j+k]=mo(x+y); a[j+k+i]=mo(x-y);
    53             }
    54         }
    55     }
    56     if(f==-1) {
    57         LL inv_n=ksm(n,p-2);
    58         For(i,0,n) a[i]=a[i]*inv_n%p;
    59     }
    60 }
    61 
    62 LL tp[N];
    63 void get_inv(LL a[],LL b[],int n,int l) {
    64     if(n==0) {
    65         b[0]=ksm(a[0],p-2);
    66         return ;
    67     }
    68     get_inv(a,b,n>>1,l-1);
    69     For(i,0,n-1) tp[i]=a[i],tp[i+n]=0;
    70     FNT(tp,n<<1,l+1,1);
    71     FNT(b,n<<1,l+1,1);
    72     For(i,0,(n<<1)) tp[i]=mo(2LL-tp[i]*b[i]%p)*b[i]%p;
    73     FNT(tp,n<<1,l+1,-1);
    74     For(i,0,n-1) b[i]=tp[i],b[i+n]=0;
    75 }
    76 
    77 int main() {
    78 #ifdef ANS
    79     freopen(".in","r",stdin);
    80     freopen(".out","w",stdout);
    81 #endif
    82     read(n);
    83     For(i,0,n-1) read(a[i]);
    84     int nn,l;
    85     for(nn=1,l=0;nn<=n;nn<<=1) l++;
    86     get_inv(a,ans,nn,l);
    87     For(i,0,n-1) printf("%lld ",ans[i]);
    88     Formylove;
    89 }
    多项式求逆
     1 //Achen
     2 #include<algorithm>
     3 #include<iostream>
     4 #include<cstring>
     5 #include<cstdlib>
     6 #include<vector>
     7 #include<cstdio>
     8 #include<queue>
     9 #include<cmath>
    10 #include<set>
    11 #include<map>
    12 #define Formylove return 0
    13 #define For(i,a,b) for(int i=(a);i<=(b);i++)
    14 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    15 const int N=262155,p=998244353,g=3,gi=332748118; 
    16 typedef long long LL;
    17 typedef double db;
    18 using namespace std;
    19 int n;
    20 LL a[N],b[N],tp[N];
    21 
    22 template<typename T>void read(T &x)  {
    23     char ch=getchar(); x=0; T f=1;
    24     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
    25     if(ch=='-') f=-1,ch=getchar();
    26     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    27 }
    28 
    29 LL ksm(LL a,LL b) {
    30     LL rs=1,bs=a;
    31     while(b) {
    32         if(b&1) rs=rs*bs%p;
    33         bs=bs*bs%p;
    34         b>>=1; 
    35     }    
    36     return rs;
    37 }
    38 
    39 LL mo(LL x) { return x<0?x+p:(x>=p?x-p:x); }
    40 int rev[N];
    41 void FNT(LL a[],int n,int l,int f) {
    42     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    43     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
    44     for(int i=1;i<n;i<<=1) {
    45         LL wi=ksm((f==1?g:gi),(p-1)/(i<<1));
    46         for(int j=0,pp=(i<<1);j<n;j+=pp) {
    47             LL w=1;
    48             for(int k=0;k<i;k++,w=w*wi%p) {
    49                 LL x=a[j+k],y=w*a[j+k+i]%p;
    50                 a[j+k]=mo(x+y); a[j+k+i]=mo(x-y);
    51             }
    52         } 
    53     }
    54     if(f==-1) {
    55         LL inv_n=ksm(n,p-2);
    56         For(i,0,n) a[i]=a[i]*inv_n%p;
    57     }
    58 }
    59 
    60 void get_inv(LL a[],LL b[],LL tp[],int n,int l) {
    61     if(n==1) { b[0]=ksm(a[0],p-2); return; }
    62     get_inv(a,b,tp,(n>>1),l-1);
    63     For(i,0,n-1) tp[i]=a[i],tp[i+n]=0;
    64     FNT(tp,(n<<1),l+1,1);
    65     FNT(b,(n<<1),l+1,1);
    66     For(i,0,(n<<1)) tp[i]=mo(2LL-tp[i]*b[i]%p)*b[i]%p;
    67     FNT(tp,(n<<1),l+1,-1);
    68     For(i,0,n-1) b[i]=tp[i],b[i+n]=0;
    69 }
    70 
    71 LL inv[N];
    72 void get_ln(LL a[],LL b[],int n,int l) {
    73     static LL a_dao[N],a_inv[N];
    74     For(i,0,(n<<1)) a_inv[i]=a_dao[i]=0;
    75     For(i,0,n-1) a_dao[i]=a[i+1]*(i+1)%p;
    76     get_inv(a,a_inv,tp,n,l);
    77     FNT(a_inv,(n<<1),l+1,1);
    78     FNT(a_dao,(n<<1),l+1,1);
    79     For(i,0,(n<<1)) a_inv[i]=a_inv[i]*a_dao[i]%p;
    80     FNT(a_inv,(n<<1),l+1,-1);
    81     b[0]=0;
    82     For(i,1,n-1) b[i]=a_inv[i-1]*inv[i]%p;
    83 }
    84 
    85 int main() {
    86 #ifdef ANS
    87     freopen(".in","r",stdin);
    88     freopen(".out","w",stdout);
    89 #endif
    90     read(n);
    91     For(i,0,n-1) read(a[i]);
    92     int nn,l;
    93     for(nn=1,l=0;nn<=n;nn<<=1) l++;
    94     inv[0]=inv[1]=1;
    95     For(i,2,n) inv[i]=mo(p-p/i*inv[p%i]%p);
    96     get_ln(a,b,nn,l);
    97     For(i,0,n-1) printf("%lld ",b[i]);
    98     Formylove;
    99 }
    多项式对数函数(ln)
      1 //Achen
      2 #include<algorithm>
      3 #include<iostream>
      4 #include<cstring>
      5 #include<cstdlib>
      6 #include<vector>
      7 #include<cstdio>
      8 #include<queue>
      9 #include<cmath>
     10 #include<set>
     11 #include<map>
     12 #define Formylove return 0
     13 #define For(i,a,b) for(int i=(a);i<=(b);i++)
     14 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
     15 const int N=262155,p=998244353,g=3,gi=332748118; 
     16 typedef long long LL;
     17 typedef double db;
     18 using namespace std;
     19 int n;
     20 LL a[N],b[N];
     21 
     22 template<typename T>void read(T &x)  {
     23     char ch=getchar(); x=0; T f=1;
     24     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
     25     if(ch=='-') f=-1,ch=getchar();
     26     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
     27 }
     28 
     29 LL ksm(LL a,LL b) {
     30     LL rs=1,bs=a;
     31     while(b) {
     32         if(b&1) rs=rs*bs%p;
     33         bs=bs*bs%p;
     34         b>>=1; 
     35     }    
     36     return rs;
     37 }
     38 
     39 LL mo(LL x) { return x<0?x+p:(x>=p?x-p:x); }
     40 int rev[N];
     41 void FNT(LL a[],int n,int l,int f) {
     42     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
     43     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
     44     for(int i=1;i<n;i<<=1) {
     45         LL wi=ksm((f==1?g:gi),(p-1)/(i<<1));
     46         for(int j=0,pp=(i<<1);j<n;j+=pp) {
     47             LL w=1;
     48             for(int k=0;k<i;k++,w=w*wi%p) {
     49                 LL x=a[j+k],y=w*a[j+k+i]%p;
     50                 a[j+k]=mo(x+y); a[j+k+i]=mo(x-y);
     51             }
     52         } 
     53     }
     54     if(f==-1) {
     55         LL inv_n=ksm(n,p-2);
     56         For(i,0,n) a[i]=a[i]*inv_n%p;
     57     }
     58 }
     59 
     60 LL tp[N];
     61 void get_inv(LL a[],LL b[],LL tp[],int n,int l) {
     62     if(n==1) { b[0]=ksm(a[0],p-2); return; }
     63     get_inv(a,b,tp,(n>>1),l-1);
     64     For(i,0,n-1) tp[i]=a[i],tp[i+n]=0;
     65     FNT(tp,(n<<1),l+1,1);
     66     FNT(b,(n<<1),l+1,1);
     67     For(i,0,(n<<1)) tp[i]=mo(2LL-tp[i]*b[i]%p)*b[i]%p;
     68     FNT(tp,(n<<1),l+1,-1);
     69     For(i,0,n-1) b[i]=tp[i],b[i+n]=0;
     70 }
     71 
     72 LL inv[N];
     73 void get_ln(LL a[],LL b[],int n,int l) {
     74     static LL a_dao[N],a_inv[N];
     75     For(i,0,(n<<1)) a_inv[i]=a_dao[i]=0;
     76     For(i,0,n-1) a_dao[i]=a[i+1]*(i+1)%p; 
     77     a_dao[n-1]=0;
     78     get_inv(a,a_inv,tp,n,l);
     79     FNT(a_inv,(n<<1),l+1,1);
     80     FNT(a_dao,(n<<1),l+1,1);
     81     For(i,0,(n<<1)) a_inv[i]=a_inv[i]*a_dao[i]%p;
     82     FNT(a_inv,(n<<1),l+1,-1);
     83     b[0]=0;
     84     For(i,1,n-1) b[i]=a_inv[i-1]*inv[i]%p;
     85 }
     86 
     87 LL ln_b[N];
     88 void get_exp(LL a[],LL b[],int n,int l) {
     89     if(n==1) { b[0]=1; return; }
     90     get_exp(a,b,(n>>1),l-1);
     91     For(i,0,(n<<1)) ln_b[i]=0;
     92     get_ln(b,ln_b,n,l); ////
     93     For(i,0,n-1) ln_b[i]=((LL)(i==0)-ln_b[i]+a[i]+p)%p;
     94     FNT(b,(n<<1),l+1,1);
     95     FNT(ln_b,(n<<1),l+1,1);
     96     For(i,0,(n<<1)) ln_b[i]=ln_b[i]*b[i]%p;
     97     FNT(ln_b,(n<<1),l+1,-1);
     98     For(i,0,n-1) b[i]=ln_b[i],b[i+n]=0; 
     99 }
    100 
    101 int main() {
    102 #ifdef ANS
    103     freopen(".in","r",stdin);
    104     freopen(".out","w",stdout);
    105 #endif
    106     read(n);
    107     For(i,0,n-1) read(a[i]);
    108     int nn,l;
    109     for(nn=1,l=0;nn<=n;nn<<=1) l++;
    110     inv[0]=inv[1]=1;
    111     For(i,2,n) inv[i]=mo(p-p/i*inv[p%i]%p);
    112     get_exp(a,b,nn,l);
    113     For(i,0,n-1) printf("%lld ",b[i]);
    114     Formylove;
    115 }
    多项式指数函数(exp)
     1 LL h_sqr[N],a_t[N],b_t[N];
     2 void get_sqr(LL a[],LL b[],int n,int l) {
     3     if(n==1) {
     4         b[0]=sqrt(a[0]);
     5         return;
     6     }
     7     get_sqr(a,b,n>>1,l-1);
     8     For(i,0,n<<1) b_t[i]=0;
     9     For(i,0,n-1) a_t[i]=a[i],a_t[i+n]=0;
    10     get_inv(b,b_t,n,l);
    11     FFT(b,n<<1,l+1,1);
    12     FFT(a_t,n<<1,l+1,1);
    13     FFT(b_t,n<<1,l+1,1);
    14     For(i,0,n<<1) b[i]=(b[i]+a_t[i]*b_t[i]%p)%p*inv_2%p;
    15     FFT(b,n<<1,l+1,-1);
    16     For(i,0,n-1) b[n+i]=0;
    17 }
    多项式开方(牛顿迭代)
      1 //Achen
      2 #include<algorithm>
      3 #include<iostream>
      4 #include<cstring>
      5 #include<cstdlib>
      6 #include<vector>
      7 #include<cstdio>
      8 #include<queue>
      9 #include<cmath>
     10 #include<set>
     11 #include<map>
     12 #define Formylove return 0
     13 #define For(i,a,b) for(int i=(a);i<=(b);i++)
     14 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
     15 const int N=262155,p=998244353,g=3,gi=332748118; 
     16 typedef long long LL;
     17 typedef double db;
     18 using namespace std;
     19 int n,m;
     20 LL a[N],b[N],D[N],R[N];
     21 
     22 template<typename T>void read(T &x)  {
     23     char ch=getchar(); x=0; T f=1;
     24     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
     25     if(ch=='-') f=-1,ch=getchar();
     26     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
     27 }
     28 
     29 LL ksm(LL a,LL b) {
     30     LL rs=1,bs=a;
     31     while(b) {
     32         if(b&1) rs=rs*bs%p;
     33         bs=bs*bs%p;
     34         b>>=1; 
     35     }    
     36     return rs;
     37 }
     38 
     39 LL mo(LL x) { return x<0?x+p:(x>=p?x-p:x); }
     40 int rev[N];
     41 void FNT(LL a[],int n,int l,int f) {
     42     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
     43     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
     44     for(int i=1;i<n;i<<=1) {
     45         LL wi=ksm((f==1?g:gi),(p-1)/(i<<1));
     46         for(int j=0,pp=(i<<1);j<n;j+=pp) {
     47             LL w=1;
     48             for(int k=0;k<i;k++,w=w*wi%p) {
     49                 LL x=a[j+k],y=w*a[j+k+i]%p;
     50                 a[j+k]=mo(x+y); a[j+k+i]=mo(x-y);
     51             }
     52         } 
     53     }
     54     if(f==-1) {
     55         LL inv_n=ksm(n,p-2);
     56         For(i,0,n) a[i]=a[i]*inv_n%p;
     57     }
     58 }
     59 
     60 LL tp[N];
     61 void get_inv(LL a[],LL b[],int n,int l) {
     62     if(n==1) { 
     63         b[0]=ksm(a[0],p-2); 
     64         return; 
     65     }
     66     get_inv(a,b,(n>>1),l-1);
     67     For(i,0,n-1) tp[i]=a[i],tp[i+n]=0;
     68     FNT(tp,(n<<1),l+1,1);
     69     FNT(b,(n<<1),l+1,1);
     70     For(i,0,(n<<1)) tp[i]=mo(2LL-tp[i]*b[i]%p)*b[i]%p;
     71     FNT(tp,(n<<1),l+1,-1);
     72     For(i,0,n-1) b[i]=tp[i],b[i+n]=0;
     73 }
     74 
     75 void get_R(LL a[],LL a_R[],int n) { For(i,0,n) a_R[i]=a[n-i]; }
     76 
     77 void mul(LL a[],LL b[],LL c[],int n,int m) {
     78     static LL A_l[N],B_l[N];
     79     int nn,l;
     80     for(nn=1,l=0;nn<=n+m;nn<<=1) l++;
     81     For(i,0,n) A_l[i]=a[i]; For(i,n+1,nn) A_l[i]=0;
     82     For(i,0,m) B_l[i]=b[i]; For(i,m+1,nn) B_l[i]=0;
     83     FNT(A_l,nn,l,1);
     84     FNT(B_l,nn,l,1);
     85     For(i,0,nn) (A_l[i]*=B_l[i])%=p;
     86     FNT(A_l,nn,l,-1);
     87     For(i,0,n+m) c[i]=A_l[i];
     88 } 
     89 
     90 void div(LL a[],int n,LL b[],int m,LL D[],LL R[]) {
     91     static LL A_R[N],B_R[N],D_R[N];
     92     int nn,l;
     93     for(nn=1,l=0;nn<=n-m+1;nn<<=1) l++;
     94     For(i,0,(nn<<1)) D_R[i]=0;
     95     get_R(a,A_R,n); 
     96     get_R(b,B_R,m);
     97     get_inv(B_R,D_R,nn,l);
     98     mul(A_R,D_R,D_R,n,n-m);
     99     get_R(D_R,D,n-m);
    100     for(nn=1,l=0;nn<=(n<<1);nn<<=1) l++;
    101     mul(b,D,R,m,n-m);
    102     For(i,0,nn) R[i]=mo(a[i]-R[i]);
    103 }
    104 
    105 //#define ANS
    106 int main() {
    107 #ifdef ANS
    108     freopen("1.in","r",stdin);
    109     //freopen(".out","w",stdout);
    110 #endif
    111     read(n); read(m);
    112     For(i,0,n) read(a[i]);
    113     For(i,0,m) read(b[i]);
    114     div(a,n,b,m,D,R);
    115     For(i,0,n-m) printf("%lld ",D[i]); puts("");
    116     For(i,0,m-1) printf("%lld ",R[i]); puts("");
    117     Formylove;
    118 }
    多项式除法/取模

    一个通往之前写的题的传送门

     

    1.bzoj3527: [Zjoi2014]力

    $F_j= sum_{i<j} q_i/(i-j)^2-sum_{i>j}q_i/(i-j)^2$

    设$E1_j=sum_{i<j} q_i/(j-i)^2$

    $E2_j=sum_{i>j}q_i/(i-j)^2$

    $B={1/i^2},B_0=0$

    $F=E1-E2$

    $E1=q*B$

    $E2_j= sum_{i=j}^{n-1}q_i/(i-j)^2$

    设$qR_{i}=q_{n-i-1}$

    $E2_j= sum_{i=j}^{n-1}qR_{n-i-1}/(i-j)^2$

    $=sum_{j+k=n-j-1}qR_j*B_k$

    $=qR*B(n-j-1)$

     1 // luogu-judger-enable-o2
     2 //Achen
     3 #include<algorithm>
     4 #include<iostream>
     5 #include<cstring>
     6 #include<cstdlib>
     7 #include<vector>
     8 #include<cstdio>
     9 #include<queue>
    10 #include<cmath>
    11 #include<set>
    12 #include<map>
    13 #define pi acos(-1)
    14 #define Formylove return 0
    15 #define For(i,a,b) for(int i=(a);i<=(b);i++)
    16 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    17 const int N=2100007;
    18 typedef long long LL;
    19 typedef double db;
    20 using namespace std;
    21 int n,nn,m,l,rev[N];
    22 
    23 template<typename T>void read(T &x)  {
    24     char ch=getchar(); x=0; T f=1;
    25     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
    26     if(ch=='-') f=-1,ch=getchar();
    27     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    28 }
    29 
    30 struct E {
    31     db a,b;
    32     E(db a=0,db b=0):a(a),b(b){}
    33 }p[N],p_R[N],b[N];
    34 E operator + (const E&A,const E&B) { return E(A.a+B.a,A.b+B.b); }
    35 E operator - (const E&A,const E&B) { return E(A.a-B.a,A.b-B.b); } 
    36 E operator * (const E&A,const E&B) { return E(A.a*B.a-A.b*B.b,A.a*B.b+A.b*B.a); } 
    37 E operator / (const E&A,const int&B) { return E(A.a/B,A.b/B); } 
    38 
    39 void FFT(E a[],int n,int f) {
    40     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
    41     for(int i=1;i<n;i<<=1) {
    42         E wi(cos(pi/i),f*sin(pi/i));
    43         for(int j=0,pp=(i<<1);j<n;j+=pp) {
    44             E w(1,0);
    45             for(int k=0;k<i;k++,w=w*wi) {
    46                 E x=a[j+k],y=w*a[j+k+i];
    47                 a[j+k]=x+y; a[j+k+i]=x-y;
    48             }
    49         }
    50     }
    51     if(f==-1) For(i,0,n) a[i]=a[i]/n;
    52 } 
    53 
    54 int main() {
    55 #ifdef ANS
    56     freopen(".in","r",stdin);
    57     freopen(".out","w",stdout);
    58 #endif
    59     read(n);
    60     For(i,0,n-1) { scanf("%lf",&p[i].a); p_R[n-i-1]=p[i]; }
    61     For(i,1,n-1) b[i].a=1.0/i/i;
    62     for(nn=1;nn<=n+n+1;nn<<=1) l++;
    63     For(i,0,nn) rev[i]=((rev[i>>1]>>1)|((i&1)<<(l-1)));
    64     FFT(p,nn,1); FFT(p_R,nn,1); FFT(b,nn,1);
    65     For(i,0,nn) p[i]=p[i]*b[i],p_R[i]=p_R[i]*b[i];
    66     FFT(p,nn,-1); FFT(p_R,nn,-1);
    67     For(i,0,n-1) printf("%lf
    ",p[i].a-p_R[n-i-1].a);
    68     Formylove;
    69 }
    View Code

     

    2.bzoj3160: 万径人踪灭

    先把序列中间加上间隔符。

    把a的位置设为1,其余为0,自己卷自己得到第2*i项就是关于i对称的位置的数量。同理处理b。然后每个位置计算答案累加。

    不允许连续,再跑一遍马拉车减去答案即可。

      1 //Achen
      2 #include<algorithm>
      3 #include<iostream>
      4 #include<cstring>
      5 #include<cstdlib>
      6 #include<vector>
      7 #include<cstdio>
      8 #include<queue>
      9 #include<cmath>
     10 #include<set>
     11 #include<map>
     12 #define pi acos(-1)
     13 #define Formylove return 0
     14 #define For(i,a,b) for(int i=(a);i<=(b);i++)
     15 #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
     16 const int N=1e6+7,mod=1e9+7;
     17 typedef long long LL;
     18 typedef double db;
     19 using namespace std;
     20 int len,n,l,cnt,rev[N];
     21 char s[N],ss[N*2];
     22 LL ans;
     23 
     24 template<typename T>void read(T &x)  {
     25     char ch=getchar(); x=0; T f=1;
     26     while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
     27     if(ch=='-') f=-1,ch=getchar();
     28     for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
     29 }
     30 
     31 struct E {
     32     db a,b;
     33     E(db a=0,db b=0):a(a),b(b){}
     34 }a[N],b[N];
     35 E operator + (const E&A,const E&B) { return E(A.a+B.a,A.b+B.b); }
     36 E operator - (const E&A,const E&B) { return E(A.a-B.a,A.b-B.b); } 
     37 E operator * (const E&A,const E&B) { return E(A.a*B.a-A.b*B.b,A.a*B.b+A.b*B.a); } 
     38 E operator / (const E&A,const int&B) { return E(A.a/B,A.b/B); } 
     39 
     40 void FFT(E a[],int n,int f) {
     41     For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
     42     for(int i=1;i<n;i<<=1) {
     43         E wi(cos(pi/i),f*sin(pi/i));
     44         for(int j=0,pp=(i<<1);j<n;j+=pp) {
     45             E w(1,0);
     46             for(int k=0;k<i;k++,w=w*wi) {
     47                 E x=a[j+k],y=w*a[j+k+i];
     48                 a[j+k]=x+y; a[j+k+i]=x-y;
     49             }
     50         }
     51     }
     52     if(f==-1) For(i,0,n) a[i]=a[i]/n;
     53 } 
     54 
     55 LL ksm(LL a,LL b) {
     56     LL rs=1,bs=a;
     57     while(b) {
     58         if(b&1) rs=rs*bs%mod;
     59         bs=bs*bs%mod;
     60         b>>=1;
     61     }
     62     return rs;
     63 }
     64 
     65 int rad[N<<1];
     66 void manacher() {
     67     LL tot=len;
     68     for(int i=1,k=0,j;i<cnt;) {
     69         while(ss[i-k-1]==ss[i+k+1]) k++;
     70         rad[i]=k;
     71         for(j=1;j<=k&&rad[i-j]!=rad[i]-j;j++)
     72             rad[i+j]=min(rad[i]-j,rad[i-j]);
     73         i+=j;
     74         k=max(0,k-j);
     75     }
     76     For(i,0,cnt-1) tot=(tot+rad[i]/2)%mod;
     77     ans=(ans-tot+mod)%mod;
     78 }
     79 
     80 int main() {
     81 #ifdef ANS
     82     freopen(".in","r",stdin);
     83     freopen(".out","w",stdout);
     84 #endif
     85     scanf("%s",s);
     86     len=strlen(s);
     87     cnt=0;
     88     ss[cnt++]='&';
     89     ss[cnt++]='#';
     90     For(i,0,len-1) {
     91         ss[cnt++]=s[i];
     92         ss[cnt++]='#';
     93     } 
     94     ss[cnt++]='*';
     95     For(i,0,cnt-1) 
     96         if(ss[i]=='a') a[i].a=1;
     97         else if(ss[i]=='b') b[i].a=1;
     98     for(n=1;n<=cnt*2;n<<=1) l++;
     99     For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
    100     FFT(a,n,1); FFT(b,n,1);
    101     For(i,0,n) a[i]=a[i]*a[i],b[i]=b[i]*b[i];
    102     FFT(a,n,-1); FFT(b,n,-1);
    103     For(i,0,cnt-1) {
    104         LL c; 
    105         if(ss[i]=='a') c=((int)(a[2*i]+0.5).a+1)/2+((int)(b[2*i].a+0.5))/2;
    106         else if(ss[i]=='b') c=((int)(a[2*i].a+0.5))/2+((int)(b[2*i].a+0.5)+1)/2;
    107         else c=((int)(a[2*i].a+0.5))/2+((int)(b[2*i].a+0.5))/2;
    108         if(c>0) ans=(ans+ksm(2,c)-1+mod)%mod;
    109     }
    110     manacher();
    111     printf("%lld
    ",ans);
    112     Formylove;
    113 }
    114 /*
    115 abaabaa
    116 
    117 aaabbbaaa
    118 
    119 aaaaaaaa
    120 */
    View Code

     

    3.bzoj4259: 残缺的字符串

    把*值设为0,定义比较函数f(i,j)=(a[i]-b[j])^2*a[i]*b[j]。

    把b取反化成卷积,拆开后FFT即可。

    把(len1-1)/2写成len1/2交换的时候多换了一位,调了一晚上。。。然后又数组开大了一会T一会MLE。。。

    //Achen
    #include<algorithm>
    #include<iostream>
    #include<cstring>
    #include<cstdlib>
    #include<vector>
    #include<cstdio>
    #include<queue>
    #include<cmath>
    #include<set>
    #include<map>
    #define pi acos(-1)
    #define Formylove return 0
    #define For(i,a,b) for(int i=(a);i<=(b);i++)
    #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    const int N=1048576;
    typedef long long LL;
    typedef double db;
    using namespace std;
    int len1,len2,n,l,rev[N],q[N],a[N],b[N];
    char s1[N],s2[N];
    
    template<typename T>void read(T &x)  {
        char ch=getchar(); x=0; T f=1;
        while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
        if(ch=='-') f=-1,ch=getchar();
        for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    }
    
    struct E {
        db a,b;
        E(db a=0,db b=0):a(a),b(b){}
    }A[N],B[N],t[N];
    E operator + (const E&A,const E&B) { return E(A.a+B.a,A.b+B.b); }
    E operator - (const E&A,const E&B) { return E(A.a-B.a,A.b-B.b); } 
    E operator * (const E&A,const E&B) { return E(A.a*B.a-A.b*B.b,A.a*B.b+A.b*B.a); } 
    E operator / (const E&A,const int&B) { return E(A.a/B,A.b/B); } 
    
    void FFT(E a[],int n,int f) {
        For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
        for(int i=1;i<n;i<<=1) {
            E wi(cos(pi/i),f*sin(pi/i));
            for(int j=0,pp=(i<<1);j<n;j+=pp) {
                E w(1,0);
                for(int k=0;k<i;k++,w=w*wi) {
                    E x=a[j+k],y=w*a[j+k+i];
                    a[j+k]=x+y; a[j+k+i]=x-y;
                }
            }
        }
        if(f==-1) For(i,0,n) a[i]=a[i]/n;
    } 
    
    int main() {
        read(len1); read(len2);
        scanf("%s%s",s1,s2);
        For(i,0,(len1-1)/2) swap(s1[i],s1[len1-i-1]);
        For(i,0,len1-1) if(s1[i]!='*') b[i]=s1[i]-'a'+1;
        For(i,0,len2-1) if(s2[i]!='*') a[i]=s2[i]-'a'+1;
        for(n=1;n<=len1+len2;n<<=1) l++;
        For(i,1,n) rev[i]=((rev[i>>1]>>1)|((i&1)<<(l-1)));
        
        For(i,0,len1-1) B[i].a=b[i]*b[i]*b[i];
        For(i,0,len2-1) A[i].a=a[i];
        FFT(A,n,1),FFT(B,n,1); 
        For(i,0,n) t[i]=A[i]*B[i];
        
        For(i,0,len1-1) B[i]=E(b[i],0); For(i,len1,n) B[i]=E(0,0);
        For(i,0,len2-1) A[i]=E(a[i]*a[i]*a[i],0); For(i,len2,n) A[i]=E(0,0);
        FFT(A,n,1),FFT(B,n,1); For(i,0,n) t[i]=t[i]+A[i]*B[i];
        
        For(i,0,len1-1) B[i]=E(b[i]*b[i],0); For(i,len1,n) B[i]=E(0,0);
        For(i,0,len2-1) A[i]=E(a[i]*a[i],0); For(i,len2,n) A[i]=E(0,0);
        FFT(A,n,1),FFT(B,n,1); For(i,0,n) t[i]=t[i]-A[i]*B[i]*E(2,0);
        
        /*For(i,0,n) { 
            t[i]=A[2][i]*B[0][i]+A[0][i]*B[2][i]-A[1][i]*B[1][i]-A[1][i]*B[1][i];
        }*/
        FFT(t,n,-1);
        int ans=0;
        For(i,len1-1,len2-1) 
            if((int)(t[i].a+0.5)==0) ans++,q[++q[0]]=i-len1+2;
        printf("%d
    ",ans);
        if(q[0]>1) For(i,1,q[0]-1) printf("%d ",q[i]);
        if(q[0]) printf("%d",q[ans]);
        Formylove;
    }
    /*
    4 30
    d*ad
    *cbacbeeae*b**debabcecdb*aaacb
    */
    View Code

     

    4.hdu4609 3-idiots

    FFT求出和为i的木棍对的数量,枚举三角形最长的一边,把和比它大的木棍对的数量减去不合法的(两根都比它小为合法),再加上等边和等腰即可。

    因为清0和开LL又wa了好久。。。

    //Achen
    #include<algorithm>
    #include<iostream>
    #include<cstring>
    #include<cstdlib>
    #include<vector>
    #include<cstdio>
    #include<queue>
    #include<cmath>
    #include<set>
    #include<map>
    #define Formylove return 0
    #define For(i,a,b) for(int i=(a);i<=(b);i++)
    #define Rep(i,a,b) for(int i=(a);i>=(b);i--)
    typedef long long LL;
    typedef double db;
    const int N=262149;
    const db pi=acos(-1);
    using namespace std;
    int T,n,m,l,rev[N],L[N];
    LL ans,suml[N],cnt[N];
    
    template<typename T>void read(T &x)  {
        char ch=getchar(); x=0; T f=1;
        while(ch!='-'&&(ch<'0'||ch>'9')) ch=getchar();
        if(ch=='-') f=-1,ch=getchar();
        for(;ch>='0'&&ch<='9';ch=getchar()) x=x*10+ch-'0'; x*=f;
    }
    
    struct E {
        db a,b;
        E(db a=0,db b=0):a(a),b(b){}
    }a[N],b[N];
    E operator +(const E&A,const E&B) { return E(A.a+B.a,A.b+B.b); }
    E operator -(const E&A,const E&B) { return E(A.a-B.a,A.b-B.b); }
    E operator *(const E&A,const E&B) { return E(A.a*B.a-A.b*B.b,A.a*B.b+A.b*B.a); }
    E operator /(const E&A,const db&B) { return E(A.a/B,A.b/B); }
    
    void FFT(E a[],int n,int f) {
        For(i,0,n) if(i<rev[i]) swap(a[i],a[rev[i]]);
        for(int i=1;i<n;i<<=1) {
            E wi(cos(pi/i),f*sin(pi/i));
            for(int j=0,pp=(i<<1);j<n;j+=pp) {
                E w(1,0);
                for(int k=0;k<i;k++,w=w*wi) {
                    E x=a[j+k],y=w*a[j+k+i];
                    a[j+k]=x+y; a[j+k+i]=x-y;
                } 
            }
        }
        if(f==-1) For(i,0,n) a[i]=a[i]/n; 
    }
    
    LL C_2(int n) { return n<2?0:(LL)n*(n-1)/2;}
    LL C_3(int n) { return n<3?0:(LL)n*(n-1)*(n-2)/6; }
    
    //#define ANS
    int main() {
    #ifdef ANS
        freopen("std.in","r",stdin);
        //freopen(".out","w",stdout);
    #endif
        read(T);
        while(T--) {
            read(m); int mxlen=0;
            For(i,1,m) { read(L[i]); mxlen=max(mxlen,L[i]); }
            for(n=1,l=0;n<=mxlen+mxlen;n<<=1) l++;
            For(i,1,n) rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1)); 
            For(i,0,mxlen*2) cnt[i]=0;
            For(i,0,n) a[i]=E(0,0);
            For(i,1,m) a[L[i]].a++,cnt[L[i]]++;
            FFT(a,n,1);
            For(i,0,n) a[i]=a[i]*a[i];
            FFT(a,n,-1);
            For(i,0,n) a[i].a=(LL)(a[i].a+0.5);
            For(i,1,mxlen) a[i*2].a-=cnt[i];
            For(i,0,n) a[i].a/=2;
            For(i,1,mxlen) suml[i]=suml[i-1]+cnt[i];
            LL sum=0; ans=0;
            Rep(i,mxlen*2,1) {
                if(cnt[i]) {
                    LL tp=suml[mxlen]-suml[i];
                    ans+=(sum-tp*suml[i-1]-C_2(tp)-C_2(cnt[i])-(LL)cnt[i]*(m-cnt[i]))*cnt[i]+C_3(cnt[i])+C_2(cnt[i])*suml[i-1];
                }
                sum+=a[i].a;
            }
            db pt=(db)ans/C_3(m);
            printf("%.7lf
    ",pt); 
        }
        Formylove;
    }
    /*
    1
    7
    2 3 2 2 3 3 4
    */
    View Code

     

  • 相关阅读:
    【oracle】约束之非空约束
    【oracle】操作表数据之修改和删除
    【oracle】操作表中的数据之添加
    HashMap和HashTable的比较
    java集合类-总纲
    性能测试工具篇(开源&商业)
    比较好的idea工具介绍
    根据mysql某一条记录生成对应实体类的初始化
    Spring+mybatis多数据源切换笔记
    Jmeter的Throughput和平均响应时间计算方法整理
  • 原文地址:https://www.cnblogs.com/Achenchen/p/9456439.html
Copyright © 2020-2023  润新知