三维几何，四面体（压纸器，LA 4795）

来源：互联网发布：德力西知乎编辑：程序博客网时间：2024/04/26 19:15
学会了求四面体重心的公式和算法。
注意讨论清楚四点共面的情况。
最后是求芯片到底的距离而不是重心到底的距离。。。
理解题意是你的硬伤哎。。其实就是没有用心看，用心想，总是胡乱看一遍然后自己脑补。
代码
#include<bits/stdc++.h>#include<limits.h>using namespace std;struct Point3{    double x,y,z;    Point3(double x=0,double y=0,double z=0):x(x),y(y),z(z){}    int Read()    {        scanf("%lf",&x);        if(scanf("%lf",&y)==-1) return -1;        else return scanf("%lf",&z);    }};typedef Point3 Vector3;const double eps = 1e-10;int dcmp(double x){    if(fabs(x)<eps) return 0;    else return x<0?-1:1;}Point3 operator + (Point3 A,Vector3 B){    return Point3(A.x+B.x,A.y+B.y,A.z+B.z);}Vector3 operator - (Point3 A,Point3 B){    return Vector3(A.x-B.x,A.y-B.y,A.z-B.z);}Vector3 operator * (Vector3 A,double t){    return Vector3(A.x*t,A.y*t,A.z*t);}Vector3 operator / (Vector3 A,double t){    return Vector3(A.x/t,A.y/t,A.z/t);}double Dot(Vector3 A,Vector3 B){    return A.x*B.x+A.y*B.y+A.z*B.z;}Vector3 Cross(Vector3 A,Vector3 B){    return Vector3(A.y*B.z-A.z*B.y,A.z*B.x-A.x*B.z,A.x*B.y-A.y*B.x);}double Len(Vector3 V){    return sqrt(Dot(V,V));}double DistToPlane(Point3 p,Point3 p0,Vector3 n){    return Dot(p-p0,n);}double Vol6(Point3 A,Point3 B,Point3 C,Point3 D){    return fabs(Dot(D-A,Cross(A-B,A-C)));}double Area2(Point3 A,Point3 B,Point3 C){    return Len(Cross(B-A,C-A));}bool PointInTri(Point3 P,Point3 A,Point3 B,Point3 C){    return dcmp(Area2(P,A,B)+Area2(P,B,C)+Area2(P,C,A)-Area2(A,B,C))==0;}Point3 GetSMTZX(Point3 A,Point3 B,Point3 C,Point3 D){    return (A+B+C+D)/4;}bool SegSegNormalIntersection(Point3 A,Point3 B,Point3 C,Point3 D){    Vector3 V1=C-A;    Vector3 V2=D-A;    Vector3 V3=A-C;    Vector3 V4=B-C;    Vector3 C1=Cross(B-A,V1);    Vector3 C2=Cross(B-A,V2);    Vector3 C3=Cross(D-C,V3);    Vector3 C4=Cross(D-C,V4);    return Dot(C1,C2)<0&&Len(Cross(C1,C2))==0&&Dot(C3,C4)<0&&Len(Cross(C3,C4))==0;}Point3 P[10];Point3 F;Point3 ZX;bool wen(Point3 A,Point3 B,Point3 C,Vector3 N){    Vector3 n;    n=Cross(N,B-A);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,A,n)-0.2)<0) return false;    n=Cross(N,C-B);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,B,n)-0.2)<0) return false;    n=Cross(N,A-C);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,C,n)-0.2)<0) return false;    return true;}bool WEN(Point3 A,Point3 B,Point3 C,Point3 D,Vector3 N){    Vector3 n;    n=Cross(N,B-A);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,A,n)-0.2)<0) return false;    n=Cross(N,C-B);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,B,n)-0.2)<0) return false;    n=Cross(N,D-C);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,C,n)-0.2)<0) return false;    n=Cross(N,A-D);    n=n/Len(n);    if(dcmp(DistToPlane(ZX,D,n)-0.2)<0) return false;    return true;}void GetZX(){    double V1=Vol6(P[0],P[1],P[2],P[3]);    double V2=Vol6(P[0],P[1],P[2],P[4]);    Point3 A=GetSMTZX(P[0],P[1],P[2],P[3]);    Point3 B=GetSMTZX(P[0],P[1],P[2],P[4]);    ZX=(A*V1+B*V2)/(V1+V2);}int kase;int main(){    while(~P[0].Read())    {        for(int i=1;i<5;i++)            P[i].Read();        F.Read();        GetZX();        double MAX=DBL_MIN;        double MIN=DBL_MAX;        for(int i=0;i<5;i++)            for(int j=0;j<i;j++)                for(int k=0;k<i;k++) if(j!=k)                {                    Vector3 n=Cross(P[j]-P[i],P[k]-P[j]);                    n=n/Len(n);                    bool ok=true;                    int si=-1;                    for(int l=0;l<5;l++)                    {                        if(l==i) continue;                        if(l==j) continue;                        if(l==k) continue;                        double DIST=DistToPlane(P[l],P[i],n);                        if(dcmp(DIST)<0)                        {                            ok=false;                            break;                        }                        if(dcmp(DIST)==0) si=l;                    }                    if(!ok) continue;                    if(si!=-1)                    {                        if(PointInTri(P[i],P[j],P[k],P[si])) continue;                        if(PointInTri(P[j],P[k],P[i],P[si])) continue;                        if(PointInTri(P[k],P[i],P[j],P[si])) continue;                        if(SegSegNormalIntersection(P[i],P[j],P[k],P[si])&&!WEN(P[i],P[si],P[j],P[k],n)) continue;                        if(SegSegNormalIntersection(P[j],P[k],P[i],P[si])&&!WEN(P[j],P[si],P[k],P[i],n)) continue;                        if(SegSegNormalIntersection(P[k],P[i],P[j],P[si])&&!WEN(P[k],P[si],P[i],P[j],n)) continue;                    }                    else if(!wen(P[i],P[j],P[k],n)) continue;                    double DIST=DistToPlane(F,P[i],n);                    MIN=min(MIN,DIST);                    MAX=max(MAX,DIST);                }        printf("Case %d: %.5lf %.5lf\n",++kase,MIN,MAX);    }    return 0;}
0 0