HDU5572(2015ACM/ICPC亚洲区上海站A题）_计算几何（积累计算几何的板子）

题目链接：An Easy Physics Problem
题意：给A、B两点，方向向量V和圆柱体的圆心和半径。问题是：在一个光滑的平面上，给出大小可忽略不计的小球，起始点是A，沿方向向量V给一个速度，问小球会不会经过B这个点。如果小球撞到圆柱体，小球会发生完全弹性碰撞。
题解：①小球不能和圆柱体碰撞，那么判断小球的运动路径上经不经过B就行。
②小球和圆柱碰撞，在撞之前经过B点
③小球和圆柱碰撞，在撞之后经过B点。求出B点的对称点，看AB是否在一条直线上。

给出两个代码，网上借鉴了一下：

法①：

#include <iostream>
#include <cstdio>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <deque>
#include <vector>
#include <queue>
#include <string>
#include <cstring>
#include <map>
#include <stack>
#include <set>
#include <bitset>
#define Max(a,b) a>b?a:b
#define Min(a,b) a>b?b:a
#define mem(a,b) memset(a,b,sizeof(a))
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
int dir[4][2]= {{1,0},{-1,0},{0,1},{0,-1}};
const double eps = 1e-8;
const double inf = 1e20;
const double pi = acos(-1.0);
const int INF=0x3f3f3f3f;//1061109567

inline int dcmp(double x){
    if(fabs(x) < eps) return 0;
    return x < 0? -1: 1;
}

struct Point{
    double x,y;
    Point(double x = 0, double y = 0): x(x), y(y){}
    void read(){
        scanf("%lf%lf",&x,&y);
    }
    Point operator + (const Point &u) const{ return Point(x+u.x, y+u.y); }
    Point operator - (const Point &u) const{ return Point(x-u.x, y-u.y); }
    Point operator * (const double &k) const{ return Point(x*k, y*k); }
};

typedef Point Vector;

struct Circle{
    Point o;
    double r;
    void read(){
        scanf("%lf%lf%lf", &o.x, &o.y, &r);
    }
};

Circle C;
Point A,B;
Vector V;

double getCross(Vector a, Vector b){ return a.x*b.y - a.y*b.x; }
Vector getNormal(Vector a){ return Vector(-a.y, a.x); }

//已知点起始位置，初始方向向量和圆的位置，求线和圆的交点个数以及交点与起始点的横坐标或者纵坐标的值之差。
int getLineCircleIntersection(Point p, Vector v, Circle O, double &dis1, double &dis2){
    double a = v.x, b = p.x - O.o.x, c = v.y, d = p.y - O.o.y;
    double e = a*a + c*c, f = 2*(a*b+c*d), g = b*b+d*d-O.r*O.r;
    double delta = f*f - 4*e*g;
    if(dcmp(delta) < 0) return 0;
    if(dcmp(delta) == 0){
        dis1 = dis2 = -f/(2*e);
        return 1;
    }
    dis1 = (-f - sqrt(delta)) / (2*e);
    dis2 = (-f + sqrt(delta)) / (2*e);
    return 2;
}

bool onLine(Point a, Vector v, Point b){
    return dcmp(getCross(v,b-a)) == 0;//叉乘等于0，表示v向量和(b-a)向量平行
}

double getPos(Vector v, Vector t){//判断两个平行的向量是正向还是反向，并且把初始方向向量作为单位长度
    if(dcmp(v.x) == 0){
        return t.y/v.y;
    }
    return t.x/v.x;
}

bool getIntersection(Point p, Vector v, Point q, Vector w, Point &o){//得到T(即o)点，用于求B点的对称点
    if(dcmp(getCross(v,w)) == 0) return false;//v,w向量既平行又垂直这是不可能的
    Vector u = p - q;
    double k = getCross(w,u) / getCross(v,w);
    o = p + v * k;
    return true;
}

bool judge(){
    double dis1,dis2,hint = inf;

    int k = getLineCircleIntersection(A,V,C,dis1,dis2);//得到线和圆的交叉点
    //dis1 <= dis2 如果射线与圆的两个交点的x相等，dis是A点的纵坐标与两交点的纵坐标之差
                //否则，dis是A点的横坐标与两交点的横坐标之差


    if(k > 1 && dcmp(dis1) >= 0) hint = dis1;
    else if(k > 1 && dcmp(dis2) >= 0) hint = dis2;

    //如果不产生碰撞，hint是无穷大

//    printf("x1:%f x2:%f\n",dis1,dis2);
//    printf("hint:%f\n",hint);

    if(onLine(A,V,B)){
        double t = getPos(V, B-A);
        if(dcmp(t) >= 0 && t < hint) return true;//t>=0保证向量A->B同向，判断到撞之前有没有经过B点
    }
    if(dcmp(hint - inf) < 0){
        Point I = A + V*hint, T;//得到交点坐标
        Vector H = getNormal(C.o - I);//向量I->C是半径方向指向圆心的法向量，
                                      //H是该法线的法向量即圆在点H处的切线向量
        getIntersection(I,H,B,(C.o - I),T);
        B = T * 2 - B;
        if(onLine(A,V,B)){
            double t = getPos(V, B - A);
            if(dcmp(t) >= 0)
                return true;//同向就行
        }

    }

    return false;
}

int main(){
    int cas;
    scanf("%d",&cas);
    for(int i = 1; i <= cas; i++){
        C.read(), A.read(), V.read(), B.read();

        printf("Case #%d: %s\n",i,judge()?"Yes":"No");
    }

    return 0;
}

法②：

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

const double eps = 1e-8;

int cmp(double x){
    return x < -eps ? -1 : (x>eps);
}

double pow2(double x){
    return x * x;
}

double mySqrt(double x){
    return sqrt(max((double)0, x));
}

struct Vec
{
    double x, y;
    Vec(){}
    Vec(double xx, double yy):x(xx), y(yy){}

    Vec& operator=(const Vec& v){
        x = v.x;
        y = v.y;
        return *this;
    }

    double norm(){
        return sqrt(pow2(x) + pow2(y));
    }
    //返回单位向量
    Vec unit(){
        return Vec(x, y) / norm();
    }
    //判等
    friend bool operator==(const Vec& v1, const Vec& v2){
        return cmp(v1.x - v2.x)==0 && cmp(v1.y - v2.y)==0;
    }

    //+, -, 数乘
    friend Vec operator+(const Vec& v1, const Vec& v2){
        return Vec(v1.x + v2.x, v1.y + v2.y);
    }
    friend Vec operator-(const Vec& v1, const Vec& v2){
        return Vec(v1.x - v2.x, v1.y - v2.y);
    }
    friend Vec operator*(const Vec& v, const double c){
        return Vec(c*v.x, c*v.y);
    }
    friend Vec operator*(const double c, const Vec& v){
        return Vec(c*v.x, c*v.y);
    }
    friend Vec operator/(const Vec& v, const double c){
        return Vec(v.x/c, v.y/c);
    }
};

int T;
int ans;
Vec O, A, V, B, C, B1;
int R;

//点乘
double dot(const Vec v1, const Vec v2){
    return v1.x*v2.x + v1.y*v2.y;
}
//叉乘的模长
double product(const Vec v1, const Vec v2){
    return v1.x*v2.y - v1.y*v2.x;
}

//点p到直线v1,v2的投影
Vec project(Vec& v1, Vec& v2, Vec& p){
    Vec v = v2 - v1;
    return v1 + v * dot(v, p-v1) / dot(v, v);
}
//点p关于直线v1,v2的对称点
Vec mirrorPoint(Vec& v1, Vec& v2, Vec& p){
    Vec c = project(v1, v2, p);
    //printf("project: %lf, %lf\n", c.x, c.y);
    return (double)2*c - p;
}

//判断点p是否在线段v1,v2上
bool onSeg(Vec& v1, Vec& v2, Vec& p){
    if(cmp(product(p-v1, v2-v1))==0 && cmp(dot(p-v1, p-v2))<=0)
        return true;
    return false;
}

bool calc_C(){
    //将AC表示为关于t的参数方程
    //则与圆的方程联立得到关于t的一元二次方程, a,b,c为一般式的三个系数
    double a = pow2(V.x) + pow2(V.y);
    double b = 2*V.x*(A.x - O.x) + 2*V.y*(A.y - O.y);
    double c = pow2(A.x - O.x) + pow2(A.y - O.y) - pow2(R);
    double delta = pow2(b) - 4*a*c;
    if(cmp(delta) <= 0) return false;
    else{
        double t1 = (-b - mySqrt(delta))/(2*a);
        double t2 = (-b + mySqrt(delta))/(2*a);
        double t;
        if(cmp(t1) >= 0) t = t1;
        if(cmp(t2) >= 0 && t2 < t1) t = t2;//取较小的那个，即为离A近的那个交点
        C.x = A.x + V.x*t;
        C.y = A.y + V.y*t;
        return true; //有交点
    }
}

int main()
{
    scanf("%d", &T);
    for(int ca = 1; ca <= T; ca++){
        scanf("%lf%lf%d", &O.x, &O.y, &R);
        scanf("%lf%lf%lf%lf", &A.x, &A.y, &V.x, &V.y);
        scanf("%lf%lf", &B.x, &B.y);
        if(calc_C()){
            if(onSeg(A, C, B)) ans = 1;
            else{
                //printf("has intersection (%lf, %lf)\n", C.x, C.y);
                Vec A1 = mirrorPoint(O, C, A);
                  // printf("A: %lf, %lf\n", A.x, A.y);
                  // printf("A1: %lf, %lf\n", A1.x, A1.y);
                //printf("B1 (%lf, %lf)\n", B1.x, B1.y);
                if(A==A1){
                    Vec u = B - O;
                    Vec v = C - O;
                    // printf("OB: %lf %lf\n", u.unit().x, u.unit().y);
                    // printf("OC: %lf %lf\n", v.unit().x, v.unit().y);
                    if(u.unit() == v.unit()){
                        ans = 1;
                    }else ans = 0;
                }else {
                    Vec u = B - C;
                    Vec v = A1 - C;
                    if(u.unit() == v.unit()){
                        ans = 1;
                    }
                    else ans = 0;
                }
            }
        }else{//运动方向不变，则AB与V同向才可碰到B
            //printf("no intersection\n");
            Vec temp = B - A;
            if(temp.unit() == V.unit())
                ans = 1;
            else ans = 0;
        }
        printf("Case #%d: %s\n", ca, ans ? "Yes" : "No");
    }
    return 0;
}