jzoj3251. 【GDOI三校联考】Zone

Description

给出平面上的n个圆，求最小的把n个圆都包围起来的轮廓线长度。

Input

输入包括多组数据，第一行一个正整数T表示数据组数。
每组数据第一行一个正整数n表示圆的个数。
接下来n行，每行三个整数Xi，Yi，Ri表示一个圆心为（Xi，Yi）半径为Ri的圆。

Output

对于每组数据，输出一个实数表示最小轮廓线的长度。（答案与标准答案的绝对误差或相对误差小于10^-7即被认为是正确的）。

Sample Input

2
3
2 2 1
8 2 1
5 6 1
2
6 4 2
2 4 1

Sample Output

22.2831853072
17.6761051635

Data Constraint

对于10%的数据，n=2。
对于10%的数据，n=3。
对于另外40%的数据，所有圆的半径相等。
对于100%的数据，n<=250，|Xi|，|Yi|<=100，Ri<=100，T<=5。

题解

这题题意大概就是在原来凸包那个橡皮筋围钉子的模型，把钉子放大成一个圆。
然后就是很多个圆，要求用橡皮筋围住。
40%：这部分分是半径相等的圆。
我们可以看到jzoj2881. 【SHTSC2012day1】信用卡凸包这道题。
与之相似。
100%：
怎么做？我们发现，由于圆的半径相等时很好处理，把圆心求一个凸包即可。
然而，半径并不相等。
怎么办？
我们一步一步来搞。
先考虑求任意两个圆之间的公切线的交点
我们发现，对于两个相离的圆，它们有两个外公切线和内公切线。

然后，我们发现，内公切线（红色）是无法包围两个圆的，于是乎，在这题中我们只需要求外公切线与圆的交点即可。
怎么求？
我们看到一道题——
http://acm.zju.edu.cn/onlinejudge/showProblem.do?problemId=199
这题求的是两条外公切线的交点

这个点怎么求？
利用相似！！！
设一个圆点为(x1,y1)半径为r1，另一个(x2,y2)半径为r2
那么要求的点坐标即为：
$x=(x2*r1-x1*r2)/(r1-r2)$
$y=(y2*r1-y1*r2)/(r1-r2)$
然而这样子求出来的点不是公切线，我们就考虑推出与圆的交点。
由于两条线是切线，所以我们可以得到这个图：

我们可以得到a,b,c的长度，以及A点与B点坐标。
现在就是要求C的坐标了。怎么求？
一个想法——
我们先求出θ的大小——
$θ=cos^{-1}(frac{b^2+c^2-a^2}{2bc})$
然后利用向量来求出C点坐标。
这里有向量的做法——https://blog.csdn.net/hjq376247328/article/details/45113563/
然而我对这个向量旋转的方法感到迷惑，怎么办？
考虑点旋转！
先把B旋转到AC线上去。
$nx1(新B点横坐标)=(x1-x0)*cos(θ)-(y1-y0)*sin(θ)+x0$
$ny1(新B点纵坐标)=(x1-x0)*sin(θ)+(y1-y0)*cos(θ)+y0$
然后，我们考虑相似把这个坐标再放大或缩小得到的结果——
$x2=frac{(nx1-x0)*c}{b}+x0$
$y2=frac{(ny1-y0)*c}{b}+y0$
应该是这样了。
这样我们还要考虑C在AB下方的情况，我们只需要把上面的θ变成-θ即可。
至此，我们就把公切线与圆的交点求出来了。

如果两个圆半径相等，怎么办？
直接利用相似即可。

求出来后怎么办？
凸包！
最后我们做完凸包后，可以得到一个类似于下面的一个图——

对于蓝色边，我们可以直接得到长度，但是如果两个红点在同圆上呢？
我们还好求弧长。
由于只需要求出相对应的圆心角即可，于是我们再利用上面求θ的方法求即可。
看起来是一道码农题啊。

半成品

（这个代码是在判断圆上弧为优弧难以判断的代码）

uses math;
type
        new=record
                x,y,jj,jl:extended;
                kind,bh:longint;
        end;
        new1=record
                x,y,r:extended;
        end;
        new2=record
                a,b,c:extended;
        end;
var
        i,j,k,l,n,m,gs,t,num,now:longint;
        cir:array[1..250] of new1;
        jd,a,b:new;
        spot:array[1..250000] of new;
        dl:array[1..250000] of longint;
        dt,nx,ny,answer:extended;
        tir:new2;
function dist(x1,y1,x2,y2:extended):extended;
begin
        exit(sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)));
end;
function pyt(a,b:extended):extended;
begin
        exit(sqrt(a*a-b*b));
end;
procedure insert(i:longint);
begin
        nx:=((cir[i].x-jd.x)*cos(dt)-(cir[i].y-jd.y)*sin(dt))*tir.c/tir.b+jd.x;
        ny:=((cir[i].x-jd.x)*sin(dt)+(cir[i].y-jd.y)*cos(dt))*tir.c/tir.b+jd.y;
        inc(num);
        spot[num].x:=nx;spot[num].y:=ny;spot[num].kind:=now;spot[num].bh:=i;
end;
function judge(ax,ay,bx,by,cx,cy:extended):extended;
begin
        exit(ax*by+cx*ay+bx*cy-cx*by-bx*ay-ax*cy);
end;
procedure con_bag;
var
        i,j,k,up:longint;
        x,y,miny,a,b,c:extended;
procedure qsort(l,r:longint);
var
        i,j:longint;
        m,m1:extended;
        k:new;
begin
        i:=l;j:=r;
        m:=spot[(l+r) div 2].jj;
        m1:=spot[(l+r) div 2].jl;
        repeat
                while (spot[i].jj<m) or ((spot[i].jj=m) and (spot[i].jl<m1)) do inc(i);
                while (spot[j].jj>m) or ((spot[j].jj=m) and (spot[j].jl>m1)) do dec(j);
                if i<=j then
                begin
                        k:=spot[i];
                        spot[i]:=spot[j];
                        spot[j]:=k;
                        inc(i);dec(j);
                end;
        until i>j;
        if l<j then qsort(l,j);
        if r>i then qsort(i,r);
end;
begin
        gs:=0;
        miny:=maxlongint;
        for i:=1 to num do
        begin
                if spot[i].y<miny then
                begin
                        miny:=spot[i].y;
                        j:=i;
                end;
        end;
        x:=spot[j].x;
        y:=spot[j].y;
        for i:=1 to num do
        begin
                spot[i].x:=spot[i].x-x;
                spot[i].y:=spot[i].y-y;
                if spot[i].y=0 then spot[i].jj:=0
                else
                if spot[i].x=0 then spot[i].jj:=90
                else
                if spot[i].x<0 then spot[i].jj:=arctan(spot[i].y/spot[i].x)*(180/pi)+180
                else spot[i].jj:=arctan(spot[i].y/spot[i].x)*(180/pi);
                spot[i].jl:=dist(0,0,spot[i].x,spot[i].y);
        end;
        qsort(1,num);
        up:=2;
        dl[1]:=1;
        dl[2]:=2;
        for i:=3 to num do
        begin
                while judge(spot[dl[up-1]].x,spot[dl[up-1]].y,spot[dl[up]].x,spot[dl[up]].y,spot[i].x,spot[i].y)<0 do dec(up);
                inc(up);
                dl[up]:=i;
        end;
        inc(up);
        dl[up]:=dl[1];
        answer:=0;
        for i:=2 to up do
        begin
                if spot[dl[i]].kind=spot[dl[i-1]].kind then
                begin
                        answer:=answer+dist(spot[dl[i]].x,spot[dl[i]].y,spot[dl[i-1]].x,spot[dl[i-1]].y);
                end
                else
                begin
                        a:=dist(spot[dl[i]].x,spot[dl[i]].y,spot[dl[i-1]].x,spot[dl[i-1]].y);
                        b:=cir[spot[dl[i]].bh].r;
                        c:=cir[spot[dl[i]].bh].r;
                        dt:=arccos((b*b-a*a+c*c)/(2*b*c))*(180/pi);
                        answer:=answer+2*c*dt*pi/360;
                end;
        end;
        writeln(answer:0:10);
end;
procedure qsortc(l,r:longint);
var
        i,j:longint;
        m,m1:extended;
        k:new1;
begin
        i:=l;j:=r;
        m:=cir[(l+r) div 2].x;
        m1:=cir[(l+r) div 2].y;
        repeat
                while (cir[i].x<m) or ((cir[i].x=m) and (cir[i].y<m1)) do inc(i);
                while (cir[j].x>m) or ((cir[j].x=m) and (cir[j].y>m1)) do dec(j);
                if i<=j then
                begin
                        k:=cir[i];
                        cir[i]:=cir[j];
                        cir[j]:=k;
                        inc(i);dec(j);
                end;
        until i>j;
        if l<j then qsortc(l,j);
        if r>i then qsortc(i,r);
end;
begin
        //assign(input,'zone.in');reset(input);
        readln(t);
        while t>0 do
        begin
                dec(t);
                readln(n);
                for i:=1 to n do
                begin
                         readln(cir[i].x,cir[i].y,cir[i].r);
                end;
                qsortc(1,n);
                num:=0;now:=0;
                for i:=1 to n do
                begin
                        for j:=i+1 to n do
                        begin
                                tir.c:=dist(cir[i].x,cir[i].y,cir[j].x,cir[j].y);
                                if (cir[i].r<>cir[j].r) and (tir.c>abs(cir[i].r-cir[j].r)) then
                                begin
                                        inc(now);
                                        jd.x:=(cir[j].x*cir[i].r-cir[i].x*cir[j].r)/(cir[i].r-cir[j].r);
                                        jd.y:=(cir[j].y*cir[i].r-cir[i].y*cir[j].r)/(cir[i].r-cir[j].r);
                                        tir.a:=cir[i].r;
                                        tir.b:=pyt(tir.c,tir.a);
                                        dt:=arccos((tir.b*tir.b-tir.a*tir.a+tir.c*tir.c)/(2*tir.b*tir.c));
                                        insert(i);insert(j);
                                        dt:=-dt;
                                        insert(i);insert(j);
                                end
                                else
                                if (cir[i].r=cir[j].r) and (tir.c>abs(cir[i].r-cir[j].r)) then
                                begin
                                        inc(now);
                                        if (cir[i].x<>cir[j].x) and (cir[i].y<>cir[j].y) then
                                        begin
                                                dt:=arctan(abs(cir[i].y-cir[j].y)/abs(cir[i].x-cir[j].x));
                                                jd.y:=cos(dt)*cir[i].r;
                                                jd.x:=sin(dt)*cir[i].r;
                                                if cir[i].y<cir[j].y then
                                                begin
                                                        inc(num);spot[num].x:=cir[i].x-jd.x;spot[num].y:=cir[i].y+jd.y;spot[num].kind:=now;spot[num].bh:=i;
                                                        inc(num);spot[num].x:=cir[j].x-jd.x;spot[num].y:=cir[j].y+jd.y;spot[num].kind:=now;spot[num].bh:=j;
                                                        inc(num);spot[num].x:=cir[i].x+jd.x;spot[num].y:=cir[i].y-jd.y;spot[num].kind:=now;spot[num].bh:=i;
                                                        inc(num);spot[num].x:=cir[j].x+jd.x;spot[num].y:=cir[j].y-jd.y;spot[num].kind:=now;spot[num].bh:=j;
                                                end
                                                else
                                                begin
                                                        inc(num);spot[num].x:=cir[i].x-jd.x;spot[num].y:=cir[i].y-jd.y;spot[num].kind:=now;spot[num].bh:=i;
                                                        inc(num);spot[num].x:=cir[j].x-jd.x;spot[num].y:=cir[j].y-jd.y;spot[num].kind:=now;spot[num].bh:=j;
                                                        inc(num);spot[num].x:=cir[i].x+jd.x;spot[num].y:=cir[i].y+jd.y;spot[num].kind:=now;spot[num].bh:=i;
                                                        inc(num);spot[num].x:=cir[j].x+jd.x;spot[num].y:=cir[j].y+jd.y;spot[num].kind:=now;spot[num].bh:=j;
                                                end;
                                        end
                                        else
                                        if (cir[i].x=cir[j].x) then
                                        begin
                                                inc(num);spot[num].x:=cir[i].x-cir[i].r;spot[num].y:=cir[i].y;spot[num].kind:=now;spot[num].bh:=i;
                                                inc(num);spot[num].x:=cir[j].x-cir[i].r;spot[num].y:=cir[j].y;spot[num].kind:=now;spot[num].bh:=j;
                                                inc(num);spot[num].x:=cir[i].x+cir[i].r;spot[num].y:=cir[i].y;spot[num].kind:=now;spot[num].bh:=i;
                                                inc(num);spot[num].x:=cir[j].x+cir[i].r;spot[num].y:=cir[j].y;spot[num].kind:=now;spot[num].bh:=j;
                                        end
                                        else
                                        if (cir[i].y=cir[j].y) then
                                        begin
                                                inc(num);spot[num].x:=cir[i].x;spot[num].y:=cir[i].y-cir[i].r;spot[num].kind:=now;spot[num].bh:=i;
                                                inc(num);spot[num].x:=cir[j].x;spot[num].y:=cir[j].y-cir[i].r;spot[num].kind:=now;spot[num].bh:=j;
                                                inc(num);spot[num].x:=cir[i].x;spot[num].y:=cir[i].y+cir[i].r;spot[num].kind:=now;spot[num].bh:=i;
                                                inc(num);spot[num].x:=cir[j].x;spot[num].y:=cir[j].y+cir[i].r;spot[num].kind:=now;spot[num].bh:=j;
                                        end;
                                end;

                        end;
                end;
                con_bag;
        end;
end.