[计算几何]平面最近点对
【算法分析】:
暴力枚举O(N^2)
分治解决O(NlogN)
传送门:http://blogold.chinaunix.net/u2/63316/showart_2334130.html
没有多想……按照算法的原理自己写了一个……效率还是蛮高的!10W的数据秒杀,暴力则跑了2min多种才出解。
正确性无从而知,因此拿暴力对拍了一下,还好自测数据都是对的。
【Code】:
const maxn=100001; type point=record x,y:extended; end; var q,q1,q2:array[0..maxn]of longint; p,pp:array[0..maxn]of point; n,i:longint; function min(x,y:extended):extended; begin if x<y then exit(x); exit(y); end; procedure swap(var a,b:point); var t:point; begin t:=a; a:=b; b:=t; end; procedure qsort(l,r:longint); var i,j:longint; x:extended; begin i:=l; j:=r; x:=p[(l+r)>>1].x; repeat while p[i].x<x do inc(i); while p[j].x>x do dec(j); if i<=j then begin swap(p[i],p[j]); inc(i); dec(j); end; until i>j; if l<j then qsort(l,j); if i<r then qsort(i,r); end; function dist(p0,p1:point):extended; begin exit(sqrt(sqr(p0.x-p1.x)+sqr(p0.y-p1.y))); end; function solve(l,r:longint):extended; var tl,tr,tx,mx:extended; m,i,j,k:longint; begin if l=r then exit(1e99); //区间大小为1 if l+1=r then begin if p[l].y>p[r].y then swap(p[l],p[r]); exit(dist(p[l],p[r])); end; //区间大小为2 m:=(l+r)>>1; mx:=(p[m].x+p[m+1].x)/2; tl:=solve(l,m); //左区间 tr:=solve(m+1,r); //右区间 tx:=min(tl,tr); q1[0]:=0; q2[0]:=0; q[0]:=0; for i:=l to m do if p[i].x>=mx-tx then begin inc(q1[0]); q1[q1[0]]:=i; end; //筛出左边可能的点 for i:=m+1 to r do if p[i].y<=mx+tx then begin inc(q2[0]); q2[q2[0]]:=i; end; //筛出右边肯能的点 j:=1; k:=1; for i:=1 to q1[0]+q2[0] do begin if (j>q1[0]) then q[i]:=q2[k] else if (k>q2[0]) then q[i]:=q1[j] else if (p[q1[j]].y<p[q2[k]].y) then q[i]:=q1[j] else q[i]:=q2[k]; if (q[i]=q1[j]) and (j<=q1[0]) then inc(j) else inc(k); end; //合并到一个线性表中 for i:=1 to q1[0]+q2[0] do for j:=1 to 6 do if (i+j<=q1[0]+q2[0]) then tx:=min(tx,dist(p[q[i]],p[q[i+j]])); //枚举周围6个的可能的解 j:=l; k:=m+1; for i:=l to r do begin if (j>m) then pp[i]:=p[k] else if (k>r) then pp[i]:=p[j] else if (p[j].y<p[k].y) then pp[i]:=p[j] else pp[i]:=p[k]; if (pp[i].x=p[j].x) and (pp[i].y=p[j].y) and (j<=m) then inc(j) else inc(k); end; for i:=l to r do p[i]:=pp[i]; //将区间按照y坐标排序 exit(tx); end; begin readln(n); for i:=1 to n do readln(p[i].x,p[i].y); qsort(1,n); //将区间按照x坐标排序 writeln(solve(1,n):0:4); end.
posted on 2011-03-17 11:22 Skywalker_Q 阅读(344) 评论(0) 收藏 举报
浙公网安备 33010602011771号