关键点是，水平的圆投影到水平面之后仍然是与原先全等的圆。
然后圆与圆之间通过曲面无缝连接，所以投影下来之后圆与圆之间通过公切线连接。
直接求有点难。把投影区域的上边界当成一个函数，然后套Simpson积分就简单了。
我的实现用了指针从左往右扫，其实没必要，暴力好像还会快一点。
代码(darkbzoj 218ms)：
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <unordered_set>
#include <queue>
using namespace std;
#define DEBUG 0
#define MAXN 511
typedef pair<double, int> pdi;
const double eps = 1e-8;
int sgn(double x) {
return x < -eps ? -1 : (x > eps ? 1 : 0);
}
const double pi = acos(-1);
double cot(double x) {
return tan(pi / 2 - x);
}
template <typename T>
T sq(T x) {
return x * x;
}
struct Func {
bool type; //0: line, 1: circle
double l, r;
double c, R;
bool operator < (const Func& b) const {
return l != b.l ? l < b.l : (r < b.r);
}
double y(double x) const {
if (type) {
return sqrt(sq(R) - sq(x - c));
} else {
return c + (R - c) / (r - l) * (x - l);
}
}
};
#define fsnow const Func fs[], const unordered_set<int>& now
double Simpson(fsnow, double l, double r) {
auto f = [&](double x) {
double y = numeric_limits<double>::min();
for (int i : now) {
y = max(y, fs[i].y(x));
}
return y;
};
return (f(l) + 4 * f((l+r)/2) + f(r)) * (r - l) / 6;
};
double Integrate(fsnow, double l, double r, double eps, double pre) {
double mid = (l+r) / 2;
double L = Simpson(fs, now, l, mid), R = Simpson(fs, now, mid, r), delta = L + R - pre;
return fabs(delta) <= 15*eps ? L + R + delta/15 :
Integrate(fs, now, l, mid, eps/2, L) + Integrate(fs, now, mid, r, eps/2, R); //If you TLE, then try not divide 2, i.e. "eps/2" -> "eps"
}
double Integrate(fsnow, double l, double r, double eps) {
return Integrate(fs, now, l, r, eps, Simpson(fs, now, l, r));
}
double Integrate(Func fs[], int n, double eps) {
struct REC {
double r;
int i;
bool operator < (const REC& b) const {
return r > b.r;
}
};
sort(fs, fs + n);
#if DEBUG
for (int i = 0; i < n; ++i) {
printf("%d %f %f %f %f\n", fs[i].type, fs[i].l, fs[i].r, fs[i].c, fs[i].R);
}
#endif
unordered_set<int> now;
priority_queue<REC> right;
double l = fs[0].l, ans = 0;
for (int i = 0; i < n; ++i) {
while (!right.empty() && right.top().r < fs[i].l) {
ans += Integrate(fs, now, l, right.top().r, eps);
l = right.top().r;
now.erase(right.top().i);
right.pop();
}
ans += Integrate(fs, now, l, fs[i].l, eps);
right.push(REC{fs[i].r, i});
now.insert(i);
l = fs[i].l;
}
while (!right.empty()) {
ans += Integrate(fs, now, l, right.top().r, eps);
l = right.top().r;
now.erase(right.top().i);
right.pop();
}
return ans;
}
#undef fsnow
int main() {
int n;
double a;
static double x[MAXN];
static double r[MAXN];
static Func fs[MAXN << 1];
scanf("%d%lf", &n, &a);
for (int i = 0; i <= n; ++i) {
scanf("%lf", x + i);
}
for (int i = 1; i <= n; ++i) {
scanf("%lf", r + i);
}
x[0] = 0;
for (int i = 1; i <= n; ++i) {
x[i] = x[i-1] + x[i] * cot(a);
}
#if DEBUG
printf("%f\n", cot(a));
for (int i = 0; i <= n; ++i) {
printf("%f ", x[i]);
}
putchar('\n');
#endif
int tot = 0;
for (int i = 1; i <= n; ++i) {
fs[tot++] = Func{true, x[i-1] - r[i], x[i-1] + r[i], x[i-1], r[i]};
}
r[n+1] = 0;
for (int i = 1; i <= n; ++i) {
double s = (r[i+1] - r[i]) / (x[i] - x[i-1]);
double d1 = r[i] * s, d2 = r[i+1] * s;
fs[tot++] = Func{false, x[i-1] - d1, x[i] - d2, sqrt(sq(r[i]) - sq(d1)), sqrt(sq(r[i+1]) - sq(d2))};
}
printf("%.2f\n", Integrate(fs, tot, 1e-6) * 2.0);
return 0;
}