BZOJ 1502 月下柠檬树 simpson积分

#include <iostream>
#include <cstring>
#include <cstdio>
#include <cstdlib>
#include <algorithm>
#include <cmath>

#define N 520
#define EPS 1e-6

using namespace std;

int n;
double a
,b
,af,h,lt,rt;

struct C
{
double x,y,p,q;
}c
;

inline int dc(double x)
{
if(x>EPS) return 1;
else if(x<-EPS) return -1;
return 0;
}

inline void read()
{
scanf("%d%lf",&n,&af);
af=1/tan(af);
for(int i=1;i<=n+1;i++)
{
scanf("%lf",&a[i]);
h+=a[i]; a[i]=h*af;
}
lt=a[1];rt=a[n+1];
for(int i=1;i<=n;i++)
{
scanf("%lf",&b[i]);
lt=min(lt,a[i]-b[i]);
rt=max(rt,a[i]+b[i]);
}
}

inline void calc()
{
for(int i=1;i<=n;i++)
if(a[i+1]-a[i]>fabs(b[i+1]-b[i]))
{
c[i].x=a[i]+b[i]*(b[i]-b[i+1])/(a[i+1]-a[i]);
c[i].y=sqrt(b[i]*b[i]-(c[i].x-a[i])*(c[i].x-a[i]));
c[i].p=a[i+1]+b[i+1]*(b[i]-b[i+1])/(a[i+1]-a[i]);
c[i].q=sqrt(b[i+1]*b[i+1]-(c[i].p-a[i+1])*(c[i].p-a[i+1]));
}
}

inline double f(double p)
{
double res=0;
for(int i=1;i<=n;i++)
{
if(fabs(a[i]-p)<b[i]) res=max(res,sqrt(b[i]*b[i]-(a[i]-p)*(a[i]-p)));
if(p>c[i].x&&p<c[i].p) res=max(res,c[i].y-(p-c[i].x)*(c[i].y-c[i].q)/(c[i].p-c[i].x));
}
return res;
}

inline double simpson(double l,double r,double fl,double fmid,double fr)
{
return (fl+fr+4*fmid)*(r-l)/6;
}

inline double rsimpson(double l,double r,double fl,double fmid,double fr)
{
double p,q,mid,x,y,z;
mid=(l+r)/2;
p=f((l+mid)/2); q=f((mid+r)/2);
x=simpson(l,r,fl,fmid,fr); y=simpson(l,mid,fl,p,fmid); z=simpson(mid,r,fmid,q,fr);
if(dc(fabs(x-y-z))==0) return y+z;
else return rsimpson(l,mid,fl,p,fmid)+rsimpson(mid,r,fmid,q,fr);
}

inline void go()
{
calc();
printf("%.2lf\n",2*rsimpson(lt,rt,0,f(lt+rt)/2,0));
}

int main()
{
read(),go();
return 0;
}