LA 4973 Ardenia (3D Geometry + Simulation)

ACM-ICPC Live Archive

　　三维几何，题意是要求求出两条空间线段的距离。题目难度在于要求用有理数的形式输出，这就要求写一个有理数类了。

　　开始的时候写出来的有理数类就各种疯狂乱套，TLE的结果是显然的。后来发现，在计算距离前都是不用用到有理数类的，所以就将开始的部分有理数改成直接用long long。其实好像可以用int来做的，不过我的方法比较残暴，中间运算过程居然爆int了。所以就只好用long long了。

代码如下，附带debug以及各种强的数据：

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

using namespace std;

template<class T> T gcd(T a, T b) { return b ? gcd(b, a % b) : a;}
template <class T> T sqr(T x) { return x * x;}
typedef long long LL;
struct Rat {
    LL a, b;
    Rat() {}
    Rat(LL x) { a = x, b = 1;}
    Rat(LL _a, LL _b) {
        LL GCD = gcd(_b, _a);
        a = _a / GCD, b = _b / GCD;
    }
    double val() { return (double) a / b;}
    bool operator < (Rat x) const { return a * x.b < b * x.a;}
    bool operator > (Rat x) const { return a * x.b > b * x.a;}
    bool operator == (Rat x) const { return a * x.b == b * x.a;}
    bool operator < (LL x) const { return Rat(a, b) < Rat(x);}
    bool operator > (LL x) const { return Rat(a, b) > Rat(x);}
    bool operator == (LL x) const { return Rat(a, b) == Rat(x);}
    Rat operator + (Rat x) {
        LL tb = b * x.b, ta = a * x.b + b * x.a;
        LL GCD = gcd(abs(tb), abs(ta));
        return Rat(ta / GCD, tb / GCD);
    }
    Rat operator - (Rat x) {
        LL tb = b * x.b, ta = a * x.b - b * x.a;
        LL GCD = gcd(abs(tb), abs(ta));
        return Rat(ta / GCD, tb / GCD);
    }
    Rat operator * (Rat x) {
        if (a * x.a == 0) return Rat(0, 1);
//        if (b * x.b == 0) { puts("..."); while (1) ;}
        LL tb = b * x.b, ta = a * x.a;
        LL GCD = gcd(abs(tb), abs(ta));
        return Rat(ta / GCD, tb / GCD);
    }
    Rat operator / (Rat x) {
        if (a * x.b == 0) return Rat(0, 1);
//        if (b * x.a == 0) { puts("!!!"); while (1) ;}
        LL GCD, tb = b * x.a, ta = a * x.b;
        GCD = gcd(abs(tb), abs(ta));
        return Rat(ta / GCD, tb / GCD);
    }
    void fix() {
        a = abs(a), b = abs(b);
        LL GCD = gcd(b, a);
        a /= GCD, b /= GCD;
    }
} ;

struct Point {
    LL x[3];
    Point operator + (Point a) {
        Point ret;
        for (int i = 0; i < 3; i++) ret.x[i] = x[i] + a.x[i];
        return ret;
    }
    Point operator - (Point a) {
        Point ret;
        for (int i = 0; i < 3; i++) ret.x[i] = x[i] - a.x[i];
        return ret;
    }
    Point operator * (LL p) {
        Point ret;
        for (int i = 0; i < 3; i++) ret.x[i] = x[i] * p;
        return ret;
    }
    bool operator == (Point a) const {
        for (int i = 0; i < 3; i++) if (!(x[i] == a.x[i])) return false;
        return true;
    }
    void print() {
        for (int i = 0; i < 3; i++) cout << x[i] << ' ';
        cout << endl;
    }
} ;
typedef Point Vec;

struct Line {
    Point s, t;
    Line() {}
    Line (Point s, Point t) : s(s), t(t) {}
    Vec vec() { return t - s;}
} ;
typedef Line Seg;

LL dotDet(Vec a, Vec b) {
    LL ret = 0;
    for (int i = 0; i < 3; i++) ret = ret + a.x[i] * b.x[i];
    return ret;
}

Vec crossDet(Vec a, Vec b) {
    Vec ret;
    for (int i = 0; i < 3; i++) {
        ret.x[i] = a.x[(i + 1) % 3] * b.x[(i + 2) % 3] - a.x[(i + 2) % 3] * b.x[(i + 1) % 3];
    }
    return ret;
}

inline LL vecLen(Vec x) { return dotDet(x, x);}
inline bool parallel(Line a, Line b) { return vecLen(crossDet(a.vec(), b.vec())) == 0;}
inline bool onSeg(Point x, Point a, Point b) { return parallel(Line(a, x), Line(b, x)) && dotDet(a - x, b - x) < 0;}
inline bool onSeg(Point x, Seg s) { return onSeg(x, s.s, s.t);}

Rat pt2Seg(Point p, Point a, Point b) {
    if (a == b) return Rat(vecLen(p - a));
    Vec v1 = b - a, v2 = p - a, v3 = p - b;
    if (dotDet(v1, v2) < 0) return Rat(vecLen(v2));
    else if (dotDet(v1, v3) > 0) return Rat(vecLen(v3));
    else return Rat(vecLen(crossDet(v1, v2)), vecLen(v1));
}
inline Rat pt2Seg(Point p, Seg s) { return pt2Seg(p, s.s, s.t);}
inline Rat pt2Plane(Point p, Point p0, Vec n) { return Rat(sqr(dotDet(p - p0, n)), vecLen(n));}
inline bool segIntersect(Line a, Line b) {
    Vec v1 = crossDet(a.s - b.s, a.t - b.s);
    Vec v2 = crossDet(a.s - b.t, a.t - b.t);
    Vec v3 = crossDet(b.s - a.s, b.t - a.s);
    Vec v4 = crossDet(b.s - a.t, b.t - a.t);
//    v1.print();
//    v2.print();
//    cout << dotDet(v1, v2).val() << "= =" << endl;
    return dotDet(v1, v2) < 0 && dotDet(v3, v4) < 0;
}//        cout << "same plane" << endl;

pair<Rat, Rat> getIntersect(Line a, Line b) {
    Point p = a.s, q = b.s;
    Vec v = a.vec(), u = b.vec();
    LL uv = dotDet(u, v), vv = dotDet(v, v), uu = dotDet(u, u);
    LL pv = dotDet(p, v), qv = dotDet(q, v), pu = dotDet(p, u), qu = dotDet(q, u);
    if (uv == 0) return make_pair(Rat(qv - pv, vv), Rat(pu - qu, uu));
//    if (vv == 0 || uv == 0 || uv / vv - uu / uv == 0) { puts("shit!"); while (1) ;}
    Rat y = (Rat(pv - qv, vv) - Rat(pu - qu, uv)) / (Rat(uv, vv) - Rat(uu, uv));
    Rat x = (y * uv - pv + qv) / vv;
//    cout << x.a << ' ' << x.b << ' ' << y.a << ' ' << y.b << endl;
    return make_pair(x, y);
}

void work(Point *pt) {
    Line a = Line(pt[0], pt[1]);
    Line b = Line(pt[2], pt[3]);
    if (parallel(a, b)) {
        if (onSeg(pt[0], b) || onSeg(pt[1], b)) { puts("0 1"); return ;}
        if (onSeg(pt[2], a) || onSeg(pt[3], a)) { puts("0 1"); return ;}
//        cout << "parallel" << endl;
        Rat tmp = min(min(pt2Seg(pt[0], b), pt2Seg(pt[1], b)), min(pt2Seg(pt[2], a), pt2Seg(pt[3], a)));
        tmp.fix();
        printf("%lld %lld
", tmp.a, tmp.b);
        return ;
    }
    Vec nor = crossDet(a.vec(), b.vec());
    Rat ans = pt2Plane(pt[0], pt[2], nor);
//    cout << "~~~" << endl;
    if (ans == 0) {
//        cout << "same plane" << endl;
        if (segIntersect(a, b)) { puts("0 1"); return ;}
        Rat tmp = min(min(pt2Seg(pt[0], b), pt2Seg(pt[1], b)), min(pt2Seg(pt[2], a), pt2Seg(pt[3], a)));
        tmp.fix();
        printf("%lld %lld
", tmp.a, tmp.b);
        return ;
    } else {
//        cout << "diff plane" << endl;
        pair<Rat, Rat> tmp = getIntersect(a, b);
//        cout << tmp.first.val() << "= =" << tmp.second.val() << endl;
//        cout << (tmp.first > 0) << endl;
        if (tmp.first > 0 && tmp.first < 1 && tmp.second > 0 && tmp.second < 1) {
//            cout << "cross" << endl;
            ans.fix();
            printf("%lld %lld
", ans.a, ans.b);
        } else {
//            cout << "not cross" << endl;
            Rat t = min(min(pt2Seg(pt[0], b), pt2Seg(pt[1], b)), min(pt2Seg(pt[2], a), pt2Seg(pt[3], a)));
            t.fix();
            printf("%lld %lld
", t.a, t.b);
        }
    }
}

int main() {
//    freopen("in", "r", stdin);
//    freopen("out", "w", stdout);
    Point pt[4];
    int T;
    cin >> T;
    while (T--) {
        for (int i = 0; i < 4; i++) {
            for (int j = 0; j < 3; j++) {
                scanf("%lld", &pt[i].x[j]);
            }
        }
        work(pt);
    }
    return 0;
}

/*
13

-20 -20 -20 20 20 19
0 0 0 1 1 1

-20 -20 -20 20 19 20
-20 -20 20 20 20 -20

0 0 0 20 20 20
0 0 10 0 20 10

0 0 0 1 1 1
2 3 4 1 2 2

0 0 0 0 0 0
0 1 1 1 2 3

0 0 0 10 10 10
11 12 13 10 11 11

0 0 0 1 1 1
1 1 1 2 2 2

1 0 0 0 1 0
1 1 0 2 2 0

1 0 0 0 1 0
0 0 0 1 1 0

0 0 0 0 0 20
20 0 10 0 20 10

0 0 0 20 20 20
1 1 2 1 1 2

0 0 0 20 20 20
0 20 20 0 20 20

0 0 0 0 0 20
20 20 0 20 20 20
*/

——written by Lyon