題面

我居然不會寫多項式卷積了，把$a[i]=a[i]*b[i]$寫成了$a[i].x=a[i].x*b[i].x$還調了好長時間，有毒

把式子寫成$C^2(x)-A^2(x)-B^2(x)=0$的形式，然後牛頓叠代即可

  1 #include<cmath>
  2 #include<cstdio>
  3 #include<cctype>
  4 #include<cstring>
  5 #include<algorithm>
  6 using namespace std;
  7 const int N=800005;
  8 const
 
 double pai=acos(-1),inf=1e9,eps=1e-6;
  9 struct cpx
 10 {
 11     double x,y;
 12 }a[N],b[N],c[N];
 13 cpx operator + (cpx a,cpx b)
 14 {
 15     return (cpx){a.x+b.x,a.y+b.y};
 16 }
 17 cpx operator - (cpx a,cpx b)
 18 {
 19     return (cpx){a.x-b.x,a.y-b.y};
 20 }
 21 cpx operator * (cpx a,cpx b)
 
 22 {
 23     double x1=a.x,y1=a.y,x2=b.x,y2=b.y;
 24     return (cpx){x1*x2-y1*y2,x1*y2+y1*x2};
 25 }
 26 int la,lb,lc,n,m,nm,rev[N]; 
 27 double ll,rr,Sin[32],Cos[32],f[N],d[N];
 28 void Read(double &x)
 29 {
 30     int ret=0; x=0;
 31     char ch=getchar();
 32     while(!isdigit(ch))
 33         ch=getchar();
 
 34     while(isdigit(ch))
 35         ret=(ret<<1)+(ret<<3)+(ch^48),ch=getchar();
 36     x=ret; return ;
 37 }
 38 void Prework()
 39 {
 40     m=2*(max(la,max(lb,lc))+1),n=1; while(n<=m) n<<=1;
 41     for(int i=0;i<=n;i++) 
 42         rev[i]=(rev[i>>1]>>1)+(i&1)*(n>>1);
 43     for(int i=1;i<=24;i++)
 44         Sin[i]=sin(2*pai/(1<<i)),Cos[i]=cos(2*pai/(1<<i));
 45 }
 46 void Trans(cpx *cop,int len,int typ)
 47 {
 48     register int i,j,k;
 49     for(i=0;i<=len;i++) 
 50         if(rev[i]>i) swap(cop[rev[i]],cop[i]);
 51     for(i=2;i<=len;i<<=1)
 52     {
 53         int lth=i>>1,lgg=log2(i);
 54         cpx omg=(cpx){Cos[lgg],typ*Sin[lgg]};
 55         for(j=0;j<len;j+=i)
 56         {
 57             cpx ori={1,0};
 58             for(k=j;k<j+lth;k++,ori=ori*omg)
 59             {
 60                 cpx tmp=ori*cop[k+lth];
 61                 cop[k+lth]=cop[k]-tmp,cop[k]=cop[k]+tmp;
 62             }
 63         }
 64     }
 65     if(typ==-1)
 66         for(int i=0;i<=len;i++) cop[i].x/=len;
 67 }
 68 void Square(cpx *cop,int len)
 69 {
 70     register int i;
 71     Trans(cop,len,1);
 72     for(i=0;i<=len;i++) 
 73         cop[i]=cop[i]*cop[i];
 74     Trans(cop,len,-1);
 75 }
 76 double Ask(double *f,double x,int len)
 77 {
 78     double ret=0,var=1;
 79     for(int i=0;i<=len;i++)
 80         ret+=f[i]*var,var*=x;
 81     return ret;
 82 }
 83 double Iterate(double x)
 84 {
 85     for(int i=1;i<=30;i++)
 86     {
 87         double y=Ask(f,x,n);
 88         if(fabs(y)<=eps) return x;
 89         x-=y/Ask(d,x,nm),x=max(x,ll),x=min(x,rr);
 90     }
 91     return inf;
 92 }
 93 int main()
 94 {
 95     register int i;
 96     scanf("%d%d%d%lf%lf",&la,&lb,&lc,&ll,&rr);
 97     for(i=0;i<=la;i++) Read(a[i].x);
 98     for(i=0;i<=lb;i++) Read(b[i].x);
 99     for(i=0;i<=lc;i++) Read(c[i].x);
100     Prework(),Square(a,n),Square(b,n),Square(c,n);
101     for(i=0;i<=n;i++) c[i]=c[i]-a[i]-b[i]; nm=n-1;
102     for(i=0;i<=n;i++) f[i]=c[i].x,d[i]=(i+1)*c[i+1].x; 
103     double ans=Iterate((ll+rr)/2);
104     ans>=inf?printf("Inconsistent!"):printf("%.8f",ans);
105     return 0;
106 }

解題：洛谷 4986 逃離