1. 程式人生 > >【BZOJ1502】[NOI2005]月下檸檬樹 Simpson積分

【BZOJ1502】[NOI2005]月下檸檬樹 Simpson積分

整數 容易 就是 分享 out width inline pri ret

【BZOJ1502】[NOI2005]月下檸檬樹

Description

李哲非常非常喜歡檸檬樹,特別是在靜靜的夜晚,當天空中有一彎明月溫柔地照亮地面上的景物時,他必會悠閑地坐在他親手植下的那棵檸檬樹旁,獨自思索著人生的哲理。李哲是一個喜愛思考的孩子,當他看到在月光的照射下檸檬樹投在地面上的影子是如此的清晰,馬上想到了一個問題:樹影的面積是多大呢?李哲知道,直接測量面積是很難的,他想用幾何的方法算,因為他對這棵檸檬樹的形狀了解得非常清楚,而且想好了簡化的方法。李哲將整棵檸檬樹分成了n 層,由下向上依次將層編號為1,2,…,n。從第1到n-1 層,每層都是一個圓臺型,第n 層(最上面一層)是圓錐型。對於圓臺型,其上下底面都是水平的圓。對於相鄰的兩個圓臺,上層的下底面和下層的上底面重合。第n 層(最上面一層)圓錐的底面就是第n-1 層圓臺的上底面。所有的底面的圓心(包括樹頂)處在同一條與地面垂直的直線上。李哲知道每一層的高度為h1,h2,…,hn,第1 層圓臺的下底面距地面的高度為h0,以及每層的下底面的圓的半徑r1,r2,…,rn。李哲用熟知的方法測出了月亮的光線與地面的夾角為alpha。 技術分享圖片
為了便於計算,假設月亮的光線是平行光,且地面是水平的,在計算時忽略樹幹所產生的影子。李哲當然會算了,但是他希望你也來練練手

Input

第1行包含一個整數n和一個實數alpha,表示檸檬樹的層數和月亮的光線與地面夾角(單位為弧度)。 第2行包含n+1個實數h0,h1,h2,…,hn,表示樹離地的高度和每層的高度。 第3行包含n個實數r1,r2,…,rn,表示檸檬樹每層下底面的圓的半徑。 上述輸入文件中的數據,同一行相鄰的兩個數之間用一個空格分隔。 輸入的所有實數的小數點後可能包含1至10位有效數字。 1≤n≤500,0.3<alpha<π/2,0<hi≤100,0<ri≤100

Output

輸出1個實數,表示樹影的面積。四舍五入保留兩位小數。

Sample Input

2 0.7853981633
10.0 10.00 10.00
4.00 5.00

Sample Output

171.97

題解:簡潔題意就是讓你求一堆圓和梯形的面積交。

Simpson積分: 技術分享圖片

相當於用一個3次函數去擬合所求的圖形,可以用於任意連續不規則圖形,但是誤差很大。自適應simpson積分呢,就是再對f(l,mid)和f(mid,r)分別算一下。如果f(l,mid)+f(mid,r)與f(l,r)誤差很小,則直接返回f(l,r),否則繼續遞歸計算。這樣誤差就比較小了(雖說也可以卡)。

剩下的問題就是如何求兩圓的公切線。比較容易的方法是設兩圓半徑為R,r,先令R‘=R-r,r‘=0,這樣就把第二個圓縮成了一個點,變成求點與圓的切線,再平移回去即可。

#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <iostream>
#define pi acos(-1.0)
using namespace std;
typedef double db;
const db eps=1e-6;
int n,m;
db alpha,h[510],O[510],R[510],ax[510],ay[510],bx[510],by[510];
inline db f(db x)
{
	db ret=0;
	int i;
	for(i=1;i<=n;i++)	if(x>=O[i]-R[i]&&x<=O[i]+R[i])
		ret=max(ret,sqrt(R[i]*R[i]-(O[i]-x)*(O[i]-x)));
	for(i=1;i<=m;i++)	if(x>=ax[i]&&x<=bx[i])
		ret=max(ret,ay[i]+(by[i]-ay[i])*(x-ax[i])/(bx[i]-ax[i]));
	return ret;
}
inline db simpson(db a,db b)
{
	return (b-a)/6*(f(a)+f(b)+f((a+b)/2)*4);
}
inline db calc(db l,db r,db val)
{
	db mid=(l+r)/2,a=simpson(l,mid),b=simpson(mid,r);
	if(fabs(a+b-val)<eps)	return val;
	return calc(l,mid,a)+calc(mid,r,b);
}
int main()
{
	scanf("%d%lf",&n,&alpha);
	int i;
	db l=1e9,r=-1e9;
	for(i=0;i<=n;i++)
	{
		scanf("%lf",&h[i]),h[i]/=tan(alpha);
		if(i)	h[i]+=h[i-1];
	}
	for(i=1;i<=n;i++)
	{
		scanf("%lf",&R[i]),O[i]=h[i-1],l=min(l,O[i]-R[i]),r=max(r,O[i]+R[i]);
		if(i!=1&&O[i]-O[i-1]>fabs(R[i]-R[i-1]))
		{
			db a=(R[i-1]-R[i])/(O[i]-O[i-1]),b=sqrt(1-a*a);
			ax[++m]=O[i-1]+a*R[i-1],ay[m]=b*R[i-1];
			bx[m]=O[i]+a*R[i],by[m]=b*R[i];
		}
	}
	if(h[n]>O[n]+R[n])
	{
		db a=R[n]/(h[n]-O[n]),b=sqrt(1-a*a);
		r=h[n];
		ax[++m]=O[n]+a*R[n],ay[m]=b*R[n];
		bx[m]=h[n],by[m]=0;
	}
	printf("%.2lf",calc(l,r,simpson(l,r))*2);
	return 0;
}//2 0.7853981633 10.0 10.00 10.00 4.00 5.00

【BZOJ1502】[NOI2005]月下檸檬樹 Simpson積分