hdu 4773 圆的反演变换

Problem of Apollonius

                                                                          Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
                                                                           Total Submission(s): 744    Accepted Submission(s): 182
Special Judge

Problem Description

　　Apollonius of Perga (ca. 262 BC - ca. 190 BC) was a Greek geometer and astronomer. In his noted work Epaphai, he posed and solved such a problem: constructing circles that are tangent to three given circles in a plane. Two tangent circles can be internally or externally tangent to each other, thus Apollonius's problem generically have eight solutions. 
　　Now considering a simplified case of Apollonius's problem: constructing circles that are externally tangent to two given circles, and touches a given point(the given point must be on the circle which you find, can't be inside the circle). In addition, two given circles have no common points, and neither of them are contained by the other, and the given point is also located strictly outside the given circles. You should be thankful that modern mathematics provides you with plenty of useful tools other than euclidean geometry that help you a lot in this problem.

 

Input

　　The first line of input contains an integer T (T ≤ 200), indicating the number of cases.
　　Each ease has eight positive integers x1, y1, r1, x2, y2, r2, x3, y3 in a single line, stating two circles whose centres are (x1, y1), (x2, y2) and radius are r1 and r2 respectively, and a point located at (x3, y3). All integers are no larger than one hundred.

 

Output

　　For each case, firstly output an integer S, indicating the number of solutions.
　　Then output S lines, each line contains three float numbers x, y and r, meaning that a circle, whose center is (x, y) and radius is r, is a solution to this case. If there are multiple solutions (S > 1), outputing them in any order is OK. Your answer will be accepted if your absolute error for each number is no more than 10 ^-4. 

 

Sample Input

  
  

   
   1
12 10 1 8 10 1 10 10
  
  

 

Sample Output

  
  

   
   2
10.00000000 8.50000000 1.50000000
10.00000000 11.50000000 1.50000000


   
   

    
    

     
     Hint
    
    
This problem is special judged.

   
   
   
    
  
  

 

可以参考一下这篇博客：http://blog.csdn.net/acdreamers/article/details/16966369#comments

给定不相交的两个圆以及圆外一点，找一个经过给定点的圆与其他两个圆相切

首先来看反演变换，首先是给定一个圆圆心为O，半径为R

1、圆外一点P与圆内一点P‘会一一对应的反演OP*OP'=R*R

2、经过O的圆，反演后成为不经过O的一条直线

3、不经过O的圆，反演后成为另一个圆，且圆心并不对应

4、不经过O的直线反演后成为一个经过O的圆

5、过O的直线反演后不变

那么这道题就把两个圆反演之后求两圆的公切线，然后反演回去，就会成为一个过O的圆，且与另两圆相切（由反演的过程可看出）

求反演后的圆的圆心，可以将过O的直线所在的那条直径两端点反演，然后再求圆心.

#include <iostream>
#include <string.h>
#include <stdio.h>
#include <math.h>
using namespace std;
double const eps = 1e-8;
struct Point
{
	double x, y;
	Point(double a = 1.0, double b = 1.0) :x(a), y(b) {}
	Point operator + (const Point &a)
	{
		return Point(x + a.x, y + a.y);
	}
	Point operator - (const Point &a)
	{
		return Point(x - a.x, y - a.y);
	}
	Point operator * (const double a)
	{
		return Point(a*x, a*y);
	}
	Point Trans()
	{
		return Point(-y, x);
	}
	void Input()
	{
		scanf("%lf%lf", &x, &y);
	}
};
struct Circle
{
	Point o;
	double r;
	Circle(Point a = Point(), double b = 1.0) :o(a), r(b) {}
	Point getPoint(double alpha)
	{
		return o + Point(r*cos(alpha), r*sin(alpha));
	}
	void Input()
	{
		o.Input();
		scanf("%lf", &r);
	}
	void Output()
	{
		printf("%.8lf %.8lf %.8lf\n", o.x, o.y, r);
	}
};
Point p;
Circle c[15];
double dist(Point A, Point B)
{
	return sqrt((A.x - B.x)*(A.x - B.x) + (A.y - B.y)*(A.y - B.y));
}
double cross(Point A, Point B, Point C)
{
	return (B.x - A.x)*(C.y - A.y) - (B.y - A.y)*(C.x - A.x);
}
int sign(double x)
{
	return (x > eps) - (x < -eps);
}
Circle Inverse(Circle C)
{
	Circle T;
	double t = dist(C.o, p);
	double x = 1.0 / (t - C.r);
	double y = 1.0 / (t + C.r);
	T.r = (x - y) / 2.0;
	double s = (x + y) / 2.0;
	T.o = p + (C.o - p) * (s / t);
	return T;
}
void add(Point a, Point b, int &k)
{
	double t = cross(a, p, b);
	if (t < 0) t = -t;
	double d = dist(a, b);
	t /= d;
	if (t > eps)
	{
		double w = 0.5 / t;
		Point dir = (b - a).Trans();
		Point a1 = p + dir * (w / d);
		Point b1 = p - dir * (w / d);
		if (fabs(cross(a, b, a1)) < fabs(cross(a, b, b1)))
			c[k++] = Circle(a1, w);
		else
			c[k++] = Circle(b1, w);
	}
}
int Work()
{
	c[0] = Inverse(c[0]);
	c[1] = Inverse(c[1]);
	if (c[1].r > c[0].r) swap(c[1], c[0]);
	Point v = c[1].o - c[0].o;
	double alpha = atan2(v.y, v.x);
	double d = dist(c[0].o, c[1].o);
	double beta = acos((c[0].r - c[1].r) / d);
	int k = 2;
	Point a = c[0].getPoint(alpha + beta);
	Point b = c[1].getPoint(alpha + beta);
	if (sign(cross(a, c[0].o, b)) == sign(cross(a, p, b)) &&	sign(cross(a, c[1].o, b)) == sign(cross(a, p, b)))
		add(a, b, k);
	a = c[0].getPoint(alpha - beta);
	b = c[1].getPoint(alpha - beta);
	if (sign(cross(a, c[0].o, b)) == sign(cross(a, p, b)) &&sign(cross(a, c[1].o, b)) == sign(cross(a, p, b)))
		add(a, b, k);
	return k - 2;
}
int main()
{
	int T;
	scanf("%d", &T);
	while (T--)
	{
		c[0].Input();
		c[1].Input();
		p.Input();
		int num = Work();
		printf("%d\n", num);
		for (int i = 0;i<num;i++)
			c[i + 2].Output();
	}
	return 0;
}