bzoj 4909 [Sdoi2017]龙与地下城

题面

https://www.lydsy.com/JudgeOnline/problem.php?id=4909

题解

目前为止仅仅在LOJ上A掉这道题(Loj真快!)

当然不是标准做法

显然我们只要求一个

然后$a^n$的系数就表示选n个的方案数

那么我们找到

然后$a^n$的系数就表示选n个的概率

FFT即可

按理说这东西只能过60分但是LOJ的评测机成功过掉...而且时限4秒最慢一个点只用3秒!!!

Code

  1 #include<bits/stdc++.h>
  2 using namespace std;
  3 typedef long long ll;
  4 
  5 ll read(){
  6     ll x=0,f=1;char c=getchar();
  7     while(c<'0' || c>'9'){if(c=='-')f=-1;c=getchar();}
  8     while(c>='0' && c<='9'){x=x*10+c-'0';c=getchar();}
  9     return x*f;
 10 }
 11 const int maxn=8000200;
 12 struct Complex{
 13     double re,im;
 14     Complex(){
 15         re=im=0;
 16     }
 17     Complex(double alpha){
 18         re=cos(alpha);
 19         im=sin(alpha);
 20     }
 21     Complex(double _re,double _im){
 22         re=_re;
 23         im=_im;
 24     }
 25     Complex operator + (const Complex &x){
 26         return Complex(re+x.re,im+x.im);
 27     }
 28     Complex operator - (const Complex &x){
 29         return Complex(re-x.re,im-x.im);
 30     }
 31     Complex operator * (const Complex &x){
 32         return Complex(re*x.re-im*x.im,re*x.im+im*x.re);
 33     }
 34     Complex operator += (const Complex &x){
 35         return *this=*this+x;
 36     }
 37     Complex operator *= (const Complex &x){
 38         return *this=*this*x;
 39     }
 40 } A[maxn],B[maxn];
 41 
 42 int fft_lst,poly_rev[maxn];
 43 inline void fft_init(int n){
 44     if(fft_lst==n) return;
 45     fft_lst=n;
 46     for(int i=1,j=n>>1;i+1<n;i++){
 47         poly_rev[i]=j;
 48         int k=n>>1;
 49         while(j>=k){
 50             j-=k;
 51             k>>=1;
 52         }
 53         j+=k;
 54     }
 55 }
 56 
 57 inline void poly_fft(Complex *a,int len,bool f){
 58     fft_init(len);
 59     for(int i=1;i+1<len;i++)
 60         if(i<poly_rev[i]) swap(a[i],a[poly_rev[i]]);
 61     for(int i=1;i<len;i<<=1){
 62         Complex off((f ? -acos(-1.0) : acos(-1.0))/i);
 63         for(int j=0;j<len;j+=i<<1){
 64             Complex cur(0);
 65             for(int k=j;k<j+i;k++,cur*=off){
 66                 Complex x=a[k+i]*cur;
 67                 a[k+i]=a[k]-x;
 68                 a[k]+=x;
 69             }
 70         }
 71     }
 72     if(f){
 73         for(int i=0;i<len;i++)
 74             a[i].re/=len;
 75     }
 76 }
 77 
 78 int tc;
 79 
 80 int main(){
 81 #ifdef LZT
 82     freopen("in","r",stdin);
 83 #endif
 84     tc=read();
 85     while(tc--){
 86         int x=read(),y=read();
 87         int mx=x*y,len=1;
 88         while(len<mx) len<<=1;
 89         for(int i=0;i<len;i++)
 90             A[i]=B[i]=Complex();
 91         for(int i=0;i<x;i++)
 92             A[i].re=1.0/x;
 93         poly_fft(A,len,false);
 94         for(int i=0;i<len;i++){
 95             int nw=y;
 96             B[i]=Complex(1,0);
 97             while(nw){
 98                 if(nw&1) B[i]*=A[i];
 99                 A[i]*=A[i];
100                 nw>>=1;
101             }
102         }
103         poly_fft(B,len,true);
104         for(int i=1;i<len;i++)
105             B[i].re+=B[i-1].re;
106         for(int i=1;i<=10;i++){
107             int l=read(),r=read();
108             printf("%.8lf
",B[r].re-(l?B[l-1].re:0));
109         }
110     }
111     return 0;
112 }
View Code
原文地址:https://www.cnblogs.com/wawawa8/p/9374853.html