计算几何 模板

来源:互联网 发布:linux打开u盘命令 编辑:程序博客网 时间:2024/06/15 07:04
//凸包+最远点对#include<iostream>#include<cstdio>#include<cstring>#include<cstdlib>#include<cmath>#include<algorithm>#include<vector>using namespace std;const double eps = 1e-8;int cmp(double x){    if(fabs(x)<eps)        return 0;    if(x>0)        return 1;    return -1;}inline double sqr(double x){    return x*x;}struct point{    double x,y;    point(){}    point(int xx,int yy)    {        x=xx;        y=yy;    }    friend point operator + (const point &a,const point &b)    {        return point(a.x+b.x,a.y+b.y);    }    friend point operator - (const point &a,const point &b)    {        return point(a.x-b.x,a.y-b.y);    }    friend bool operator == (const point &a,const point &b)    {        return cmp(a.x-b.x)==0 && cmp(a.y-b.y)==0;    }    void input()    {        scanf("%lf%lf",&x,&y);    }    double norm()    {        return sqrt(sqr(x)+sqr(y));    }};struct polygon_convex{    vector<point>P;    polygon_convex(int Size=0)    {        P.resize(Size);    }};bool comp_less(const point &a,const point &b){    return cmp(a.x-b.x)<0 || (cmp(a.x-b.x==0) && cmp(a.y-b.y)<0);}double det(point a,point b){    return a.x*b.y-a.y*b.x;}double dist(const point &a,const point &b){    return (a-b).norm();}long long dd(const point &a,const point &b){    return sqr(a.x-b.x)+sqr(a.y-b.y);}polygon_convex convex_hull(vector<point> a){    polygon_convex res(2*a.size()+5);    sort(a.begin(),a.end(),comp_less);    a.erase(unique(a.begin(),a.end()),a.end());    int m=0;    for(int i=0;i<a.size();++i)    {        while(m>1              &&              cmp(det( res.P[m-1]-res.P[m-2],a[i]-res.P[m-2] ))<=0)            --m;        res.P[m++]=a[i];    }    int k=m;    for(int i=int(a.size())-2;i>=0;--i)    {        while(m>k && cmp(det(res.P[m-1]-res.P[m-2],a[i]-res.P[m-2]))<=0)--m;        res.P[m++]=a[i];    }    res.P.resize(m);    if(a.size()>1)        res.P.resize(m-1);    return res;}long long convex_diameter(polygon_convex &a,int &First,int &Second){    vector<point> &p = a.P;    int n=p.size();    long long maxd=0;    if(n==1)    {        First=Second=0;        return maxd;    }    #define next(i) ((i+1)%n)    for(int i=0,j=1;i<n;i++)    {        while(cmp(det(p[next(i)]-p[i],p[j]-p[i])-det(p[next(i)]-p[i],p[next(j)]-p[i]))<0)        {            j=next(j);        }        long long d=dd(p[i],p[j]);        if(d>maxd)        {            maxd=d;            First=i,Second=j;        }        d=dd(p[next(i)],p[next(j)]);        if(d>maxd)        {            maxd=d;            First=i,Second=j;        }    }    return maxd;}int n;vector<point> p;polygon_convex pc;point pp;int main(){    scanf("%d",&n);    for(int i=1;i<=n;i++)    {        pp.input();        p.push_back(pp);    }    int x,y;    pc = convex_hull(p);    long long dis = convex_diameter(pc,x,y);    cout<<dis;}


//圆与多边形面积交#include<iostream>#include<cmath>#include<cstdlib>#include<cstdio>#include<cmath>using namespace std;const int Max = 1011;const double PI =acos(-1.0);const double pi =acos(-1.0);double ax,ay,bx,by,K;double cx,cy,tx,ty,R;int n;const double eps=1e-8;int cmp(double x){    if(fabs(x)<eps)        return 0;    if(x>0)        return 1;    return -1;}inline double sqr(double x){    return x*x;}struct point{    double x,y;    point(){};    point(double a,double b):x(a),y(b){}    friend point operator +(const point &a,const point &b)    {        return point(a.x+b.x,a.y+b.y);    }    friend point operator -(const point &a,const point &b)    {        return point(a.x-b.x,a.y-b.y);    }    friend bool operator ==(const point &a,const point &b)    {        return (cmp(a.x-b.x)==0)&&(cmp(a.y-b.y)==0);    }    friend point operator *(const point &a,const double &b)    {        return point(a.x*b,a.y*b);    }    friend point operator *(const double &a,const point &b)    {        return point(a*b.x,a*b.y);    }    friend point operator /(const point &a,const double &b)    {        return point(a.x/b,a.y/b);    }    double norm()    {        return sqrt(x*x+y*y);    }};double det(const point &a,const point &b){    return a.x*b.y-a.y*b.x;}double dot(const point &a,const point &b){    return a.x*b.x+a.y*b.y;}double dist(const point &a,const point &b){    return (a-b).norm();}point rotate_point(const point &a,const double A){    double tx=a.x,ty=a.y;    return point(tx*cos(A)-ty*sin(A),tx*sin(A)+ty*cos(A));}int dcmp(double k){    return k<-eps?-1:k>eps?1:0;}double cross(const point &a,const point &b){    return a.x*b.y-a.y*b.x;}double abs(const point &o){    return sqrt(dot(o,o));}point crosspt(const point &a,const point &b,const point &p,const point &q){    double a1=cross(b-a,p-a);    double a2=cross(b-a,q-a);    return (p*a2-q*a1)/(a2-a1);}point res[Max];double r;double mysqrt(double n){    return sqrt(max(0.0,n));}void circle_cross_line(const point &a,const point &b,const point &o,double r,point ret[],int &num){    double x0=o.x,y0=o.y;    double x1=a.x,y1=a.y;    double x2=b.x,y2=b.y;    double dx=x2-x1,dy=y2-y1;    double A=dx*dx+dy*dy;    double B=2*dx*(x1-x0)+2*dy*(y1-y0);    double C=(x1-x0)*(x1-x0)+(y1-y0)*(y1-y0)-r*r;    double delta=B*B-4*A*C;    num=0;    if(dcmp(delta)>=0)    {        double t1=(-B-mysqrt(delta))/(2*A);        double t2=(-B+mysqrt(delta))/(2*A);        if(dcmp(t1-1)<=0&&dcmp(t1)>=0) ret[num++]=point(x1+t1*dx,y1+t1*dy);        if(dcmp(t2-1)<=0&&dcmp(t2)>=0) ret[num++]=point(x1+t2*dx,y1+t2*dy);    }}double sector_area(const point &a,const point &b){    double theta=atan2(a.y,a.x)-atan2(b.y,b.x);    while(theta<=0) theta+=2*PI;    while(theta>2*PI) theta-=2*PI;    theta=min(theta,2*PI-theta);    return r*r*theta/2;}double calc(const point &a,const point &b){    point p[2];    int num =0;    int ina=dcmp(abs(a)-r)<0;    int inb=dcmp(abs(b)-r)<0;    if(ina)    {        if(inb)        {            return fabs(cross(a,b))/2.0;        }        else        {            circle_cross_line(a,b,point(0,0),r,p,num);            return sector_area(b,p[0])+fabs(cross(a,p[0]))/2.0;        }    }    else    {        if(inb)        {            circle_cross_line(a,b,point(0,0),r,p,num);            return sector_area(p[0],a)+fabs(cross(p[0],b))/2.0;        }        else        {            circle_cross_line(a,b,point(0,0),r,p,num);            //cout<<"num="<<num<<endl;            if(num==2)            {                //cout<<"num=2"<<endl;                return sector_area(a,p[0])+sector_area(p[1],b)+fabs(cross(p[0],p[1]))/2.0;            }            else            {                //cout<<" not in "<<a.x<<" "<<a.y<<" "<<b.x<<" "<<b.y<<endl;                return sector_area(a,b);            }        }    }}double area(){    double ret=0;    for(int i=0;i<n;i++){        int sgn=dcmp(cross(res[i],res[i+1]));        if(sgn!=0){            ret+=sgn*calc(res[i],res[i+1]);            //cout<<"sgn="<<sgn<<" calc="<<calc(res[i],res[i+1])<<endl;            //cout<<"point "<<res[i].x<<" "<<res[i].y<<" "<<res[i+1].x<<" "<<res[i+1].y<<endl;            //cout<<"ret="<<ret<<endl;        }    }    return ret;}int main(){    int it=0;    while(scanf("%d",&n)==1)    {        it++;        //cout<<"PI="<<PI<<endl;        scanf("%lf",&K);        //K=1/K;        for(int i=0;i<n;i++)            scanf("%lf%lf",&res[i].x,&res[i].y);        res[n].x=res[0].x;        res[n].y=res[0].y;        scanf("%lf%lf%lf%lf",&ax,&ay,&bx,&by);        tx=ax+(bx-ax)/(1+K);        ty=ay+(by-ay)/(1+K);        cx=ax+(bx-ax)+(bx-ax)*K/(1-K);        cy=ay+(by-ay)+(by-ay)*K/(1-K);        cx=(cx+tx)/2;        cy=(cy+ty)/2;        R=sqrt((tx-cx)*(tx-cx)+(ty-cy)*(ty-cy));        r=R;        for(int i=0;i<=n;i++)        {            res[i].x-=cx;            res[i].y-=cy;            //cout<<"res["<<i<<"].x="<<res[i].x<<endl;            //cout<<"res["<<i<<"].y="<<res[i].y<<endl;        }        printf("Case %d: %.8f\n",it,fabs(area()));    }}


0 0