[代码]ural 1913 Titan Ruins: Old Generators Are Fine Too

Abstract

ural 1913 Titan Ruins: Old Generators Are Fine Too

implementation 几何 圆交

Source

http://acm.timus.ru/problem.aspx?space=1&num=1913

Solution

主要思路就是判断距离到某两个点是r，另一个点是2r的点。

比赛时候用三分乱搞，大概是被卡了精度还是什么的，没有通过。

赛后用圆交重写了一遍，一开始还是没通过。

然后发现自己的圆交模板是针对面积问题的，单点的情况没有考虑清楚，乱改的时候改错了。大致就是处理圆覆盖的情况时，在预处理统计时算被覆盖了一次，在算两圆交点的时候又被统计了一次。还是没弄熟悉圆交的原理所致。

Code

#include <iostream>
#include <cassert>
#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

const double EPS = 1e-8;
const double PI = acos(-1.0);

inline int sgn(double x) {
    if (fabs(x)<EPS) return 0;
    return x>0?1:-1;
}

struct sp {
    double x, y;
    double pa;
    int cnt;
    sp() {}
    sp(double a, double b): x(a), y(b) {}
    sp(double a, double b, double c, int d): x(a), y(b), pa(c), cnt(d) {}
    bool operator<(const sp &rhs) const {
        if (sgn(pa-rhs.pa)==0) return cnt>rhs.cnt;
        return pa<rhs.pa;
    }
    void read() {scanf("%lf%lf", &x, &y);}
    void write() {printf("%.10f %.10f\n", x, y);}
}t, s1, s2, p[5], e[20];
double r, rad[5];
int cover[5];

sp operator*(double d, const sp &v) {
    return sp(d*v.x, d*v.y);
}

sp operator-(const sp &u, const sp &v) {
    return sp(u.x-v.x, u.y-v.y);
}

sp operator+(const sp &u, const sp &v) {
    return sp(u.x+v.x, u.y+v.y);
}

double dot(const sp &u, const sp &v) {
    return u.x*v.x+u.y*v.y;
}

double det(const sp &u, const sp &v) {
    return u.x*v.y-u.y*v.x;
}

double dis(const sp &u, const sp &v) {
    double dx = u.x-v.x;
    double dy = u.y-v.y;
    return sqrt(dx*dx+dy*dy);
}

double dissqr(const sp &u, const sp &v) {
    double dx = u.x-v.x;
    double dy = u.y-v.y;
    return dx*dx+dy*dy;
}

void write(sp u, sp v) {
    puts("Now we have enough power");
    u.write(); v.write();
}

bool cirint(const sp &u, double ru, const sp &v, double rv, sp &a, sp &b) {
    double d = dis(u, v);
    if (sgn(d-(ru+rv))>0 || sgn(d-fabs(ru-rv))<=0) return 0;
    sp c = u-v;
    double ca, sa, cb, sb, csum, ssum;

    ca = c.x/d, sa = c.y/d;
    cb = (rv*rv+d*d-ru*ru)/(2*rv*d), sb = sqrt(1-cb*cb);
    csum = ca*cb-sa*sb;
    ssum = sa*cb+sb*ca;
    a = sp(rv*csum, rv*ssum);
    a = a+v;

    sb = -sb;
    csum = ca*cb-sa*sb;
    ssum = sa*cb+sb*ca;
    b = sp(rv*csum, rv*ssum);
    b = b+v;

    return 1;
}

bool cu(int N, sp p[], double r[], sp &res) {
    int i, j, k;
    memset(cover, 0, sizeof cover);
    for (i = 0; i < N; ++i)
        for (j = 0; j < N; ++j) {
            double rd = r[i]-r[j];
            if (i!=j && sgn(rd)>0 && sgn(dis(p[i], p[j])-rd)<=0)
                cover[j]++;
        }
    sp a, b;
    for (i = 0; i < N; ++i) {
        int ecnt = 0;
        e[ecnt++] = sp(p[i].x-r[i], p[i].y, -PI, 1);
        e[ecnt++] = sp(p[i].x-r[i], p[i].y, PI, -1);
        for (j = 0; j < N; ++j) {
            if (j==i) continue;
            if (cirint(p[i], r[i], p[j], r[j], a, b)) {
                e[ecnt++] = sp(a.x, a.y, atan2(a.y-p[i].y, a.x-p[i].x), 1);
                e[ecnt++] = sp(b.x, b.y, atan2(b.y-p[i].y, b.x-p[i].x), -1);
                if (sgn(e[ecnt-2].pa-e[ecnt-1].pa)>0) {
                    e[0].cnt++;
                    e[1].cnt--;
                }
            }
        }
        sort(e, e+ecnt);
        int cnt = e[0].cnt;
        for (j = 1; j < ecnt; ++j) {
            if (cover[i]+cnt==N) {
                double a = (e[j].pa+e[j-1].pa)/2;
                res = p[i]+sp(r[i]*cos(a), r[i]*sin(a));
                return 1;
            }
            cnt += e[j].cnt;
        }
    }
    return 0;
}

bool solve(const sp &o, const sp &u, const sp &v, double r) {
    p[0] = o, rad[0] = r+r;
    p[1] = u; p[2] = v; rad[1] = rad[2] = r;
    sp a, b;
    if (cu(3, p, rad, a)) {
        b = 0.5*(a-o);
        b = o+b;
        write(a, b);
        return 1;
    }
    else return 0;
}

int main() {
    cin>>r;
    t.read(); s1.read(); s2.read();
    sp a, b, c, d;
    if (cirint(s1, r, s2, r, a, b)) {
        if (solve(t, s1, s2, r)) return 0;
    }
    else {
        if (cirint(s1, r, t, r, a, b) && cirint(s2, r, t, r, c, d)) {
            write(a, c);
            return 0;
        }
        else {
            if (solve(s1, t, s2, r)) return 0;
            if (solve(s2, t, s1, r)) return 0;
        }
    }
    puts("Death");
    return 0;
}

顺便贴圆交模板，这回的应该没问题了。

#include <cmath>
#include <algorithm>
using namespace std;

const double PI = acos(-1.0);
const double EPS = 1e-8;

int sgn(double d) {
    if (fabs(d)<EPS) return 0;
    return d>0?1:-1;
}

struct sp {
    double x, y;
    double pa;
    int cnt;
    sp() {}
    sp(double a, double b): x(a), y(b) {}
    sp(double a, double b, double c, int d): x(a), y(b), pa(c), cnt(d) {}
    bool operator<(const sp &rhs) const {
        /* 需要处理单点的情况时
        if (sgn(pa-rhs.pa)==0) return cnt>rhs.cnt;
        */
        return pa<rhs.pa;
    }
    void read() {scanf("%lf%lf", &x, &y);}
    void write() {printf("%lf %lf\n", x, y);}
};

sp operator+(const sp &u, const sp &v) {
    return sp(u.x+v.x, u.y+v.y);
}

sp operator-(const sp &u, const sp &v) {
    return sp(u.x-v.x, u.y-v.y);
}

double det(const sp &u, const sp &v) {
    return u.x*v.y-v.x*u.y;
}

double dir(const sp &p, const sp &u, const sp &v) {
    return det(u-p, v-p);
}

double dis(const sp &u, const sp &v) {
    double dx = u.x-v.x;
    double dy = u.y-v.y;
    return sqrt(dx*dx+dy*dy);
}

//计算两圆交点
//输入圆心坐标和半径
//返回两圆是否相交(相切不算)
//如果相交，交点返回在a和b(从a到b为u的交弧)
bool cirint(const sp &u, double ru, const sp &v, double rv, sp &a, sp &b) {
    double d = dis(u, v);
    /* 需要处理单点的情况时 覆盖的情况在这里不统计
    if (sgn(d-(ru+rv))>0 || sgn(d-fabs(ru-rv))<=0) return 0;
    */
    if (sgn(d-(ru+rv))>=0 || sgn(d-fabs(ru-rv))<=0) return 0;
    sp c = u-v;
    double ca, sa, cb, sb, csum, ssum;

    ca = c.x/d, sa = c.y/d;
    cb = (rv*rv+d*d-ru*ru)/(2*rv*d), sb = sqrt(1-cb*cb);
    csum = ca*cb-sa*sb;
    ssum = sa*cb+sb*ca;
    a = sp(rv*csum, rv*ssum);
    a = a+v;

    sb = -sb;
    csum = ca*cb-sa*sb;
    ssum = sa*cb+sb*ca;
    b = sp(rv*csum, rv*ssum);
    b = b+v;

    return 1;
}

sp e[222];
int cover[111];

//求圆并
//输入点数N，圆心数组p，半径数组r，答案数组s
//s[i]表示被至少i个圆覆盖的面积(最普通的圆并就是s[1])
void circle_union(int N, sp p[], double r[], double s[]) {
    int i, j, k;
    memset(cover, 0, sizeof cover);
    for (i = 0; i < N; ++i)
        for (j = 0; j < N; ++j) {
            double rd = r[i]-r[j];
            //覆盖的情况在这里统计
            if (i!=j && sgn(rd)>0 && sgn(dis(p[i], p[j])-rd)<=0)
                cover[j]++;
        }
    for (i = 0; i < N; ++i) {
        int ecnt = 0;
        e[ecnt++] = sp(p[i].x-r[i], p[i].y, -PI, 1);
        e[ecnt++] = sp(p[i].x-r[i], p[i].y, PI, -1);
        for (j = 0; j < N; ++j) {
            if (j==i) continue;
            if (cirint(p[i], r[i], p[j], r[j], a, b)) {
                e[ecnt++] = sp(a.x, a.y, atan2(a.y-p[i].y, a.x-p[i].x), 1);
                e[ecnt++] = sp(b.x, b.y, atan2(b.y-p[i].y, b.x-p[i].x), -1);
                if (sgn(e[ecnt-2].pa-e[ecnt-1].pa)>0) {
                    e[0].cnt++;
                    e[1].cnt--;
                }
            }
        }
        sort(e, e+ecnt);
        int cnt = e[0].cnt;
        for (j = 1; j < ecnt; ++j) {
            double pad = e[j].pa-e[j-1].pa;
            s[cover[i]+cnt] += (det(e[j-1], e[j])+r[i]*r[i]*(pad-sin(pad)))/2;
            cnt += e[j].cnt;
        }
    }
}