bzoj2823

最小圆覆盖

有个东西叫作随机增量法，具体可以baidu

这里来说说怎么求三点共圆

这其实就是求两条线段的交点

在编程中，我们解方程是比较麻烦的一个比较好的方法是利用相似三角形

设线段AB，CD交P，则PC：PD=Sabc:Sabd

然后用定比分点就可以求的交点坐标了

 1 const eps=1e-6;
 2
 3 type point=record
 4        x,y:double;
 5      end;
 6
 7 var p:array[0..500010] of point;
 8     i,n,j,k:longint;
 9     r:double;
10
11 procedure swap(var a,b:point);
12   var c:point;
13   begin
14     c:=a;
15     a:=b;
16     b:=c;
17   end;
18
19 function dis(i,j:longint):double;
20   begin
21     exit(sqrt(sqr(p[i].x-p[j].x)+sqr(p[i].y-p[j].y)));
22   end;
23
24 function cov(i:longint):boolean;
25   begin
26     if dis(i,0)-r>eps then exit(false)
27     else exit(true);
28   end;
29
30 procedure little(i,j:longint);
31   begin
32     p[0].x:=(p[i].x+p[j].x)/2;
33     p[0].y:=(p[i].y+p[j].y)/2;
34     r:=dis(i,j)/2;
35   end;
36
37 function cross(x1,y1,x2,y2:double):double;
38   begin
39     exit(x1*y2-x2*y1);
40   end;
41
42 procedure circle(i,j,k:longint);
43   var ret,xa,ya,xb,yb,xc,yc,xd,yd,t1,t2:double;
44   begin
45     ret:=100/dis(i,j);
46     xa:=(p[i].x+p[j].x)/2; ya:=(p[i].y+p[j].y)/2;
47     xb:=xa-ret*(p[i].y-p[j].y); yb:=ya+ret*(p[i].x-p[j].x);
48     ret:=100/dis(j,k);
49     xc:=(p[j].x+p[k].x)/2; yc:=(p[j].y+p[k].y)/2;
50     xd:=xc-ret*(p[j].y-p[k].y); yd:=yc+ret*(p[j].x-p[k].x);
51     t1:=cross(xc-xa,yc-ya,xb-xa,yb-ya);
52     t2:=cross(xb-xa,yb-ya,xd-xa,yd-ya);
53     p[0].x:=(t1*xd+t2*xc)/(t1+t2);
54     p[0].y:=(t1*yd+t2*yc)/(t1+t2);
55     r:=dis(0,i);
56  end;
57
58 begin
59   readln(n);
60   for i:=1 to n do
61   begin
62     readln(p[i].x,p[i].y);
63     swap(p[i],p[trunc(random*i)+1]);
64   end;
65   little(1,2);
66   for i:=3 to n do
67     if not cov(i) then
68     begin
69       little(1,i);
70       for j:=2 to i-1 do
71         if not cov(j) then
72         begin
73           little(i,j);
74           for k:=1 to j-1 do
75             if not cov(k) then circle(i,j,k);
76         end;
77     end;
78   writeln(p[0].x:0:2,‘ ‘,p[0].y:0:2,‘ ‘,r:0:2);
79 end.