sgu217：Two Cylinders(simpson积分)

题目大意：
      给两个半径分别为R1,R2的无限长的圆柱，满足它们的轴相交且垂直，求公共体积。

分析：
      三维的simpson积分直接上，注意一下精度。

AC code：

#include <cstdio>#include <cmath>#include <algorithm>#define ONLINE_JUDGEtypedef double DB;using namespace std;const DB eps = 1e-10;DB r1, r2;DB simpson(DB l, DB r, DB fl, DB fm, DB fr){    return (r-l)*(fl+fm*4+fr)/6;}DB sqr(DB x) {return x*x;}DB f2(DB x, DB y){    if(sqr(x)+sqr(y)-sqr(r1) >= eps) return 0;    return sqrt(sqr(r2)-sqr(x))*2;}DB sol2(DB l, DB m, DB r, DB p, DB fl, DB fm, DB fr, DB ans){    DB m1 = (l+m)*0.5, m2 = (m+r)*0.5, fm1 = f2(m1, p), fm2 = f2(m2, p);    DB left = simpson(l, m, fl, fm1, fm), right = simpson(m, r, fm, fm2, fr);    return fabs(left+right-ans)<=eps ? ans : sol2(l, m1, m, p, fl, fm1, fm, left)+sol2(m, m2, r, p, fm, fm2, fr, right);}DB f1(DB p){    DB fl = f2(-r2, p), fm = f2(0, p), fr = f2(r2, p);    return sol2(-r2, 0, r2, p, fl, fm, fr, simpson(-r2, r2, fl, fm, fr));}DB sol1(DB l, DB m, DB r, DB fl, DB fm, DB fr, DB ans){    DB m1 = (l+m)*0.5, m2 = (m+r)*0.5, fm1 = f1(m1), fm2 = f1(m2);    DB left = simpson(l, m, fl, fm1, fm), right = simpson(m, r, fm, fm2, fr);    return fabs(left+right-ans)<=eps ? ans : sol1(l, m1, m, fl, fm1, fm, left)+sol1(m, m2, r, fm, fm2, fr, right);}int main(){    #ifndef ONLINE_JUDGE    freopen("sgu217.in", "r", stdin);    freopen("sgu217.out", "w", stdout);    #endif    scanf("%lf%lf", &r1, &r2);    if(r1 > r2) swap(r1, r2);       DB fl = f1(-r1), fm = f1(0), fr = f1(r1);    printf("%.4lf", sol1(-r1, 0, r1, fl, fm, fr, simpson(-r1, r1, fl, fm, fr)));    #ifndef ONLINE_JUDGE    fclose(stdin);    fclose(stdout);    #endif    return 0;}