标签：平面 online inter return space -- clu src include
傻逼题。我从来没见过eps这样的。。。
打破了我对计算几何美好的幻想。
eps=1e-6=>wa3 eps=1e-8->wa2 eps 1e-4->AC
真的自闭，真的猜不到eps要设成1e-4，最后看了先杭电的代码，发现他们是1e-4,我就改了下，然后？？？？过了？？？？
nmsl
nmsl
nmsl
思路正确无比。赛后发现板子错了，把板子删了，tle3，把三分套三分换成了三分+板子，就开始无限WA制了
精简一下就是一句话：求线段到三角形的距离。
前置技能：三分，点到三角形的距离，点到线段的距离，点到平面的投影，点是否在三角形内
然后就没了。赛时-13赛后+37
1 #include <bits/stdc++.h> 2 using namespace std; 3 typedef double db; 4 const db eps = 1e-4; 5 const db pi = acos(-1); 6 int sign(db k){if(k>eps)return 1;else if(k<-eps)return -1;return 0;} 7 int cmp(db k1,db k2){return sign(k1-k2);} 8 int inmid(db k1,db k2,db k3){return sign(k1-k3)*sign(k2-k3)<=0;}// k3 在 [k1,k2] 内 9 struct point{ 10 db x,y,z; 11 point operator + (point k1){return (point){x+k1.x,y+k1.y,z+k1.z};} 12 point operator - (point k1){return (point){x-k1.x,y-k1.y,z-k1.z};} 13 point operator * (db k1){return (point){x*k1,y*k1,z*k1};} 14 point operator / (db k1){return (point){x/k1,y/k1,z/k1};} 15 db length()const {return sqrt(x*x+y*y+z*z);} 16 db len2()const {return x*x+y*y+z*z;} 17 }; 18 point det(const point &a,const point &b){ 19 return point{a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x}; 20 } 21 db dis(const point &a,const point &b){//距离 22 return sqrt(pow(a.x-b.x,2)+pow(a.y-b.y,2)+pow(a.z-b.z,2)); 23 } 24 struct line{ 25 point p[2]; 26 line (point k1,point k2){p[0]=k1; p[1]=k2;} 27 point& operator [] (int k){return p[k];} 28 }; 29 db dot(point k1,point k2){return k1.x*k2.x+k1.y*k2.y+k1.z*k2.z;} 30 int dot_online_in(point p,line l){//点在线段上，包括端点 31 return sign(det(p-l[0],p-l[1]).length())==0&&(l[0].x-p.x)*(l[1].x-p.x)<eps&& 32 (l[0].y-p.y)*(l[1].y-p.y)<eps&&(l[0].z-p.z)*(l[1].z-p.z)<eps; 33 } 34 struct plane{ 35 point p[3]; 36 plane(point k1,point k2,point k3){p[0]=k1,p[1]=k2,p[2]=k3;} 37 point& operator [] (int k){return p[k];} 38 }; 39 point pvec(point s1,point s2,point s3){ 40 return det(s1-s2,s2-s3); 41 } 42 point pvec(plane s){return pvec(s[0],s[1],s[2]);} 43 int dot_inplane_in(point p,plane s){ 44 return sign(det(s[0]-s[1],s[0]-s[2]).length()-det(p-s[0],p-s[1]).length()- 45 det(p-s[1],p-s[2]).length()-det(p-s[2],p-s[0]).length())==0; 46 } 47 int sameside(point p1,point p2,plane s){//两点在平面同侧 48 // point tmp = pvec(s); 49 return dot(pvec(s),p1-s[0])*dot(pvec(s),p2-s[0])>eps; 50 } 51 int intersect_in(line l,plane s){//线段和三角形是否有交点包含边界 52 return !sameside(l[0],l[1],s)&& 53 !sameside(s[2],s[0],plane{l[0],l[1],s[1]})&& 54 !sameside(s[0],s[1],plane{l[0],l[1],s[2]})&& 55 !sameside(s[1],s[2],plane{l[0],l[1],s[0]}); 56 } 57 point intersection(line l,plane s){//直线与平面的交点 58 point ret = pvec(s); 59 db t = ((ret.x*(s[0].x-l[0].x)+ret.y*(s[0].y-l[0].y))+ret.z*(s[0].z-l[0].z))/ 60 ((ret.x*(l[1].x-l[0].x)+ret.y*(l[1].y-l[0].y))+ret.z*(l[1].z-l[0].z)); 61 ret = l[0]+(l[1]-l[0])*t;return ret; 62 } 63 int t; 64 point p[22]; 65 int check(plane a,plane b){//check 0; 66 bool f=0; 67 if(intersect_in({a[0],a[1]},b))f=1; 68 if(intersect_in({a[1],a[2]},b))f=1; 69 if(intersect_in({a[0],a[2]},b))f=1; 70 if(intersect_in({b[0],b[1]},a))f=1; 71 if(intersect_in({b[0],b[2]},a))f=1; 72 if(intersect_in({b[1],b[2]},a))f=1; 73 return f; 74 } 75 db slove(plane a,plane b){//b三个点到a的距离 76 db ans = 1e18; 77 point x = pvec(a); 78 line s0 = {b[0],b[0]+x}; 79 point xx = intersection(s0,a); 80 if(dot_inplane_in(xx,a)) 81 ans = min(ans,(xx-b[0]).length()); 82 line s1 = {b[1],b[1]+x}; 83 xx = intersection(s1,a); 84 if(dot_inplane_in(xx,a)) 85 ans = min(ans,(xx-b[1]).length()); 86 line s2 = {b[2],b[2]+x}; 87 xx = intersection(s2,a); 88 if(dot_inplane_in(xx,a)) 89 ans = min(ans,(xx-b[2]).length()); 90 return ans; 91 } 92 db ptoline(point p,line l){//点到直线的距离 93 return (det(p-l[0],l[1]-l[0])/dis(l[0],l[1])).length(); 94 } 95 point proj(point k1,point k2,point q){ // q 到直线 k1,k2 的投影 96 point k=k2-k1; return k1+k*(dot(q-k1,k)/k.len2()); 97 } 98 int inmid(point k1,point k2,point k3){ 99 return inmid(k1.x,k2.x,k3.x)&&inmid(k1.y,k2.y,k3.y)&&inmid(k1.z,k2.z,k3.z); 100 } 101 db slove2(point a,line b){//点到线段的距离 102 point v1=b[1]-b[0],v2=a-b[0],v3=a-b[1]; 103 if(sign(dot(v1,v2)<0))return v2.length(); 104 else if(sign(dot(v1,v3))>0)return v3.length(); 105 else return det(v1,v2).length()/v1.length(); 106 // point x = proj(b[0],b[1],a); 107 // if(dot_online_in(x,b)) 108 // return (a-x).length(); 109 // return min((a-b[0]).length(),(a-b[1]).length()); 110 } 111 db dis_point_plane(point a,plane p){ 112 point xx = pvec(p); 113 line l = {a,a+xx}; 114 point xxx = intersection(l,p); 115 if(dot_inplane_in(xxx,p)) 116 return (a-xxx).length(); 117 return min(min(slove2(a,{p[0],p[1]}),slove2(a,{p[1],p[2]})),slove2(a,{p[2],p[0]})); 118 } 119 db slove3(line a,plane b){ 120 point l=a[0],r=a[1]; 121 for(int i=0;i<90;i++){ 122 point mid = (l+r)/2,lm=(l+mid)/2,rm=(r+mid)/2; 123 if(dis_point_plane(lm,b)<=dis_point_plane(rm,b)){ 124 r=rm; 125 } else 126 l=lm; 127 } 128 return dis_point_plane(l,b); 129 } 130 double x,y,z; 131 int main(){ 132 scanf("%d",&t); 133 while (t--){ 134 for(int i=1;i<=6;i++){ 135 scanf("%lf%lf%lf",&x,&y,&z); 136 p[i]={x,y,z}; 137 } 138 plane a = {p[1],p[2],p[3]},b = {p[4],p[5],p[6]}; 139 // if(check(a,b)){ 140 // printf("%.11f\n",0.0); 141 // }else{ 142 db ans = 1e18; 143 // ans = min(ans,slove(a,b)); 144 // ans = min(ans,slove(b,a)); 145 // for(int i=0;i<3;i++) { 146 // for(int k=0;k<3;k++) { 147 // point l=a[i],r=a[(i+1)%3],lm,rm; 148 // line t=(line){b[k],b[(k+1)%3]}; 149 // for(int j=0;j<90;j++) { 150 // lm=l+(r-l)/3,rm=r-(r-l)/3; 151 // if(slove2(lm,t)<=slove2(rm,t)) r=rm; 152 // else l=lm; 153 // } 154 // ans=std::min(ans,slove2(r,t)); 155 // } 156 // } 157 for(int i=0;i<3;i++){ 158 ans = min(ans,slove3({a[i],a[(i+1)%3]},b)); 159 ans = min(ans,slove3({b[i],b[(i+1)%3]},a)); 160 } 161 printf("%.11f\n",ans); 162 // } 163 } 164 } 165 /** 166 1 167 0 0 0 1 0 0 0 1 0 168 2 2 1 3 2 1 2 3 1 169 */
标签：平面 online inter return space -- clu src include
原文地址：https://www.cnblogs.com/MXang/p/10746442.html
