hdu 4773 Problem of Apollonius—2013年杭州赛区D题—反演变换

最新推荐文章于 2020-05-04 06:38:49 发布

dechi5213

最新推荐文章于 2020-05-04 06:38:49 发布

阅读量139

点赞数

原文链接：http://www.cnblogs.com/wangsouc/articles/3649201.html

版权

题意：给定两个圆的圆心（x1,y1,x2,y2）和半径(r1,r2)以及一个点(x0,y0)（两圆不相交且不包含，点不在圆内），求满足下面条件的圆：1.给定点在圆周上。2.和给出的两个圆外切。输出满足条件的圆的个数以及圆心和半径。

解题思路：先补充一下这道题需要用到的反演变换的知识：

1.定义：设在平面内给定一点O和常数k(k不等于零)，对于平面内任意一点A，确定A′，使A′在直线OA上一点，并且有向线段OA与OA′满足OA·OA′=k，我们称这种变换是以O为的反演中心，以k为反演幂的反演变换，简称反演。称A′为A关于O（r）的互为反演点.

2.不过反演中心的圆经过反演变换仍然是一个不过反演中心的圆.

3.不过反演中心的直线经过反演变换是一个经过反演中心的圆.

4.反演变换不改变图形的相切性.

本题就是先把给定的两个圆（O1, O2）经过一次反演变换得到新的圆(O3, O4)，然后求O3,O4的外公切线(L)，显然此时若把L和O3,O4都经过一次反演变换，L变成了过反演中心的圆(O5)，O3,O4分别变成了原来的圆O1,O2，由相切性不变可知，O5和O1,O2相切。故我们只需要求出L反演后的圆即可。

  1 #include <stdio.h>
  2 #include <string.h>
  3 #include <math.h>
  4 #define eps 0.00000001
  5 double x_0, y_0;
  6 double a[5], b[5], c[5];
  7 double ansx[5], ansy[5], ansr[5];
  8 int flag[5];
  9 // 以x_0,y_0做为反演中心，1作为反演幂
 10 void exchange(double &x, double &y)
 11 {
 12      double tmp = x;
 13      x = y;
 14      y = tmp;
 15      return ;
 16 }
 17 
 18 double dis(double x_1, double y_1, double x_2, double y_2)
 19 {
 20      return sqrt((x_1 - x_2) * (x_1 - x_2) + (y_1 - y_2) * (y_1 - y_2));
 21 }
 22 
 23 void changepoint(double x, double y, double &ex, double &ey) //点反演变换
 24 {
 25      double d = dis(x, y, x_0, y_0);
 26      ey = y_0 + (y - y_0) / (d * d);
 27      ex = x_0 + (x - x_0) / (d * d);
 28      return ;
 29 }
 30 
 31 void changecircle(double x, double y, double r, double &ex, double &ey, double &er) // 圆反演变换
 32 {
 33      double x_1, y_1, x_2, y_2;
 34      if(fabs(x - x_0) < eps)
 35      {
 36           changepoint(x, y + r, x_1, y_1);
 37           changepoint(x, y- r, x_2, y_2);
 38      }
 39      else
 40      {
 41           double k = (y - y_0) / (x - x_0);
 42           double tmp = sqrt(r * r / (1 + k * k));
 43           changepoint(tmp + x, k * tmp + y, x_1, y_1);
 44           changepoint(-tmp + x, -k * tmp + y, x_2, y_2);
 45      }
 46      ex = (x_1 + x_2) / 2;
 47      ey = (y_1 + y_2) / 2;
 48      er = dis(x_1, y_1, x_2, y_2) / 2;
 49      return ;
 50 }
 51 
 52 void rotate(double s, double t, double &aa, double &bb, double &cc, double sinct, double cosct) // 直线aax + bby = cc绕s,t点旋转ct角度，其中sin(ct) = sinct,cos(ct) = cosct
 53 {
 54      double tmpa = aa, tmpb = bb, tmpc = cc;
 55      aa = tmpa * cosct - tmpb * sinct;
 56      bb = tmpa * sinct + tmpb * cosct;
 57      cc = tmpc + tmpa * s * cosct + tmpa * t * sinct - tmpa * s
 58           - tmpb * s * sinct + tmpb * t * cosct - tmpb * t;
 59      return ;
 60 }
 61 
 62 void getline(double x_1, double y_1, double r_1, double x_2, double y_2, double r_2) // 求两个圆的外公切线
 63 {
 64      double tmpa, tmpb, tmpc;
 65      tmpa = y_2 - y_1;
 66      tmpb = x_1 - x_2;
 67      tmpc = y_1 * (x_1 - x_2) - x_1 * (y_1 - y_2);
 68      double d, detr;
 69      d = dis(x_1, y_1, x_2, y_2);
 70      detr = r_2 - r_1;
 71      double sinct = fabs(detr) / d;
 72      double cosct = sqrt(d * d - detr * detr) / d;
 73      double aa = tmpa, bb = tmpb, cc = tmpc;
 74      rotate(x_1, y_1, aa, bb, cc, sinct, cosct);
 75      a[0] = aa;
 76      b[0] = bb;
 77      c[0] = cc - r_1 * sqrt(a[0] * a[0] + b[0] * b[0]);
 78      aa = tmpa, bb = tmpb, cc = tmpc;
 79      rotate(x_1, y_1, aa, bb, cc, -sinct, cosct);
 80      a[1] = aa;
 81      b[1] = bb;
 82      c[1] = cc + r_1 * sqrt(a[1] * a[1] + b[1] * b[1]);
 83      return ;
 84 }
 85 
 86 void getcircle(double aa, double bb, double cc, int k) // 求直线aax + bby = cc经过反演变换后的圆
 87 {
 88      if(fabs(aa * x_0 + bb * y_0 - cc) < eps)
 89      {
 90           flag[k] = 0;
 91           return ;
 92      }
 93      double s, t, ex, ey;
 94      s = (aa * cc + x_0 * bb * bb - aa * bb * y_0) / (aa * aa + bb * bb);
 95      t = (bb * cc - x_0 * aa * bb + aa * aa * y_0) / (aa * aa + bb * bb);
 96      changepoint(s, t, ex, ey);
 97      ansx[k] = (x_0 + ex) / 2;
 98      ansy[k] = (y_0 + ey) / 2;
 99      ansr[k] = dis(x_0, y_0, ex, ey) / 2;
100      return ;
101 }
102 
103 double x_00,y_00,x_11,y_11,x_22,y_22,r_11,r_22;
104 bool check(double x, double y, double r)
105 {
106      if (fabs(dis(x,y,x_00,y_00) - r) < eps)
107      if (fabs(dis(x,y,x_11,y_11) - r - r_11) < eps)
108      if (fabs(dis(x,y,x_22,y_22) - r - r_22) < eps)
109           return true;
110      return false;
111 }
112 
113 double get_dis(double aa, double bb, double cc,double x, double y)
114 {
115      return fabs((aa*x + bb*y - cc)/sqrt(aa*aa + bb*bb));
116 }
117 
118 int main()
119 {
120 
121      int T = 100;
122      double x_1, y_1, r_1, x_2, y_2, r_2;
123      scanf("%d", &T);
124      while(T--)
125      {
126           scanf("%lf%lf%lf", &x_1, &y_1, &r_1);
127           scanf("%lf%lf%lf", &x_2, &y_2, &r_2);
128           scanf("%lf%lf", &x_0, &y_0);
129           x_00 = x_0;
130           y_00 = y_0;
131           x_11 = x_1;
132           y_11 = y_1;
133           x_22 = x_2;
134           y_22 = y_2;
135           r_11 = r_1;
136           r_22 = r_2;
137 
138           double ex_1, ex_2, ey_1, ey_2, er_1, er_2;
139 
140           changecircle(x_1, y_1, r_1, ex_1, ey_1, er_1);
141           changecircle(x_2, y_2, r_2, ex_2, ey_2, er_2);
142           if(er_1 > er_2)
143           {
144                exchange(ex_1, ex_2);
145                exchange(ey_1, ey_2);
146                exchange(er_1, er_2);
147           }
148           getline(ex_1, ey_1, er_1, ex_2, ey_2, er_2);
149           memset(flag, -1, sizeof(flag));
150           getcircle(a[0], b[0], c[0], 0);
151           getcircle(a[1], b[1], c[1], 1);
152           int ans = 0;
153           for (int i=0; i < 2;i++)
154           {
155                if (check(ansx[i],ansy[i],ansr[i])){
156                     ans++;
157                     flag[i] = 1;
158                }
159           }
160 
161           printf("%d\n", ans);
162 
163           for(int i = 0; i < 2; i++)
164           {
165                if(flag[i] == 1){
166                     printf("%.8lf %.8lf %.8lf\n", ansx[i], ansy[i], ansr[i]);
167                }
168           }
169      }
170      return 0;
171 }
172 
173 /*
174 10
175 12 10 1 8 10 1 10 10
176 20 40 5 -20 -40 6 5 5
177 2 4 8 33 -55 9 50 50
178 1 2 5 100 200 8 30 50
179 2 5 10 -1000 -2000 19 -10 -20
180 2 10 3 2 30 2 4 4
181 10 100 5 -100 100 8 5 5
182 */