【HDU】4773 Problem of Apollonius

题意

给定相离的两个圆(圆心坐标以及半径)以及圆外的一个定点\(P\),求出过点\(P\)的且与已知的两个圆外切的全部圆(输出总数+圆心、半径)。spa

分析

若是强行解方程,反正我是不会。
本题用到新姿式:圆的反演。
二维上的圆的反演一般是指定一个圆\(C\)为基础,其圆心\(O\)为反演中心,其半径\(r\)为反演半径。对于平面上任意一个非反演中心的点\(P\),都有且仅有一个反演点\(P'\),知足\(OP \cdot OP' = r^2\)\(P'\)\(OP\)射线上。对于任意一个二维上的图形,其反形就是图形上的全部点的反演点组成的。
因而能够证实:code

  1. \(C\)内的点的反演点在\(C\)外,\(C\)上的点的反演点就是自身,\(C\)外的点的反演点在\(C\)内。
  2. 任意一条不过反演中心的直线,其反形为过反演中心的圆。(至于过反演中心,我只能理解为这是极限意义下有点已经无限逼近了反演中心因此就算通过)
  3. 任意一个不过反演中心的圆,其反形为不过反演中心的圆,且反演中心为两圆的位似中心。

因而对于本题,交\(3\)个点的圆,能够看作一条直线关于\(P\)点反演获得的圆。因为反演点惟一,因此这条直线确定与给出的两个圆的反演圆各相交一个点。因此就是这两个圆的反演圆的公切线!咱们发现,若是是内公切线,反演成圆后会把一个圆包含,所以不符合题意。因此咱们只须要考虑外公切线便可。
可是外公切线的反形圆也可能会出现把给出的两个圆包含的状况,画一下图就能发现这种状况只出如今\(p\)和另外两个反演圆不在公切线的同一侧。ip

题解

本题要注意反演半径不要开过小,不然会有精度问题。(我是设成10)it

\((C_1, r_1)\)关于圆\((P, r)\)反演获得反形\((C_2, r_2)\),咱们来求\(r_2\)
根据反演式子能够获得:io

$$ \begin{align} (PC_1+r_1) \cdot (PC_2-r_2) = & r^2 \\ (PC_1-r_1) \cdot (PC_2+r_2) = & r^2 \\ \end{align} $$ class

消去\(PC_2\)能够解得:\[r_2 = \frac{r^2}{2} \left( \frac{1}{PC_1-r_1} - \frac{1}{PC_1+r_1} \right)\]那么再根据\(C_1\)的反演点是\(C_2\),咱们也能解出\(C_2\)来。基础

而后咱们须要求切线的反演圆\((O, r_3)\)了,假设切线与圆\((C_2, r_2)\)交于\(A\)
首先\(r_3\)能够由\(P\)到切线的距离\(t\)经过反演定义式求出:\[ 2r_3 \cdot t = r^2 \Rightarrow r_3 = \frac{r^2}{2t} \]而后因为直线\(OP\)\(AC_1\)是平行的(都与切线垂直)
因此根据类似来求出圆心坐标:\[ O = P + \frac{r_3}{r_2} \overrightarrow{C_2 A} \]di

#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
typedef double lf;
const lf eps=1e-8;
inline int dcmp(lf x) {
    return x>-eps?x>=eps:-1;
}
struct ip {
    lf x, y;
    ip(lf _x=0, lf _y=0) : x(_x), y(_y) {}
    void scan() {
        scanf("%lf%lf", &x, &y);
    }
};
typedef ip iv;
ip operator + (ip a, iv b) {
    return ip(a.x+b.x, a.y+b.y);
}
iv operator - (ip a, ip b) {
    return iv(a.x-b.x, a.y-b.y);
}
lf operator ^ (iv a, iv b) {
    return a.x*b.y-a.y*b.x;
}
iv operator * (lf k, iv a) {
    return iv(a.x*k, a.y*k);
}
inline lf sqr(lf x) {
    return x*x;
}
lf dis(ip a, ip b) {
    return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
int onleft(ip a, ip b, iv v) {
    return dcmp(v^(a-b))==1;
}
struct ic {
    ip o;
    lf r;
    void scan() {
        o.scan();
        scanf("%lf", &r);
    }
    ip gen(lf a) {
        return ip(o.x+r*cos(a), o.y+r*sin(a));
    }
    void print() {
        printf("%.8f %.8f %.8f\n", o.x, o.y, r);
    }
}o[2], ans[2], p;
ic inv(ic a) {
    ic b;
    lf d=dis(a.o, p.o), r2=sqr(p.r),
       k1=r2/(d-a.r), k2=r2/(d+a.r);
    b.r=0.5*(k1-k2);
    b.o=p.o+0.5*(k1+k2)/d*(a.o-p.o);
    return b;
}
ic inv(ip a, ip b) {
    ic c;
    c.r=sqr(p.r)/(abs(((b-a)^(p.o-a))/dis(a, b))*2);
    c.o=p.o+(c.r/o[0].r)*(a-o[0].o);
    return c;
}
int cal() {
    int tot=0;
    o[0]=inv(o[0]);
    o[1]=inv(o[1]);
    if(o[0].r<o[1].r) {
        swap(o[0], o[1]);
    }
    iv t=o[1].o-o[0].o;
    lf k1=atan2(t.y, t.x),
       k2=acos((o[0].r-o[1].r)/dis(o[0].o, o[1].o));
    ip a=o[0].gen(k1+k2), b=o[1].gen(k1+k2);
    t=b-a;
    if(onleft(o[0].o, a, t)==onleft(p.o, a, t)) {
        ans[tot++]=inv(a, b);
    }
    a=o[0].gen(k1-k2), b=o[1].gen(k1-k2);
    t=b-a;
    if(onleft(o[0].o, a, t)==onleft(p.o, a, t)) {
        ans[tot++]=inv(a, b);
    }
    return tot;
}
int main() {
    int T;
    scanf("%d", &T);
    while(T--) {
        o[0].scan();
        o[1].scan();
        p.o.scan();
        p.r=10;
        int tot=cal();
        printf("%d\n", tot);
        for(int i=0; i<tot; ++i) {
            ans[i].print();
        }
    }
    return 0;
}
相关文章
相关标签/搜索