计算几何——三角剖分 uva1331+poj3675

Posted 2020-12-13 zsben991126

看了一整天三角剖分。。能找到的题只有求最大面积最小的三角剖分。。

然而这题除了用些三角剖分的性质外。。其实是个dp

#include <bits/stdc++.h>
using namespace std;

const int INF = 1e9;
const int maxn = 50 + 5;
const double eps = 0.00001;
int n, m;
struct Point {
    int x, y;
    friend Point operator-(const Point &lhs, const Point &rhs) {
        return {rhs.x - lhs.x, rhs.y - lhs.y};
    }
    // 公式: a(x1, y1) * b(x2, y2) = x1 * y2 - x2 * y1;
    friend double operator*(const Point &lhs, const Point &rhs) {
        return std::abs(lhs.x * rhs.y - lhs.y * rhs.x) * 1.0;
    }
}A[maxn];

// 判断顶点i, 顶点j, 顶点k, 组成的三角形是否合法(即内部没有其他顶点)
// 合法就返回面积, 非法返回0
// 判断是否合法:
//  判断是否有一个点在三角形内部, 如果一个点在内部, 三角形的三个端点连接内部这个点, 将面积分为了3部分,
//  判断这3部分相加是否等于原三角形面积. 相等说明有点在内部.

// 叉积算三角形的面积
double area(int i, int j, int k) {
    Point a = A[k] - A[i];
    Point b = A[k] - A[j];
    // 四边形面积除以2便是三角形的面积
    return a * b * 1.0 / 2;
}

double is_valid(int i, int j, int k) {
    for (int p = 1; p <= m; ++p) {
        if (!(p == i || p == j || p == k)) {
            double area_S = area(i, j, k);
            double area_1 = area(i, j, p);
            double area_2 = area(i, k, p);
            double area_3 = area(k, j, p);
            // 非法
            if (std::fabs(area_S - (area_1 + area_2 + area_3)) < eps)
                return 0;
        }
    }
    // 合法
    return area(i, j, k);
}

// dp[i][j]: 表示从i, i + 1, i + 2, ... , j,(看成一个圈, j的下一个元素为i, 即i和j相邻), 这些顶点中所剖分的三角形中, 最大三角形面积的最小值
double dp[maxn][maxn];
double dfs(int L, int R) {
    // 记忆化
    if (dp[L][R] != -1)
        return dp[L][R];
    // 递归的边界
    if (std::abs(R- L) <= 2)
        return dp[L][R] = area(L, L + 1, R);
    double ret = INF;
    for (int k = L + 1; k < R; ++k) {
        if (is_valid(L, R, k)) {
            ret = std::min(ret, std::max(dfs(L, k), std::max(area(L, R, k), dfs(k, R))));
        }
    }
    return dp[L][R] = ret;
}

int main() {
    scanf("%d", &n);
    while (n--) {
        // 初始化
        for (int i = 0; i < maxn; ++i)
            std::fill(dp[i], dp[i] + maxn, -1);
        scanf("%d", &m);
        for (int i = 1; i <= m; ++i) {
            scanf("%d %d", &A[i].x, &A[i].y);
        }
        printf("%.1lf
", dfs(1, m));
    }
    return 0;
}

 

poj3675这个用三角剖分能做，但是我感觉直接圆+三角形面积交应该也行（poj真的毒。。）

#include<cstdio>
#include<vector>
#include<cmath>
#include<math.h>
#include<string>
#include<string.h>
#include<iostream>
#include<algorithm>
#include<map>
#define PI acos(-1.0)
#define pb push_back
#define F first
#define S second
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
const int N=1005;
const int MOD=1e9+7;
const double EPS=1e-10;
int sign(double x) {  //三态函数，减少精度问题
    return abs(x) < EPS ? 0 : x < 0 ? -1 : 1;
}
struct Point { //点的定义
    double x, y;
    Point(double x=0.0, double y=0.0) : x(x), y(y) {}
    Point operator + (const Point &rhs) const {  //向量加法
        return Point(x + rhs.x, y + rhs.y);
    }
    Point operator - (const Point &rhs) const {  //向量减法
        return Point(x - rhs.x, y - rhs.y);
    }
    Point operator * (double p) const {  //向量乘以标量
        return Point(x * p, y * p);
    }
    Point operator / (double p) const {  //向量除以标量
        return Point(x / p, y / p);
    }
    bool operator < (const Point &rhs) const {  //点的坐标排序
        return x < rhs.x || (x == rhs.x && y < rhs.y);
    }
    bool operator == (const Point &rhs) const {  //判断同一个点
        return sign(x - rhs.x) == 0 && sign(y - rhs.y) == 0;
    }
};
struct Circle {  //圆的定义
    Point c;  //圆心
    double r;  //半径
    Circle() {}
    Circle(Point c, double r) : c(c), r(r) {}
    Point point(double a) {  //圆上的一点
        return Point(c.x+cos(a)*r, c.y+sin(a)*r);
    }
};
typedef Point Vector;
Point p[N];
 
double dot(Vector A,Vector B){
    return A.x*B.x+A.y*B.y;
}
double length(Vector A){
    return sqrt(dot(A,A));
}
double cross(Vector A,Vector B){
    return A.x*B.y-A.y*B.x;
}
double Angle(Vector A,Vector B){
    double t=dot(A,B)/length(A)/length(B);
    return acos(t);
}
bool on_seg(Point p, Point a, Point b) {  //判断点在线段上（不包含端点），两点式
    return sign(cross(a-p, b-p)) == 0 && sign(dot(a-p, b-p)) <= 0;  //点p是否在线段ab上
}
 
double CulArea( Point A, Point B,Circle C)
{
    Vector OA = A-C.c, OB = B-C.c;
    Vector BA = A-B, BC = C.c-B;
    Vector AB = B-A, AC = C.c-A;
    double DOA = length(OA), DOB = length(OB),DAB = length(AB), r = C.r;
    if(sign(cross(OA,OB)) == 0) return 0;
    if(sign(DOA-C.r) < 0 && sign(DOB-C.r) < 0) return cross(OA,OB)*0.5;
    else if(DOB < r && DOA >= r) {
        double x = (dot(BA,BC) + sqrt(r*r*DAB*DAB-cross(BA,BC)*cross(BA,BC)))/DAB;
        double TS = cross(OA,OB)*0.5;
        return asin(TS*(1-x/DAB)*2/r/DOA)*r*r*0.5+TS*x/DAB;
    }
    else if(DOB >= r && DOA < r) {
        double y = (dot(AB,AC)+sqrt(r*r*DAB*DAB-cross(AB,AC)*cross(AB,AC)))/DAB;
        double TS = cross(OA,OB)*0.5;
        return asin(TS*(1-y/DAB)*2/r/DOB)*r*r*0.5+TS*y/DAB;
    }
    else if(fabs(cross(OA,OB)) >= r*DAB || dot(AB,AC) <= 0 || dot(BA,BC) <= 0) {
        if(dot(OA,OB) < 0){
            if(cross(OA,OB) < 0) return (-acos(-1.0)-asin(cross(OA,OB)/DOA/DOB))*r*r*0.5;
            else return ( acos(-1.0)-asin(cross(OA,OB)/DOA/DOB))*r*r*0.5;
        }
        else return asin(cross(OA,OB)/DOA/DOB)*r*r*0.5;
    }
    else {
        double x = (dot(BA,BC)+sqrt(r*r*DAB*DAB-cross(BA,BC)*cross(BA,BC)))/DAB;
        double y = (dot(AB,AC)+sqrt(r*r*DAB*DAB-cross(AB,AC)*cross(AB,AC)))/DAB;
        double TS = cross(OA,OB)*0.5;
        return(asin(TS*(1-x/DAB)*2/r/DOA)+asin(TS*(1-y/DAB)*2/r/DOB))*r*r*0.5 + TS*((x+y)/DAB-1);
    }
}
 
int main(){
    double R;
    while(scanf("%lf",&R)==1){
        Circle C1=Circle({0,0},R);//圆
 
        int n;
        scanf("%d",&n);
        double ans=0.0;
        for(int i=0;i<n;i++){
            double x,y;
            scanf("%lf%lf",&x,&y);
            p[i]={x,y};
        }
        for(int i=0;i<n;i++)    ans+=CulArea(p[i],p[(i+1)%n],C1);
        if(ans<0)   ans=-ans;
        printf("%.2f
",ans);
    }
 
    return 0;