_皎月半洒花 的博客

_皎月半洒花 的博客

反抗命运的镇魂曲停歇之前,你绝对不可以说我懦弱,绝对不可以说我没有殊死一搏。

题解 CF2C 【Commentator problem】

posted on 2019-02-12 16:31:43 | under 题解 |

一道被严重低估难度的计算几何。

喂喂凭什么CF1C是个蓝题,这个题这么麻烦是个黄题啊……

首先我们要知道本题问的是什么:

给出三个互不相交的圆,求一个点使得到这三个圆的切线夹角相同。

事实上,我们设这点为 $T$ , 三个圆心分别为 $A, B,C$ 。而事实上,圆 $A$ 的半径 $r_A$ 与 $dis(A,T)$ 的比值,就是 $sin(\frac{1}{2}\angle A_1TA_2)$ ,其中 $A_1$ 和 $A_2$ 是过T的圆A的两条切线与圆的交点。

那么也就是说,我们如果有 $\angle A_1TA_2 = \angle B_1TB_2= \angle C_1TC_2$ ,那么一定有 $$\frac{r_A}{dis(A,T)} = \frac{r_B}{dis(B,T)} = \frac{r_C}{dis(C,T)}$$

稍微移一下项,就会有 $$\frac{r_A}{r_B} = \frac{dis(A,T)}{dis(B,T)} $$

那么我们就可以发现,对于两个点而言,我们要找的目标点 $T$ 满足到两个点的距离等于一个给定的比例( $r_A$ 和 $r_B$ 给定)。

事实上,这样的点的轨迹是可以刻画的。我们列一个方程即可:

设比例系数为 $k(k \geq 1)$ , 那么:

$$\frac{\sqrt{(x_T - x_A)^2 + (y_T - y_A)^2 }}{\sqrt{(x_T - x_B)^2 + (y_T - y_B)^2 }} = k$$

$$\frac{(x_T - x_A)^2 + (y_T - y_A)^2 }{(x_T - x_B)^2 + (y_T - y_B)^2 } = k^2$$

稍微移一下项就会得到 $$(k^2-1)x_T^2 + (k^2-1)y_T^2 - 2(k^2y_B - y_A)y_T - 2(k^2x_B - x_A)x_T+k^2x_B^2 - x_A^2 + k^2y_B^2 - y_A^2 = 0$$

看起来有点儿长……我们令 $A = k^2-1, C = - 2(k^2x_B - x_A), D = -2(k^2y_B - y_A), E = k^2x_B^2 - x_A^2 + k^2y_B^2 - y_A^2$

那么就会变成 $Ax^2 + Ay^2+ Cx + Dy+E = 0$ ,由于 $A,C,D,E$ 都是常数,所以这是一个圆的一般方程。

我们其实也可以发现,当 $k=1$ 时,此时为一条直线(即中垂线),换句话说当且仅当两个圆半径相等时,点 $T$ 的轨迹是一条直线。其余的情况则是一个圆

我们不妨先记这种到两个圆的圆心的距离成定比例的轨迹为两个圆的生成曲线

那么之后呢,我们发现,圆 $A$ 和圆 $C$ 的生成曲线,与圆 $A$ 和圆 $B$ 的生成曲线,至多有两个交点。那么我们只需要:

  • $(1)~~~$ 判断三组圆的生成曲线是否都相交且交于一点,不是则无解。
  • $(2)~~~$ 对于其中两个圆的生成曲线的交点,判断是否满足条件,即是我们已经找到了符合 $$\frac{r_A}{r_B} = \frac{dis(A,T)}{dis(B,T)}$$ 的点,我们需要判断对于圆 $C$ 是否也满足 $$\frac{r_A}{r_C} = \frac{dis(A,T)}{dis(C,T)}$$
  • $(3)~~~$ 如果选取的生成曲线恰好有 $2$ 个交点且两个交点 $T',T''$ 都满足 $(2)$ 中的条件,那么我们选 $sin$ 值最大的(对于 $\leq \frac{\pi}{2}$ 的角, $sin$ 值与角的大小成正相关)。

然后算法就结束了。中间还有好多好多好多问题,比如圆与圆的交点怎么求,直线与直线的交点怎么求,圆与直线的交点怎么求……不说了,看必修二去吧phhh

ps:我求直线与直线的交点是先转化成斜截式再求的,不是因为别的,只是因为比较熟悉并且更重要的是,我不久前写完一道用斜截式做的计算几何……于是就挪过来用了。

代码很繁琐,轻喷(含注释):

$\color{red}{C}\color{cyan}{o}\color{gold}{d}\color{green}{e}$


#include <cmath>
#include <cstdio>
#include <iostream>
#include <algorithm>

using namespace std ;
const double Eps = 1e-3 ; int i ;//以下的mark都是记录状态
struct Node{ int mark ; double xa, ya, xb, yb ; } I[5] ; // 0 = inexist, 1 = exist*1, 2 = exist*2 ;
//此处我的Node存的实际上是两个点,即一个一元二次方程的两个解。
struct Line{ int mark ; double k, b ; double x, y ; }L[12] ; //0 : // x-axis, 1: // y-axis, 2: // normal ;
struct Circle{
    int mark ; // 1 : circle ; 0 : Line ;
    double x, y, r ;
    double A, B, C, D, E ;
    Circle friend operator -(const Circle &A, const Circle &B){
        return (Circle){0, A.x - B.x, A.y - B.y, A.r - B.r, A.A - B.A, A.B - B.B, A.C - B.C, A.D - B.D, A.E - B.E} ;
    }
}C[10] ;
double ansx, ansy ; bool check ;

inline bool Comp(Node A, Node B){ return A.mark < B.mark ; }
inline double get_x(Line A, Line B){ return A.mark == 0 ? A.x : B.x ; } //which is (x = k) ;
inline double get_y(Line A, Line B){ return A.mark == 1 ? A.y : B.y ; } //which is (y = k) ;
inline bool equal(double x, double y){ return (x - y <= Eps) && (x - y >= -Eps) ; }
inline double disa(Node A, Node B){ return sqrt((A.xa - B.xa) * (A.xa - B.xa) + (A.ya - B.ya) * (A.ya - B.ya)); }//第一个点之间的距离
inline double disb(Node A, Node B){ return sqrt((A.xb - B.xb) * (A.xb - B.xb) + (A.yb - B.yb) * (A.yb - B.yb)); }//第二个点之间的距离
//呃……我承认两个dis写的很麻烦……但是好像也没什么很简单的法子
inline Node Line_inter(Line A, Line B){//斜截式直线求交点(之前写的直接copy过来的)
    if (A.mark == B.mark && (A.mark == 1 || A.mark == 0) ) return (Node){0, 0, 0, 0, 0} ;
    if ((A.mark == 1 && B.mark == 0) || (A.mark == 0 && B.mark == 1)) return (Node){1, get_x(A, B), get_y(A, B), 0, 0} ;
    if (A.mark == 1 && B.mark == 2) return (Node){1, (A.y - B.b) / B.k, A.y, 0, 0} ;
    if (A.mark == 2 && B.mark == 1) return (Node){1, (B.y - A.b) / A.k, B.y, 0, 0} ;
    if (A.mark == 2 && B.mark == 0) return (Node){1, B.x, B.x * A.k + A.b, 0, 0} ;
    if (A.mark == 0 && B.mark == 2) return (Node){1, A.x, A.x * B.k + B.b, 0, 0} ;
    return (Node){1, (A.b - B.b) / (B.k - A.k), (A.b - B.b) / (B.k - A.k) * A.k + A.b, 0, 0} ;  
}
inline Node get_inter (Circle A, Circle B){//“生成曲线”求交点
    if ((A.mark == 0 && B.mark) || (A.mark && B.mark == 0)){//一条是直线,一个是圆
        if (!A.mark) {Circle C ; C = A, A = B, B = C ;} // B is a line ;
        double a = 1 + (B.C / B.D) * (B.C / B.D), del ;
        double c = A.E - B.E * A.D / B.D + B.E * B.E /((B.D) * (B.D)) ;
        double b = (A.C - B.C * A.D / B.D + 2 * B.C * B.E /((B.D) * (B.D)) ) ; 
        if ((del = (b * b - 4 * a * c)) < -Eps) return (Node){0, 0, 0, 0, 0} ; 
        // printf("%lf %lf %lf %lf\n", a, b, c, del) ;
        double xa =  (-b + sqrt(del)) / (2 * a), xb = (-b - sqrt(del)) / (2 * a) ;
        double ya = -B.C / B.D * xa - B.E / B.D, yb = -B.C / B.D * xb - B.E / B.D ; 
        // cout << "-----------------" << xa << " " << ya << " " << xb << " " << yb << endl ;
        return (Node){2, xa, ya, xb, yb} ;//此处由于误差等原因,不容易判断是否delta=0的情况,所如果delta=0直接记录两遍,不影响结果
    }
    if (!A.mark && !B.mark){
        Line La, Lb ; //两条都是直线,那么就直接转化成斜截式求。
        if (!A.C) La = (Line){1, 0, 0, 0, - A.E / A.D} ; 
        else if (!A.D) La = (Line){0, 0, 0, -A.E / A.C, 0} ; else La = (Line){2, -A.C / A.D, -A.E / A.D, 0, 0} ;
        if (!B.C) Lb = (Line){1, 0, 0, 0, - B.E / B.D} ; 
        else if (!B.D) Lb = (Line){0, 0, 0, -B.E / B.C, 0} ; else Lb = (Line){2, -B.C / B.D, -B.E / B.D, 0, 0} ;  
        return Line_inter(La, Lb) ; 
    } 
    if (A.mark && B.mark){
        Circle C = A - B ; return get_inter(C, A) ;
        //此处需要用到一点小知识,就是两个圆的交点很难求,但是我们可以通过相减求出交线来(必修二知识点),那么就直接把这条线代会第一个if里就好。
    }
}
inline Circle make_rat(Circle A, Circle B){//rat = ratio[n.]比例;比率,用来求生成曲线的函数
    double _k2 = (A.r / B.r) * (A.r / B.r) ; Circle Ans ; double t ; 
    Ans.A = Ans.B = (_k2 - 1), Ans.C = -2 * (_k2 * B.x - A.x), Ans.D = -2 * (_k2 * B.y - A.y), 
    Ans.E = (_k2 * B.x * B.x - A.x * A.x) + (_k2 * B.y * B.y - A.y * A.y), Ans.x = Ans.y = Ans.r = 0 ; 
    if (Ans.A != 0) Ans.mark = 1, t = Ans.A, Ans.A /= t, Ans.B /= t, Ans.C /= t, Ans.D /= t, Ans.E /= t ; else Ans.mark = 0 ; return Ans ;
}
inline void make_for_Ans(){//最后的结果,判断选哪个交点
    sort(I + 1, I + 3, Comp) ;//我闲的,方便一点
    if (I[1].mark <= 1) ansx = I[1].xa, ansy = I[1].ya ;
    else {
        double A1, A11, B1, B11 ;
        I[1] = get_inter(C[4], C[5]) ;
        A1 = disa(I[1], (Node){0, C[1].x, C[1].y, 0, 0}) / C[1].r ; 
        A11 = disa(I[1], (Node){0, C[3].x, C[3].y, 0, 0}) / C[3].r ;
        B1 = disb(I[1], (Node){0, 0, 0, C[1].x, C[1].y}) / C[1].r ; 
        B11 = disb(I[1], (Node){0, 0, 0, C[3].x, C[3].y}) / C[3].r ;
        if (equal(A1, A11) && !equal(B1, B11)) ansx = I[1].xa, ansy = I[1].ya ; 
        else if (!equal(A1, A11) && equal(B1, B11)) ansx = I[1].xb, ansy = I[1].yb ;
        else if (!equal(A1, A11) && !equal(B1, B11)) check = 1 ;//如果在误差范围内都不相等就说明无解。
        else {
            double Ja = sin(1 / A1), Jb = sin(1 / B1) ;//比较角的大小,通过sin来搞
            if (Ja > Jb) ansx = I[1].xa, ansy = I[1].ya ; else ansx = I[1].xb, ansy = I[1].yb ;
        }
    }
}
int main(){
    for (i = 1 ; i <= 3 ; ++ i) cin >> C[i].x >> C[i].y >> C[i].r ;
    C[4] = make_rat(C[1], C[2]), C[5] = make_rat(C[2], C[3]), C[6] = make_rat(C[3], C[1]), 
    I[1] = get_inter(C[4], C[5]), I[2] = get_inter(C[5], C[6]), I[3] = get_inter(C[4], C[6]) ; 
    /*cout << I[1].xa << " " << I[1].xb << " " << I[1].ya << " " << I[1].yb << " " << I[1].mark << endl ;
    cout << I[2].xa << " " << I[2].xb << " " << I[2].ya << " " << I[2].yb << " " << I[2].mark << endl ;
    cout << I[3].xa << " " << I[3].xb << " " << I[3].ya << " " << I[3].yb << " " << I[3].mark << endl ;*/
    if (!I[1].mark || !I[2].mark || !I[3].mark) return putchar('\n'), 0 ; make_for_Ans() ; (!check) ? printf("%.5lf %.5lf", ansx, ansy) : 1 ; return 0 ; 
}

PS:

  • 这道题我写了好久…大概3h?因为边做变走神+公式难推……
  • 233此处向我们最可爱的数学老师lxa致以敬意qaq……前不久刚学的必修二就用上了qwq