【bzoj1502】 NOI2005—月下柠檬树

http://www.lydsy.com/JudgeOnline/problem.php?id=1502 (题目链接)

今天考试题，从来没写过圆的面积之类的东西。。GG

题意：一颗树由n个圆台组成，现在有倾斜角为alpha的光，不计树干阴影，光线沿直线传播，求这个树在水平地面投影的面积。

Solution
　　Simpson积分法。
　　直接蒯一份题解算了（http://blog.csdn.net/ww140142/article/details/48296273）

　　给出一颗树。。求这棵树的阴影面积。。。
　　什么鬼题！= =
　　当然这是一道计算几何；
　　首先我们把树的尖看成半径为0的圆；
　　三维图形到二维的投影貌似很难啊，不过这题也有很多特殊性；
　　因为圆对于平行的面的投影都是相同的正圆，所以可以直接投在地面上；
　　而圆心距则是除了一个tan函数的样子，所以圆已经被扔到待求的平面上了；
　　这些圆的公切线！(这样就从母的变成公的啦233)
　　公切线的求法就是画图上勾股定理相似乱搞，调两组数据改改就好了；
　　然后正解似乎是讨论了一堆东西然后分别求面积；
　　然后我们似乎可以直接上Simpson积分！
　　具体细节不说了，代码里都有；
　　这题精度不是特别卡，EPS=1e5就够了？

　　当然，我写的是1e7= =。

代码：

// bzoj1502
#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdlib>
#include<cstdio>
#include<cmath>
#include<vector>
#define eps 1e-7
#define inf 2147483640
#define LL long long
#define free(a) freopen(a".in","r",stdin);freopen(a".out","w",stdout);
using namespace std;
inline LL getint() {
    LL x=0,f=1;char ch=getchar();
    while (ch>'9' || ch<'0') {if (ch=='-') f=-1;ch=getchar();}
    while (ch>='0' && ch<='9') {x=x*10+ch-'0';ch=getchar();}
    return x*f;
}

const int maxn=510;
struct S {double x,y,p,q;}c[maxn];
double h[maxn],s[maxn],r[maxn],alpha,ll,rr;
int n,size=0;

double sqr(double x) {return x*x;}
double f(double l) {
    double t=0.0;
    for (int i=1;i<=n;i++)
        if (fabs(s[i]-l)<r[i]) t=max(t,sqrt(sqr(r[i])-sqr(s[i]-l)));
    for (int i=1;i<=size;i++)
        if (c[i].x<l && l<c[i].p) t=max(t,c[i].y+(c[i].q-c[i].y)*(l-c[i].x)/(c[i].p-c[i].x));
    return t;
}
double Simpson(double l,double r,double fl,double fmid,double fr) {
    double m=(l+r)/2;
    double p=f((l+m)/2),q=f((m+r)/2);
    double x=(fl+4*fmid+fr)*(r-l)/6,y=(fl+4*p+fmid)*(m-l)/6,z=(fmid+4*q+fr)*(r-m)/6;
    if (fabs(x-y-z)<eps) return y+z;
    return Simpson(l,m,fl,p,fmid)+Simpson(m,r,fmid,q,fr);
}
int main() {
    scanf("%d%lf",&n,&alpha);
    alpha=1.0/tan(alpha);
    for (int i=1;i<=n+1;i++) {
        scanf("%lf",&h[i]);
        h[i]+=h[i-1];
        s[i]=h[i]*alpha;
    }
    ll=0;rr=0;
    for (int i=1;i<=n;i++) {
        scanf("%lf",&r[i]);
        ll=min(ll,s[i]-r[i]);
        rr=max(rr,s[i]+r[i]);
    }
    r[n+1]=0;
    for (int i=1;i<=n;i++) {
        double d=s[i+1]-s[i];
        if (d>fabs(r[i]-r[i+1])) {
            c[++size].x=s[i]-r[i]*(r[i+1]-r[i])/d;
            c[size].y=sqrt(sqr(r[i])-sqr(c[size].x-s[i]));
            c[size].p=s[i+1]-r[i+1]*(r[i+1]-r[i])/d;
            c[size].q=sqrt(sqr(r[i+1])-sqr(c[size].p-s[i+1]));
        }
    }
    rr=max(rr,s[n+1]);
    printf("%.2lf",2*Simpson(ll,rr,0,f((ll+rr)/2),0));
    return 0;
}