@codeforces - 607E@ Cross Sum


@description@

给定 n 条直线,构造可重点集合 I 为 n 条直线两两的交点集合。

再给定一点 (p, q)。求该点到 I 中的点的距离中前 m 小的距离之和。

原题链接。

@solution@

求前 m 小,而 m 又这么大,不难想到二分答案。

二分答案 r,可以得到一个以 (p, q) 为圆心、以 r 为半径的圆。
我们想要求圆内直线的交点个数 cnt 并判断是否有 cnt >= m。

注意到此时直线在圆中呈现弦的形式。而圆中弦的相交,可以通过将端点转成极角序,然后判断两个区间严格相交(没有包含的情况)来实现。这个直接树状数组随便做。
然而上面只能解决严格在圆内的交点,恰好在圆上的交点,需要判断一个端点连接了多少条不同的直线来求解。

最后我们得到第 m 大的距离在 [l, l + eps] 之间,且这个区间之间只含这一种距离。
怎么求前 m 大之和呢?我们可以一一找出以 l 为半径的这个圆中的直线交点。先找出严格在圆内的:直接作一个二维偏序即可。
此时在圆上的交点个数可能非常多,所以我们不能取 (l + eps) 为半径,不然就会把这些交点算进严格在圆内的。
圆上的交点个数 = (m - 严格在圆内的交点个数)。

@accepted code@

#include <cmath>
#include <cstdio>
#include <vector>
#include <iomanip>
#include <iostream>
#include <algorithm>
using namespace std;

#define double long double

const int MAXN = 100000;
const double EPS = 1E-16;
const double PI = acos(-1);

int dcmp(double x) {
    return (fabs(x) <= EPS ? 0 : (x < 0 ? -1 : 1));
}

int n, m; double K[MAXN + 5], B[MAXN + 5];

double t[MAXN + 5]; int cnt;
bool cmp(double a, double b) {return dcmp(a - b) < 0;}
bool cmp2(double a, double b) {return dcmp(a - b) == 0;}

int tree[MAXN + 5], le[MAXN + 5], ri[MAXN + 5];
int lowbit(int x) {
    return x & -x;
}
int sum(int x) {
    int ret = 0;
    for(int i=x;i;i-=lowbit(i))
        ret += tree[i];
    return ret;
}
void update(int x, int k) {
    for(int i=x;i<=cnt;i+=lowbit(i))
        tree[i] += k;
}

int tmp[MAXN + 5];
vector<int>v1[MAXN + 5], v2[MAXN + 5];
bool check(double r) {
    cnt = 0;
    for(int i=1;i<=n;i++) {
        double del = K[i]*K[i]*r*r + r*r - B[i]*B[i];
        if( dcmp(del) >= 0 ) {
            double frac = K[i]*K[i] + 1; del = sqrt(del);
            double x1 = -(del + B[i]*K[i])/frac, y1 = -(K[i]*del - B[i])/frac;
            t[++cnt] = atan2(y1, x1);
            double x2 = -(-del + B[i]*K[i])/frac, y2 = (K[i]*del + B[i])/frac;
            t[++cnt] = atan2(y2, x2);
        }
    }
    sort(t + 1, t + cnt + 1, cmp);
    cnt = unique(t + 1, t + cnt + 1, cmp2) - t - 1;
    for(int i=1;i<=cnt;i++)
        tree[i] = tmp[i] = 0, v1[i].clear(), v2[i].clear();
    for(int i=1;i<=n;i++) {
        double del = K[i]*K[i]*r*r + r*r - B[i]*B[i];
        if( dcmp(del) >= 0 ) {
            double frac = K[i]*K[i] + 1; del = sqrt(del);
            double x1 = -(del + B[i]*K[i])/frac, y1 = -(K[i]*del - B[i])/frac;
            int p = lower_bound(t + 1, t + cnt + 1, atan2(y1, x1), cmp) - t;
            double x2 = -(-del + B[i]*K[i])/frac, y2 = (K[i]*del + B[i])/frac;
            int q = lower_bound(t + 1, t + cnt + 1, atan2(y2, x2), cmp) - t;
            le[i] = min(p, q), ri[i] = max(p, q);
            v1[le[i]].push_back(i), v2[ri[i]].push_back(i);
            tmp[le[i]]++;
            if( le[i] != ri[i] ) tmp[ri[i]]++;
        }
    }
    int tot = 0;
    for(int i=1;i<=cnt;i++) {
        for(int j=0;j<v1[i].size();j++) update(i, 1);
        for(int j=0;j<v2[i].size();j++)
            update(le[v2[i][j]], -1), tot += sum(i-1) - sum(le[v2[i][j]]);
        tot += 1LL*tmp[i]*(tmp[i] - 1)/2;
    }
    return tot >= m;
}

struct node{
    int le, ri, id;
    friend bool operator < (node a, node b) {return a.ri < b.ri;}
}b[MAXN + 5], c[MAXN + 5];
double ans; int nw;

double intersect(int i, int j) {
    double x = (B[j] - B[i]) / (K[i] - K[j]);
    double y = K[i] * x + B[i];
    return sqrt(x*x + y*y);
}

void solve(int le, int ri) {
    if( le == ri ) return ;
    int m = (le + ri) >> 1;
    solve(le, m), solve(m + 1, ri);
    
    int p = le, q = m + 1, t = le;
    while( p <= m && q <= ri ) {
        if( b[p].le < b[q].le ) {
            for(int i=q;i<=ri;i++) {
                if( b[i].le >= b[p].ri ) break;
                if( b[p].ri < b[i].ri ) nw++, ans += intersect(b[p].id, b[i].id);
            }
            c[t++] = b[p++];
        }
        else c[t++] = b[q++];
    }
    while( p <= m ) c[t++] = b[p++];
    while( q <= ri ) c[t++] = b[q++];
    for(int i=le;i<=ri;i++) b[i] = c[i];
}

double get(double r) {
    cnt = 0;
    for(int i=1;i<=n;i++) {
        double del = K[i]*K[i]*r*r + r*r - B[i]*B[i];
        if( dcmp(del) >= 0 ) {
            double frac = K[i]*K[i] + 1; del = sqrt(del);
            double x1 = -(del + B[i]*K[i])/frac, y1 = -(K[i]*del - B[i])/frac;
            t[++cnt] = atan2(y1, x1);
            double x2 = -(-del + B[i]*K[i])/frac, y2 = (K[i]*del + B[i])/frac;
            t[++cnt] = atan2(y2, x2);
        }
    }
    sort(t + 1, t + cnt + 1, cmp);
    cnt = unique(t + 1, t + cnt + 1, cmp2) - t - 1;
    int tot = 0;
    for(int i=1;i<=n;i++) {
        double del = K[i]*K[i]*r*r + r*r - B[i]*B[i];
        if( dcmp(del) >= 0 ) {
            double frac = K[i]*K[i] + 1; del = sqrt(del);
            double x1 = -(del + B[i]*K[i])/frac, y1 = -(K[i]*del - B[i])/frac;
            int p = lower_bound(t + 1, t + cnt + 1, atan2(y1, x1), cmp) - t;
            double x2 = -(-del + B[i]*K[i])/frac, y2 = (K[i]*del + B[i])/frac;
            int q = lower_bound(t + 1, t + cnt + 1, atan2(y2, x2), cmp) - t;
            tot++, b[tot].le = min(p, q), b[tot].ri = max(p, q), b[tot].id = i;
        }
    }
    sort(b + 1, b + tot + 1), solve(1, tot);
    return ans + 1.0*(m - nw)*r;
}

int read() {
    int x = 0, f = 1; char ch = getchar();
    while( (ch != '-') && (ch > '9' || ch < '0') ) ch = getchar();
    if( ch == '-' ) f = -1, ch = getchar();
    while( '0' <= ch && ch <= '9' ) x = 10*x + ch - '0', ch = getchar();
    return x * f;
}

int main() {
    double x, y;
    n = read(), x = 1.0*read()/1000, y = 1.0*read()/1000, m = read();
    for(int i=1;i<=n;i++) {
        K[i] = 1.0*read()/1000, B[i] = 1.0*read()/1000;
        B[i] = K[i]*x + B[i] - y;
    }
/*
    for(int i=1;i<=n;i++) {
        for(int j=i+1;j<=n;j++)
            printf("%d %d : ", i, j), cout << dist(intersect(l[i], l[j]), pnt)  << endl;
        cout << i << " " << (pnt - a[i]) * l[i].ab << endl;
    }
*/
    double le = 0, ri = 2E9 + 5000;
    for(int i=0;i<=50;i++) {
        double mid = (le + ri) / 2;
        if( check(mid) ) ri = mid;
        else le = mid;
    }
    cout << fixed << setprecision(9) << get(le) << endl;
}

@details@

本来思路很简单,结果它先是卡我精度然后又卡我常。。。

写几点做题时遇到的问题:
(1)可以先把直线平移,并把 (p, q) 平移成坐标系原点
(2)不要觉得向量表达就是万能的,这道题恰巧向量表达精度损失非常大。。。反而用直线与圆的方程联立来得更直接。
(3)二分次数不要太多小心 TLE。cin 是真的慢,能换掉就换。

posted @ 2020-01-10 20:59  Tiw_Air_OAO  阅读(...)  评论(...编辑  收藏