计算几何/平面和空间

Posted 2020-12-09 jaszzz

总结：

注意几点：

• 二维向量的操作可以直接使用STL中的complex
• 问题中的数值往往是浮点数，此时应该注意误差的问题，不考虑误差会WA掉的，这是非常重要的一点

1、计算几何基础

POJ 1127 Jack Straws

题意：判断给定的若干组线段是否有交点

#include<cstdio>
#include<cmath>
#include<iostream>
using namespace std;
const int maxn = 25;
bool g[maxn][maxn];
double eps = 1e-10;
double add(double a, double b)
{
    if (abs(a + b) < eps * (abs(a) + abs(b)))  return 0;
    else return a + b;
}

struct P {
    double x, y;
    P() {}
    P(double x, double y) : x(x), y(y) {}
    P operator + (P p) {
        return P(add(x, p.x), add(y, p.y));
    }
    P operator -(P p) {
        return P(add(x, -p.x), add(y, -p.y));
    }
    double dot(P p) {
        return add(x * p.x, y *p.y);
    }
    double det(P p) {
        return add(x * p.y, -y * p.x);
    }
    P operator *(double d) {
        return P(x * d, y * d);
    }
};
P p[maxn], q[maxn];
bool onseg(P p1, P p2, P q)
{
    return (p1 - q).det(p2 - q) == 0 && (p1 - q).dot(p2 - q) <= 0;
}
P intersection(P p1, P p2, P q1, P q2)
{
    return p1 + (p2 - p1) * ((q2 - q1).det(q1 - p1) / (q2 - q1).det(p2 - p1));
}
int main(void)
{
    int n;
    while (scanf("%d", &n) && n) {
        for (int i = 0; i < n; i++) {
            scanf("%lf%lf%lf%lf", &p[i].x, &p[i].y, &q[i].x, &q[i].y);
        }
        for (int i = 0; i < n; i++) {
            g[i][i] = true;
            for (int j = 0; j < i; j++) {
                //如果两个直线平行，只需要他们重合且线段有重合部分
                if ((p[i] - q[i]).det(p[j] - q[j]) == 0) {
                    g[i][j] = g[j][i] = onseg(p[i], q[i], p[j]) 
                        || onseg(p[j], q[j], p[i]) 
                        || onseg(p[i], q[i], q[j]) 
                        || onseg(p[j], q[j], q[i]);
                }
                else {
                    //不平行的话先求交点，然后再判断交点是否在两条线段上
                    P tmp = intersection(p[i], q[i], p[j], q[j]);
                    g[i][j] = g[j][i] = onseg(p[i], q[i], tmp) && onseg(p[j], q[j], tmp);
                }
            }
        }
        //最后用Floyd-Warshall算法求任意两个棍子是否通过其他的棍子相连
        for (int k = 0; k < n; k++) {
            for (int i = 0; i < n; i++) {
                for (int j = 0; j < n; j++)
                    g[i][j] |= g[i][k] && g[k][j];
            }
        }
        int a, b;
        while (scanf("%d%d", &a, &b) && a + b) {
            if (g[a - 1][b - 1])
                printf("CONNECTED
");
            else
                printf("NOT CONNECTED
");
        }
    }
}

2、极限情况

1、AOJ 2308

 

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

const double g = 9.8;   //重力加速度
const double EPS = 1e-10;
int N;
double V, X, Y;
double L[55], B[55], R[55], T[55];

// 计算以 vy 的速度竖直向上射出 t 秒后的高度
double calc(double vy,double t)
{
    return vy*t-g*t*t/2;
}
// a 相对 lb 和 ub 的位置
int cmp(double lb, double ub, double a)
{ 
   return a < lb + EPS ? -1 : a > ub - EPS ? 1 : 0;
}
// 判断当射出路径经过 （qx, qy) 时，蛋能否击中猪
bool check(double qx, double qy)
{
    // 设初速度在 x 和 y 方向上的分量分别为 vx 和 vy，过（qx，qy）的时间为 t
    // 求解联立方程 vx^2 + vy^2 = v^2, vx * t = qx, vy * t - 1/2 * g * t^2 = qy
    double a = g * g / 4, b = g * qy - V * V, c = qx * qx + qy * qy;
    double D = b * b - 4 * a * c;
    if(D < 0 && D > -EPS) D = 0;
    if(D < 0) return false;
    for(int d = -1; d <= 1; d += 2)  // 验证联立方程的两个解
    {
        double t2 = (-b + d * sqrt(D)) / (2 * a);
        if(t2 <= 0) continue;
        double t = sqrt(t2);
        double vx = qx / t, vy = (qy + g * t * t / 2) / t;
        
        // 判断是否通过猪的正上方
        double yt = calc(vy, X / vx);
        if(yt < Y - EPS) continue;
        
        bool ok = true;
        for(int i = 0; i < N; i++)
        {
            if(L[i] >= X) continue;
            // 判断在猪正上方的鸟和猪之间是否有障碍物
            if(R[i] == X && Y <= T[i] && B[i] <= yt)
                ok = false;
            //判断在飞到猪正上方之前是否会撞到障碍物
            int yL = cmp(B[i], T[i], calc(vy, L[i] / vx));  //左侧的相对位置
            int yR = cmp(B[i], T[i], calc(vy, R[i] / vx));  //右侧的相对位置
            int xH = cmp(L[i], R[i], vx * (vy / g));   //最高点的相对位置
            int yH = cmp(B[i], T[i], calc(vy, vy / g));
            if(xH == 0 && yH >= 0 && yL < 0)
                ok = false;
            if(yL * yR  <= 0)
                ok = false;
        }
        if(ok)
            return true;
    }
    return false;
}
void solve(void)
{
    // 截掉猪以右的障碍物
    for(int i = 0; i < N; i++)
        R[i] = min(R[i], X);
    bool ok = check(X, Y);   // 直接撞上猪
    for(int i = 0; i < N; i++)
    {
        ok |= check(L[i], T[i]);
        ok |= check(R[i], T[i]);
    }
    puts(ok ? "Yes" : "No");
}
int main()
{
    scanf("%d %lf %lf %lf", &N, &V, &X, &Y);
    for(int i = 0; i < N; i++)
        scanf("%lf %lf %lf %lf", &L[i], &B[i], &R[i], &T[i]);
    solve();
    return 0;
}

2、POJ 1981 Circle and Points

 

题意：平面上有N个点，用单位圆去套，最多能套几个

思路：当然是先想到枚举两个，再跑一边所有点的O(n^3)的优秀算法，// 滑稽一下，不过这数据看起来也不一定会T掉，不管辣

point getmid(point p1, point p2)
{
    point mid,center;
    mid.x = (p1.x + p2.x) / 2.0;
    mid.y = (p1.y + p2.y) / 2.0;
    double angle = atan2(p1.x - p2.x, p2.y - p1.y); // 从平行四边形的角度去看，求经过中点的那条边的角度
    double dcm = sqrt(1-dis(p1, mid) * dis(p1, mid)); // 因为是单位圆，我们要计算出以二者中心为圆心，距离一半为半径的圆的半径和单位圆半径之差
    center.x = mid.x + dcm * cos(angle);//三角形思路。。。
    center.y = mid.y + dcm * sin(angle);
    return center;
}

计算两点之间圆心的函数

// #include<bits/stdc++.h>
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;

const int INF    = 0x3f3f3f3f;
const int MAX_N  = 500 + 50;
const int MOD    = 1000000007;
const double PI  = acos(-1.0);
const double EPS = 1e-10;

struct point
{
    double x,y;
    point(double x=0,double y=0):x(x),y(y) {}
}p[310];

struct node
{
    double ang;
    int in;
}arc[310*300];

double dis(point a,point b){return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));}

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

bool cmp(node a,node b)
{
    if(dcmp(a.ang-b.ang)==0) return a.in>b.in;
    return dcmp(a.ang-b.ang)<0;
}

int main(void)
{
    int i,j,n;
    while(scanf("%d",&n)&&n)
    {
        for(i = 1 ; i <= n; i++)
            scanf("%lf%lf",&p[i].x,&p[i].y);
        int g = 0;
        int ans = 0,maxz = 1;
        for(i = 1; i <= n ; i++)
        {
            ans = 0; g = 0;
            for(j = 1; j <= n ; j++)
            {
                if(dis(p[i],p[j])>2.0) continue;
                double ang1 = atan2(p[j].y-p[i].y,p[j].x-p[i].x);
                double ang2 = acos(dis(p[i],p[j])/2);
                arc[++g].ang = ang1-ang2;//这里角度的算法很巧妙
                arc[g].in = 1;
                arc[++g].ang = ang1+ang2;
                arc[g].in = -1;
            }
            sort(arc+1,arc+g+1,cmp);
            
            for(j = 1 ; j <= g;j++)
            {
                ans+=arc[j].in;
                maxz = max(maxz,ans);
            }
        }
        printf("%d
",maxz);
    }
    return 0;
}

3、POJ 1418

题意：给定Confetti的尺寸和位置以及它们的叠放次序，计算出有多少Confetti是可以看见的

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

const int  MAX_N = 128;
const double EPS = 5e-13;
const double PI = acos(-1.0);

typedef struct
{
    double x, y;
} point;

double Distance(const point & p1, const point & p2){return sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));}

double MainAngle(double a)
{
    while (a > 2 * PI) a -= 2 * PI;
    while (a < 0) a += 2 * PI;
    return a;
}

int  n;
point o[MAX_N];  //圆心
double r[MAX_N];  //圆的弧度

int  pan;   //与这个圆的交点数目
double pa[2 * MAX_N]; //存放与这个圆所有交点对应的弧度
int  visible[MAX_N];
int  ans;

int main()
{
    int i, j, k, t;
    point tp;
    double a, b, d;
    while (scanf("%d", &n)&&n)
    {
        for (i = 0; i < n; ++i)
        {
            scanf("%lf %lf %lf", &o[i].x, &o[i].y, &r[i]);
            visible[i] = 0;
        }

        for (i = 0; i < n; ++i)
        {
            pan = 0;
            pa[pan++] = 0;
            pa[pan++] = 2 * PI;

            for (j = 0; j < n; ++j)
            {
                if (j == i) continue;

                d = Distance(o[i], o[j]); //判断两个圆心距离

                if (r[i] + r[j] < d || r[i] + d < r[j] || r[j] + d < r[i]) //包含或不相交的
                    continue;

                a = atan2(o[j].y - o[i].y, o[j].x - o[i].x);//atan2(),是求这个点和x轴正方形夹角，*PI/180  得到度数
                b = acos((r[i] * r[i] + d * d - r[j] * r[j]) / (2 * r[i] * d));

                pa[pan] = MainAngle(a + b);
                pan++;
                pa[pan] = MainAngle(a - b);
                pan++;
            }

            sort(pa, pa + pan);

            for (j = 0; j < pan - 1; ++j)
            {
                a = (pa[j] + pa[j + 1]) / 2;

                for (t = -1; t <= 1; t += 2) //t = -1 或 1
                {
                    //将每段圆弧中点往内外各移动很小距离
                    tp.x = o[i].x + (r[i] + t * EPS) * cos(a);
                    tp.y = o[i].y + (r[i] + t * EPS) * sin(a);

                    for (k = n - 1; k >= 0; --k)
                        //如果找到第一个cover point i 的 Arc j 的圆，break
                        if (Distance(tp, o[k]) < r[k])
                            break;
                    visible[k] = 1;
                }
            }
        }
        ans = 0;
        for (i = 0; i < n; ++i)
            if (visible[i] == 1)
                ans++;

        printf("%d
", ans);
    }

    return 0;
}

4、AOJ 2201 Immortal Jewels

 题意：有n个宝石，宝石为圆形。给出现在有根金属棒，靠近宝石距离m之内就会吸附上去。现在给出每颗宝石的信息，求用这个棒子一次能钓出最多多少个宝石？

看起来简单但是感觉着实挺复杂，先空着吧，之后再来