hdu 4667（几何）

最新推荐文章于 2017-11-01 20:14:23 发布

laziercs

最新推荐文章于 2017-11-01 20:14:23 发布

阅读量1.1k

点赞数

CC 4.0 BY-SA版权

分类专栏：计算几何

本文链接：https://blog.youkuaiyun.com/zkfzkfzkfzkfzkfzkfzk/article/details/9981729

计算几何专栏收录该内容

4 篇文章

订阅专栏

本文深入探讨了如何运用解析几何和算法解决计算机科学领域的复杂问题，包括但不限于圆与圆之间的关系、凸包生成、三角形与圆的交互作用等。通过详细解析这些数学概念的应用实例，揭示了其在实际编程和算法设计中的价值。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

解题思路：求出所有可能在外围的点。也就是将两圆两两求外公切线，得到所有的切点；将圆和三角形的三个点依次做切线，

也得到切点，再加上所有三角形的三个点。最后对这些点求凸包，对于凸包上的每条边，如果它们在同一个圆上，用相应的圆弧替代。

#include<stdio.h>
#include<string.h>
#include<stdlib.h>
#include<time.h>
#include<math.h>
#include<bitset>
#include<string>
#include<queue>
#include<deque>
#include<ctype.h>
#include<vector>
#include<set>
#include<map>
#include<list>
#include<stack>
#include<functional>
#include<algorithm>

using namespace std;

#define N 100010

const double eps=1e-12;
const double pi=acos(-1.0);

struct Point{
    double x,y;
    int id;
    Point(double a=0,double b=0):x(a),y(b){}
};
typedef Point Vector;

Point operator-(Point a,Point b){ return Point(a.x-b.x,a.y-b.y); }
Point operator+(Point a,Point b){ return Point(a.x+b.x,a.y+b.y); }
double operator*(Vector a,Vector b){ return a.x*b.x+a.y*b.y; }
int operator==(Point a,Point b){ return a.x==b.x&&a.y==b.y; }
double Len(Vector a){  return sqrt(a.x*a.x+a.y*a.y); }

Vector Rotate(Vector a,double rad){ //向量逆时针旋转rad弧度
    return Vector(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));
}

int dcmp(double x){
    if(x>-eps&&x<eps) return 0;
    if(x>eps) return 1;
    return -1;
}
double cross(Point a,Point b){ //叉积
    return a.x*b.y-a.y*b.x;
}
double angle(Vector a,Vector b){ //两个向量的夹角
    return acos((a*b)/(Len(a)*Len(b)));
}
int online(Point a,Point b,Point c){ //判断点c是否在线段ab上
    return  (a==c)||(b==c)||(cross(a-c,b-c)==0&&((a-c)*(b-c)<0));
}

int segmentIntersection(Point a1,Point a2,Point b1,Point b2) //判断线段a1a2和b1b2相交（不包含端点）
{
    double c1=cross((a2-a1),(b1-a1)),c2=cross((a2-a1),(b2-a1));
    double c3=cross((b2-b1),(a1-b1)),c4=cross((b2-b1),(a2-b1));
    return dcmp(c1)*dcmp(c2)<0&&dcmp(c3)*dcmp(c4)<0;
}

int Polygon_contains(Point point,Point *points,int n) {  //判断一个点在多边形内
    int i, j, status=0;
    for (i=0,j=n-1;i<n;j=i++){
        if(online(points[i],points[j],point)) return 0;
        if((((points[i].y<=point.y)&&(point.y<points[j].y))
        ||((points[j].y<=point.y)&&(point.y<points[i].y)))
        &&(point.x <(points[j].x - points[i].x)*(point.y - points[i].y)
        /(points[j].y - points[i].y) + points[i].x))
            status = !status;
    }
    return status;
}

Point P[N];
bool cmp(Point a,Point b){   //极坐标排序
    if(cross(a-P[0],b-P[0])==0)
        return dcmp(Len(a-P[0])-Len(b-P[0]))<=0;
    else
        return cross(a-P[0],b-P[0])>0;
}

void Graham(int n,vector<Point> &sta){  //求解二维凸包
    for(int i=1;i<n;i++){
        if(dcmp(P[i].y-P[0].y)<0||(P[i].y==P[0].y&&P[i].x<P[0].x))
            swap(P[i],P[0]);
    }
    sort(P+1,P+n,cmp);
    sta.push_back(P[0]);
    sta.push_back(P[1]);
    for(int i=2;i<n;i++){
        while(sta.size()>1&&dcmp(cross(sta[sta.size()-1]-sta[sta.size()-2],P[i]-sta[sta.size()-1]))<=0)
            sta.pop_back();
        sta.push_back(P[i]);
    }
}
struct circle{
    double x,y,r;
    Point make_Point(double a){
        return Point(x+cos(a)*r,y+sin(a)*r);
    }
}C[60];
int getTangents(Point p,circle cc,Vector *t){ //求一个点p到圆cc的切线向量
    Vector u=Point(cc.x,cc.y)-p;
    double dist=Len(u);
    if(dcmp(dist-cc.r)<0) return 0;
    else if(dcmp(dist-cc.r)==0){
        t[0]=Rotate(u,pi/2);
        return 1;
    }else {
        double ang=asin(cc.r/dist);
        t[0]=Rotate(u,-ang);
        t[1]=Rotate(u,+ang);
        return 2;
    }
}
void meet(Point p,circle cc,Point *t){ //求一个点p到圆cc的切点
    double dis,l;
    Vector u=Point(cc.x,cc.y)-p;
    getTangents(p,cc,t);l=Len(t[0]);
    dis=sqrt(Len(u)*Len(u)-cc.r*cc.r);
    t[0].x=t[0].x*dis/l;t[0].y=t[0].y*dis/l;
    t[1].x=t[1].x*dis/l;t[1].y=t[1].y*dis/l;
    t[0]=t[0]+p;
    t[1]=t[1]+p;
}

double D(Point a,Point b,int id){ //求a点到b点在的圆上的圆弧长度
    double ang1,ang2;
    Vector v1,v2;
    v1=a-Point(C[id].x,C[id].y);
    v2=b-Point(C[id].x,C[id].y);
    ang1=atan2(v1.y,v1.x);
    ang2=atan2(v2.y,v2.x);
    if(ang2<ang1) ang2+=2*pi;
    return C[id].r*(ang2-ang1);
}

int getgqx(circle A,circle B,Point* sa,Point* sb) //求圆A和B的公切线
{
    int cnt=0;
    Point *a=sa,*b=sb;
    if(A.r<B.r){ swap(A,B);swap(a,b);}
    double d2=(A.x-B.x)*(A.x-B.x)+(A.y-B.y)*(A.y-B.y);
    double rdiff=A.r-B.r;
    double rsum=A.r+B.r;
    if(dcmp(d2-rdiff*rdiff)<0) return 0;    //内含

    double base=atan2(B.y-A.y,B.x-A.x);
    if(dcmp(d2)==0&&dcmp(A.r-B.r)==0) return -1; //无限多条切线
    if(dcmp(d2-rdiff*rdiff)==0)
    {
        a[cnt]=A.make_Point(base);b[cnt]=B.make_Point(base);cnt++;
        return 1;
    }
    //有外共切线
    double ang=acos((A.r-B.r)/sqrt(d2));
    a[cnt]=A.make_Point(base+ang);b[cnt]=B.make_Point(base+ang);cnt++;
    a[cnt]=A.make_Point(base-ang);b[cnt]=B.make_Point(base-ang);cnt++;
    if(dcmp(d2-rsum*rsum)==0)//一条内共切线
    {
        a[cnt]=A.make_Point(base);b[cnt]=B.make_Point(pi+base);cnt++;
    }
    else if(dcmp(d2-rsum*rsum)>0)
    {
        double ang=acos((A.r+B.r)/sqrt(d2));
        a[cnt]=A.make_Point(base+ang);b[cnt]=B.make_Point(pi+base+ang);cnt++;
        a[cnt]=A.make_Point(base-ang);b[cnt]=B.make_Point(pi+base-ang);cnt++;
    }
    return cnt;
}

int main(){
    int n,m,i,M,j;
    while(scanf("%d%d",&n,&m)!=EOF){
        for(i=0;i<n;i++){
            scanf("%lf%lf%lf",&C[i].x,&C[i].y,&C[i].r);
        }
        for(i=0;i<m;i++){
            scanf("%lf%lf",&P[3*i].x,&P[3*i].y);
            scanf("%lf%lf",&P[3*i+1].x,&P[3*i+1].y);
            scanf("%lf%lf",&P[3*i+2].x,&P[3*i+2].y);
        }
        if(n==1&&m==0){
            printf("%lf\n",pi*2*C[0].r);
            continue;
        }
        M=3*m;m*=3;
        for(i=0;i<M;i++) P[i].id=i-M-M;
        double ans=0;
        vector<Point> sta;
        if(M>0){
            Graham(M,sta);
            for(i=0;i<sta.size();i++) P[i]=sta[i];
            M=sta.size();m=M;
        }
        Point t[4],T[4];
        for(i=0;i<n;i++){
            for(j=i+1;j<n;j++){
                getgqx(C[i],C[j],t,T);
                t[0].id=t[1].id=i;
                T[0].id=T[1].id=j;
                P[m++]=t[0];P[m++]=t[1];
                P[m++]=T[0];P[m++]=T[1];
            }
        }
        for(i=0;i<M;i++){
            for(j=0;j<n;j++){
                meet(P[i],C[j],t);
                t[0].id=t[1].id=j;
                P[m++]=t[0];P[m++]=t[1];
            }
        }
        M=m;
        sta.clear();
        Graham(M,sta);
        if(sta[0].id==sta[sta.size()-1].id)
            ans+=D(sta[sta.size()-1],sta[0],sta[0].id);
        else
            ans+=Len(sta[0]-sta[sta.size()-1]);
        for(i=1;i<int(sta.size());i++){
            if(sta[i].id==sta[i-1].id)
                ans+=D(sta[i-1],sta[i],sta[i].id);
            else
                ans+=Len(sta[i]-sta[i-1]);
        }
        printf("%.8lf\n",ans);
    }
    return 0;
}