LA 4642 马尔法蒂问题 (二分或推公式)

最新推荐文章于 2020-11-05 14:40:01 发布

Tczxw

最新推荐文章于 2020-11-05 14:40:01 发布

阅读量1.6k

点赞数

CC 4.0 BY-SA版权

本文链接：https://blog.youkuaiyun.com/Tczxw/article/details/50807551

文章详细介绍了如何使用二分法解决几何问题，特别针对马尔法蒂问题给出了公式推导，并最终通过公式快速解决。文章还提到了在几何基础部分的二分法技巧提升。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

#include <bits/stdc++.h>
using namespace std;
struct Point
{
	double x, y;
	Point(double x = 0, double y = 0): x(x), y(y) {}
};
typedef Point Vector;
typedef vector<Point> Polygon;
Vector operator +(Vector A, Vector B)//
{
	return Vector(A.x + B.x, A.y + B.y);
}
Vector operator -(Point A, Point B)//
{
	return Vector(A.x - B.x , A.y - B.y);
}
Vector operator *(Vector A, double p)//
{
	return Vector(A.x * p, A.y * p);
}
Vector operator /(Vector A, double p)//
{
	return Vector(A.x / p, A.y / p);
}
bool operator <(const Point &a, const Point &b)//
{
	return a.x < b.x || (a.x == b.x && a.y < b.y);
}
const double eps = 1e-10;
int dcmp(double x)//
{
	if (fabs(x) < eps) return 0;
	else return x < 0 ? -1 : 1;
}
bool operator ==(const Point &a, const Point &b)//
{
	return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}
int t2(Point A, Point B)
{
	return (A.x - B.x) * (A.x - B.x) + (A.y - B.y) * (A.y - B.y);
}
double Dot(Vector A, Vector B)//
{
	return A.x * B.x + A.y * B.y;
}
double Length(Vector A)//
{
	return sqrt(Dot(A, A));
}
double Angle(Vector A, Vector B)//
{
	return acos(Dot(A, B) / Length(A) / Length(B));
}
double Cross(Vector A, Vector B)//
{
	return A.x * B.y - A.y * B.x;
}
bool read_point(Point &p)
{
	double x, y;
	scanf("%lf%lf", &x, &y);
	p = Point(x, y);
	return (x == 0 && y == 0);
}
Point A, B, C;
double lAB, lBC, lAC, r1, r2, r3, angA, angB, angC, mid, l, r, sqr1, sqr3, tlAC;
int main(int argc, char const *argv[])
{
	while (!(read_point(A) & read_point(B) & read_point(C)))
	{
		lAB = Length(A - B); lBC = Length(B - C); lAC = Length(A - C);
		angA = Angle(A - B, A - C), angB = Angle(B - A, B - C), angC = Angle(C - A, C - B);
		l = 0, r = min(lAB * tan(angB / 2), lBC * tan(angB / 2)), mid = r / 2;
		while (dcmp(r - l))
		{
			mid = (l + r) / 2;
			sqr1 = (sqrt(4 * mid  - 4 * (1 / tan(angA / 2) * (mid / tan(angB / 2) - lAB))) - 2 * sqrt(mid)) / (2 * 1 / tan(angA / 2));
			sqr3 = (sqrt(4 * mid  - 4 * (1 / tan(angC / 2) * (mid / tan(angB / 2) - lBC))) - 2 * sqrt(mid)) / (2 * 1 / tan(angC / 2));
			tlAC = sqr1 * sqr1 / tan(angA / 2) + sqr3 * sqr3 / tan(angC / 2) + 2 * sqr1 * sqr3;
			if (dcmp(tlAC - lAC) < 0)r = mid; else l = mid;
		}
		r2 = mid, r1 = sqr1 * sqr1, r3 = sqr3 * sqr3;
		printf("%.6lf %.6lf %.6lf\n", r1, r2, r3);
	}
	return 0;
}

大白书说是要推公式，恐怕这个公式相当不好推，湮没在无数的开平方中了，那么神器就是二分啦，

那么，假设一个半径，可以通过半径与底边的关系推出三元方程组，解方程就好。

值得注意的是，三个圆相切，暗示这有一个半径的约束范围，是什么呢？显然是要sqr1和sqr3有解啦。（不推出范围会导致sqr1和sqr3无解，得不出答案）。

终于做完了大白书几何的基础部分，终于圆满结束了，通过这一部分的联系，二分的水平有所提高。

等等！这还没完！这货是马尔法蒂问题，确实有公式可以推！

其中r为三角形内切圆的半径，s为半周长。

然后就可以用公式秒杀啦。

#include <bits/stdc++.h>
using namespace std;
struct Point
{
    double x, y;
};
double dis(Point a, Point b)
{
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
Point A, B, C, I;
double a, b, c, s, r, ia, ib, ic, r1, r2, r3;
int main()
{
    int i, j, k;
    while (~scanf("%lf%lf%lf%lf%lf%lf", &A.x, &A.y, &B.x, &B.y, &C.x, &C.y))
    {
        if (A.x == 0 && A.y == 0 && B.x == 0 && B.y == 0 && C.x == 0 && C.y == 0)
            break;
        a = dis(B, C);
        b = dis(A, C);
        c = dis(A, B);
        s = (a + b + c) / 2;
        I.x = (a * A.x + b * B.x + c * C.x) / (a + b + c);
        I.y = (a * A.y + b * B.y + c * C.y) / (a + b + c);
        r = sqrt(s * (s - a) * (s - b) * (s - c)) / s;
        ia = dis(I, A);
        ib = dis(I, B);
        ic = dis(I, C);
        r1 = r / (2 * (s - a)) * (s - r - (ib + ic - ia));
        r2 = r / (2 * (s - b)) * (s - r - (ia + ic - ib));
        r3 = r / (2 * (s - c)) * (s - r - (ia + ib - ic));
        printf("%.6f %.6f %.6f\n", r1, r2, r3);
    }
}

终于结束了。