【数值积分】 ZOJ 3832 Tilt Cylinder

转化成对矩形的积分。。。然后用辛普森积分就行了。。。注意控制精度。。。

#include <iostream>
#include <queue> 
#include <stack> 
#include <map> 
#include <set> 
#include <bitset> 
#include <cstdio> 
#include <algorithm> 
#include <cstring> 
#include <climits>
#include <cstdlib>
#include <cmath>
#include <time.h>
#define maxn 205
#define maxm 2000005
//#define eps 1e-10
#define mod 1000000007
#define INF 0x3f3f3f3f
#define PI (acos(-1.0))
#define lowbit(x) (x&(-x))
#define mp make_pair
#define ls o<<1
#define rs o<<1 | 1
#define lson o<<1, L, mid 
#define rson o<<1 | 1, mid+1, R
#define pii pair<int, int>
//#pragma comment(linker, "/STACK:16777216")
typedef long long LL;
typedef unsigned long long ULL;
//typedef int LL;
using namespace std;
LL qpow(LL a, LL b){LL res=1,base=a;while(b){if(b%2)res=res*base;base=base*base;b/=2;}return res;}
LL powmod(LL a, LL b){LL res=1,base=a;while(b){if(b%2)res=res*base%mod;base=base*base%mod;b/=2;}return res;}
// head

double R, H, h, a;

double f(double x)
{
	double t1 = 2.0 * sqrt(2 * R * x - x * x);
	double t2;
	if(fabs(a) < 1e-7) {
		if(x < h) t2 = H;
		else t2 = 0;
	}
	else t2 = h / sin(a) + x / tan(PI - a);
	t2 = min(H, t2);
	t2 = max(0.0, t2);
	return t1 * t2;
}

double simpson(double a, double b)
{
	double c = (a + b) / 2;
	return (f(a) + 4 * f(c) + f(b)) * (b - a) / 6;
}

double calc(double a, double b, double eps, double res)
{
	double c = (a + b) / 2;
	double L = simpson(a, c), R = simpson(c, b);
	if(fabs(R + L - res) < 15.0 * eps) return R + L + (R + L - res) / 15.0;
	return calc(a, c, eps / 2, L) + calc(c, b, eps / 2, R);
}
double calc(double a, double b, double eps)
{
	return calc(a, b, eps, simpson(a, b));
}

void work(void)
{
	if(fabs(90 - a) < 1e-7) {
		printf("%.9f\n", PI * R * R * h);
		return;
	}
	a = a / 180 * PI;
	double ans = calc(0.000001, 2 * R - 0.000001, 1e-5);
	printf("%.18f\n", ans);
}

int main(void)
{
	while(scanf("%lf%lf%lf%lf", &R, &H, &h, &a)!=EOF) {
		work();
	}

	return 0;
}