LA 3485 自适应辛普森公式

#include <cstdio>
#include <cmath>
double a;
double asr(double a, double b, double eps, double A);
double asr(double a, double b, double eps);
double F(double x)
{
	return sqrt(1 + 4 * a * a * x * x);
}
double parabola_arc_length(double w, double h)
{
	a = 4.0 * h / (w * w); // 修改全局变量a，从而改变全局函数F的行为
	return asr(0, w / 2, 1e-5) * 2;
}
double simpson(double a, double b)
{
	double c = a + (b - a) / 2;
	return (F(a) + 4 * F(c) + F(b)) * (b - a) / 6;
}
double asr(double a, double b, double eps, double A)
{
	double c = a + (b - a) / 2;
	double L = simpson(a, c), R = simpson(c, b);
	if (fabs(L + R - A) <= 15 * eps) return L + R + (L + R - A) / 15.0;
	return asr(a, c, eps / 2, L) + asr(c, b, eps / 2, R);
}
double asr(double a, double b, double eps)
{
	return asr(a, b, eps, simpson(a, b));
}

int main(int argc, char const *argv[])
{
	int T;
	scanf("%d", &T);
	for (int kase = 1; kase <= T; kase++)
	{
		int D, H, B, L;
		scanf("%d%d%d%d", &D, &H, &B, &L);
		int n = (B + D - 1) / D;
		double D1 = (double) B / n;
		double L1 = (double) L / n;
		double x = 0, y = H;
		while (y - x > 1e-5)
		{
			double m = x + (y - x) / 2;
			if (parabola_arc_length(D1, m) < L1) x = m; else y = m;
		}
		if (kase > 1) printf("\n");
		printf("Case %d:\n%.2lf\n", kase, H - x);
	}
	return 0;
}

首先介绍下辛普森公式（没法帖图片，只能口述），对于一段函数，我们要求面积，除了求反导之外，有一种近似方法，即用辛普森公式。
a,b分别指函数最左端和最右端，这个公式推导可以去维基上看，他有个直观的应用便是求拟柱体（比如棱柱）体积，f函数便代表截面面积，套一套求体积的公式或许会好理解一点。
仅仅是套公式是不行的，或许这三个点均是特殊点（比如极值），误差那不是一般的大。怎么办呢？我们可是学信息的，一边不行多套几边公式计算机还是办得到的。
下面介绍自适应辛普森
将一段区间二分，再将两端区间求出的面积合并，如此递归下去。为了保证误差足够小，如果两端用上述公式算出面积之和近似等于该区间用上述公式求出面积，便不再递归，返回值。
如此一来，即缩小误差，又节省了在函数图象较平整区域的计算时间。