注意:这是一篇个人学习笔记,如果有人因为某些原因点了进来并且要看一下,请一定谨慎地阅读,因为可能存在各种奇怪的错误,如果有人发现错误请指出谢谢!
资料待看:
https://www.cnblogs.com/yangsongyi/p/10000756.html
http://nfile.oss-cn-hangzhou.aliyuncs.com/KDTree.pdf
kdtree
将一些点组织成二叉树,对于每一个节点,有一个特定的划分维度,其左子树中点都在这一维上小于等于它,其右子树中点都在这一维上大于等于它。(等于的归于哪一边没有办法规定,解决问题的时候可能需要注意)这个划分维度一般是按层轮流划分(比如二维平面,树的第1,3,5,7,..层按x轴,第2,4,6,8,..层按y轴)
最近点对
查询离某一个点p最近的点:在kdtree上dfs,进行一些剪枝。
到某个节点时,优先向它的两个子节点中,"估计"更优的一个点走。
到任何一个节点都可以用来更新当前最优解。
然后就不要走那些一定不会产生比当前最优解更优的解的区域(也就是那些离p的最小可能距离也要大于等于当前最近距离的区域)。
怎么样”估计“更优呢?网上有方法1,每次优先走包含p点的那个半空间(https://blog.csdn.net/u013534123/article/details/80952174,https://codeforces.com/blog/entry/54080),然后如果p到用于划分的“面”的距离<=当前最短距离,才去走另外一个半空间(实际实现见代码)
1 #include<cstdio> 2 #include<ctime> 3 #include<algorithm> 4 #include<cstring> 5 #include<vector> 6 #include<cmath> 7 using namespace std; 8 #define fi first 9 #define se second 10 #define mp make_pair 11 #define pb push_back 12 typedef long long ll; 13 typedef unsigned long long ull; 14 #define cmin(a,b) (((b)<(a))?void((a)=(b)):void()) 15 template<class T> 16 inline T sqr(T x){return x*x;} 17 namespace K 18 { 19 struct N 20 { 21 double d[2]; 22 int lc,rc; 23 }d[200011],p[200011]; 24 int W;//W:当前维度 25 bool c1(const N &a,const N &b) 26 { 27 return a.d[W]<b.d[W]; 28 } 29 double dis2(const N &a,const N &b) 30 { 31 return sqr(a.d[0]-b.d[0])+sqr(a.d[1]-b.d[1]); 32 } 33 #define LC (d[u].lc) 34 #define RC (d[u].rc) 35 int build(int l,int r,int w)//w:当前维度 36 { 37 W=w; 38 int u=(l+r)>>1; 39 nth_element(p+l,p+u,p+r+1,c1); 40 d[u]=p[u]; 41 if(l!=u) LC=build(l,u-1,w^1); 42 if(u!=r) RC=build(u+1,r,w^1); 43 return u; 44 } 45 int x; 46 const double inf=1e18; 47 double ans; 48 void querymn(int u,int w) 49 { 50 if(!u) return; 51 double tmp=dis2(d[u],d[x]); 52 if(u!=x) cmin(ans,tmp); 53 int lc=LC,rc=RC; 54 if(d[x].d[w]>d[u].d[w]) swap(lc,rc); 55 querymn(lc,w^1); 56 if(sqr(d[x].d[w]-d[u].d[w])<=ans) querymn(rc,w^1); 57 } 58 } 59 int n,rt; 60 double anss; 61 int main() 62 { 63 srand(time(0)); 64 int i; 65 scanf("%d",&n); 66 { 67 using namespace K; 68 for(i=1;i<=n;++i) 69 scanf("%lf%lf",&p[i].d[0],&p[i].d[1]); 70 anss=inf; 71 rt=build(1,n,0); 72 for(i=1;i<=n;++i) 73 { 74 x=i;ans=inf; 75 querymn(rt,0); 76 anss=min(anss,ans); 77 } 78 printf("%.4f ",sqrt(anss)); 79 } 80 return 0; 81 }
还有方法2(http://hzwer.com/6954.html),维护每个区域的各维最小坐标和最大坐标,然后可以分别求出p到左右子节点的区域的某种特征(这个版本是将各个维度的最小差值直接加了起来,因为那个题是曼哈顿距离;如果是欧几里得距离,那么应该是各个维度的最小差值平方后加起来),先走较优的(这个版本就是指较小的)
可能需要特判同一个点/重合点
1 #include<cstdio> 2 #include<algorithm> 3 #include<cstring> 4 #include<vector> 5 #include<cmath> 6 using namespace std; 7 #define fi first 8 #define se second 9 #define mp make_pair 10 #define pb push_back 11 typedef long long ll; 12 typedef unsigned long long ull; 13 template<class T1,class T2> 14 inline T1 min1(T1 a,T2 b){return a<b?a:b;} 15 template<class T1,class T2> 16 inline T1 max1(T1 a,T2 b){return a<b?b:a;} 17 #define min min1 18 #define max max1 19 #define cmin(a,b) (((b)<(a))?void((a)=(b)):void()) 20 #define cmax(a,b) (((a)<(b))?void((a)=(b)):void()) 21 template<class T> 22 inline T sqr(T x){return x*x;} 23 namespace K 24 { 25 struct N 26 { 27 double d[2],mx[2],mn[2];//自身,子树min,max 28 int lc,rc; 29 }d[200011],p[200011]; 30 int W;//W:当前维度 31 bool c1(const N &a,const N &b) 32 { 33 return a.d[W]<b.d[W]; 34 } 35 double dis2(const N &a,const N &b) 36 { 37 return sqr(a.d[0]-b.d[0])+sqr(a.d[1]-b.d[1]); 38 } 39 #define LC (d[u].lc) 40 #define RC (d[u].rc) 41 void upd(int u) 42 { 43 for(int i=0;i<2;++i) 44 { 45 d[u].mx[i]=d[u].mn[i]=d[u].d[i]; 46 if(LC) 47 { 48 cmax(d[u].mx[i],d[LC].mx[i]); 49 cmin(d[u].mn[i],d[LC].mn[i]); 50 } 51 if(RC) 52 { 53 cmax(d[u].mx[i],d[RC].mx[i]); 54 cmin(d[u].mn[i],d[RC].mn[i]); 55 } 56 } 57 } 58 int build(int l,int r,int w)//w:当前维度 59 { 60 W=w; 61 int u=(l+r)>>1; 62 nth_element(p+l,p+u,p+r+1,c1); 63 d[u]=p[u]; 64 //for(int i=0;i<2;++i) 65 // d[u].mx[i]=d[u].mn[i]=d[u].d[i]; 66 if(l!=u) LC=build(l,u-1,w^1); 67 if(u!=r) RC=build(u+1,r,w^1); 68 upd(u); 69 return u; 70 } 71 int x; 72 const double inf=1e18; 73 double ans; 74 double getmn(int u)//x到u的子树区域的最小可能距离 75 { 76 double ans=0; 77 for(int i=0;i<2;++i) 78 ans+=sqr(max(d[x].d[i]-d[u].mx[i],0) 79 +max(d[u].mn[i]-d[x].d[i],0)); 80 return ans; 81 } 82 void querymn(int u) 83 { 84 double tmp=dis2(d[u],d[x]); 85 if(u!=x) cmin(ans,tmp); 86 double dl=inf,dr=inf; 87 if(LC) dl=getmn(LC); 88 if(RC) dr=getmn(RC); 89 if(dl<dr) 90 { 91 if(dl<ans) querymn(LC); 92 if(dr<ans) querymn(RC); 93 } 94 else 95 { 96 if(dr<ans) querymn(RC); 97 if(dl<ans) querymn(LC); 98 } 99 } 100 } 101 int n,rt; 102 double anss; 103 int main() 104 { 105 int i; 106 scanf("%d",&n); 107 { 108 using namespace K; 109 for(i=1;i<=n;++i) 110 scanf("%lf%lf",&p[i].d[0],&p[i].d[1]); 111 anss=inf; 112 rt=build(1,n,0); 113 for(i=1;i<=n;++i) 114 { 115 x=i;ans=inf; 116 querymn(rt); 117 anss=min(anss,ans); 118 } 119 printf("%.4f ",sqrt(anss)); 120 } 121 return 0; 122 }
(实际测试结果:方法1比方法2快很多(接近1/2),因为维护的东西少,且估价函数计算更简单;不清楚估价函数本身的优劣,应该是差不多的)
最远点对
这个一定要维护信息了...方法几乎和最近点对的方法2是一样的,只需要改一下估价函数(见下面K远点对)
K远点对
例题:[CQOI2016]K远点对 洛谷P4357 bzoj4520
几乎和最近/远点对是一样的,就是用一个堆维护一下当前的前K远点对
(对于此题由于点对是无序的,实际要求出第2*K大的点对)
1 #include<cstdio> 2 #include<algorithm> 3 #include<cstring> 4 #include<vector> 5 #include<cmath> 6 #include<queue> 7 using namespace std; 8 #define fi first 9 #define se second 10 #define mp make_pair 11 #define pb push_back 12 typedef long long ll; 13 typedef unsigned long long ull; 14 template<class T1,class T2> 15 inline T1 min1(T1 a,T2 b){return a<b?a:b;} 16 template<class T1,class T2> 17 inline T1 max1(T1 a,T2 b){return a<b?b:a;} 18 template<class T> 19 inline T abs(T a){return a>0?a:-a;} 20 #define min min1 21 #define max max1 22 #define cmin(a,b) (((b)<(a))?void((a)=(b)):void()) 23 #define cmax(a,b) (((a)<(b))?void((a)=(b)):void()) 24 template<class T> 25 inline T sqr(T x){return x*x;} 26 int n,rt;unsigned K; 27 namespace KK 28 { 29 struct N 30 { 31 ll d[2],mx[2],mn[2];//自身,子树min,max 32 int lc,rc; 33 }d[100011],p[100011]; 34 int W;//W:当前维度 35 bool c1(const N &a,const N &b) 36 { 37 return a.d[W]<b.d[W]; 38 } 39 ll dis2(const N &a,const N &b) 40 { 41 return sqr(a.d[0]-b.d[0])+sqr(a.d[1]-b.d[1]); 42 } 43 #define LC (d[u].lc) 44 #define RC (d[u].rc) 45 void upd(int u) 46 { 47 for(int i=0;i<2;++i) 48 { 49 d[u].mx[i]=d[u].mn[i]=d[u].d[i]; 50 if(LC) 51 { 52 cmax(d[u].mx[i],d[LC].mx[i]); 53 cmin(d[u].mn[i],d[LC].mn[i]); 54 } 55 if(RC) 56 { 57 cmax(d[u].mx[i],d[RC].mx[i]); 58 cmin(d[u].mn[i],d[RC].mn[i]); 59 } 60 } 61 } 62 int build(int l,int r,int w)//w:当前维度 63 { 64 W=w; 65 int u=(l+r)>>1; 66 nth_element(p+l,p+u,p+r+1,c1); 67 d[u]=p[u]; 68 if(l!=u) LC=build(l,u-1,w^1); 69 if(u!=r) RC=build(u+1,r,w^1); 70 upd(u); 71 return u; 72 } 73 int x; 74 const ll inf=0x3f3f3f3f3f3f3f3f; 75 ll getmx(int u)//x到u的子树区域的最大可能距离的平方 76 { 77 ll ans=0; 78 for(int i=0;i<2;++i) 79 ans+=sqr(max(abs(d[x].d[i]-d[u].mn[i]), 80 abs(d[u].mx[i]-d[x].d[i]))); 81 return ans; 82 } 83 priority_queue<ll,vector<ll>,greater<ll> > ans; 84 #define OK(x) ( (ans.size()<K) || (ans.top()<(x)) ) 85 void querymx(int u) 86 { 87 ll tmp=dis2(d[u],d[x]); 88 if(u!=x) 89 { 90 if(ans.size()<K) 91 ans.push(tmp); 92 else if(ans.top()<tmp) 93 ans.pop(),ans.push(tmp); 94 //printf("2t%d %lld %lld ",u,tmp,ans.top()); 95 } 96 ll dl=-inf,dr=-inf; 97 if(LC) dl=getmx(LC); 98 if(RC) dr=getmx(RC); 99 if(dl>dr) 100 { 101 if(LC && OK(dl)) querymx(LC); 102 if(RC && OK(dr)) querymx(RC); 103 } 104 else 105 { 106 if(RC && OK(dr)) querymx(RC); 107 if(LC && OK(dl)) querymx(LC); 108 } 109 } 110 /* 111 void out() 112 { 113 printf("size:%d ",int(ans.size())); 114 auto qqq=ans; 115 while(!qqq.empty()) 116 { 117 auto t=qqq.top();qqq.pop(); 118 printf("%lld ",t); 119 } 120 { 121 char s[233]; 122 scanf("%s",s); 123 } 124 } 125 */ 126 } 127 int main() 128 { 129 int i; 130 scanf("%d%u",&n,&K); 131 K*=2; 132 { 133 using namespace KK; 134 for(i=1;i<=n;++i) 135 scanf("%lld%lld",&p[i].d[0],&p[i].d[1]); 136 rt=build(1,n,0); 137 for(i=1;i<=n;++i) 138 { 139 x=i; 140 querymx(rt); 141 //out(); 142 } 143 printf("%lld ",ans.top()); 144 } 145 return 0; 146 }
改一下可以过:JZPFAR 清橙A1302 bzoj2626 洛谷P2093
(要求输出编号,注意细节)
1 #include<cstdio> 2 #include<algorithm> 3 #include<cstring> 4 #include<vector> 5 #include<cmath> 6 #include<queue> 7 using namespace std; 8 #define fi first 9 #define se second 10 #define mp make_pair 11 #define pb push_back 12 typedef long long ll; 13 typedef unsigned long long ull; 14 struct pli 15 { 16 ll fi;int se; 17 }; 18 bool operator<(const pli &a,const pli &b) 19 { 20 return (a.fi<b.fi) || (a.fi==b.fi && a.se>b.se); 21 } 22 template<class T1,class T2> 23 inline T1 min1(T1 a,T2 b){return a<b?a:b;} 24 template<class T1,class T2> 25 inline T1 max1(T1 a,T2 b){return a<b?b:a;} 26 template<class T> 27 inline T abs(T a){return a>0?a:-a;} 28 #define min min1 29 #define max max1 30 #define cmin(a,b) (((b)<(a))?void((a)=(b)):void()) 31 #define cmax(a,b) (((a)<(b))?void((a)=(b)):void()) 32 template<class T> 33 inline T sqr(T x){return x*x;} 34 int n,m,rt;unsigned K; 35 namespace KK 36 { 37 struct N 38 { 39 ll d[2],mx[2],mn[2];//自身,子树min,max 40 int lc,rc,n; 41 }d[100011],p[100011]; 42 int W;//W:当前维度 43 bool c1(const N &a,const N &b) 44 { 45 return a.d[W]<b.d[W]; 46 } 47 ll dis2(const N &a,const N &b) 48 { 49 return sqr(a.d[0]-b.d[0])+sqr(a.d[1]-b.d[1]); 50 } 51 #define LC (d[u].lc) 52 #define RC (d[u].rc) 53 void upd(int u) 54 { 55 for(int i=0;i<2;++i) 56 { 57 d[u].mx[i]=d[u].mn[i]=d[u].d[i]; 58 if(LC) 59 { 60 cmax(d[u].mx[i],d[LC].mx[i]); 61 cmin(d[u].mn[i],d[LC].mn[i]); 62 } 63 if(RC) 64 { 65 cmax(d[u].mx[i],d[RC].mx[i]); 66 cmin(d[u].mn[i],d[RC].mn[i]); 67 } 68 } 69 } 70 int build(int l,int r,int w)//w:当前维度 71 { 72 W=w; 73 int u=(l+r)>>1; 74 nth_element(p+l,p+u,p+r+1,c1); 75 d[u]=p[u]; 76 if(l!=u) LC=build(l,u-1,w^1); 77 if(u!=r) RC=build(u+1,r,w^1); 78 upd(u); 79 return u; 80 } 81 N x; 82 const ll inf=0x3f3f3f3f3f3f3f3f; 83 ll getmx(int u)//x到u的子树区域的最大可能距离的平方 84 { 85 ll ans=0; 86 for(int i=0;i<2;++i) 87 ans+=sqr(max(abs(x.d[i]-d[u].mn[i]), 88 abs(d[u].mx[i]-x.d[i]))); 89 return ans; 90 } 91 bool c2(const pli &a,const pli &b) 92 { 93 return (a.fi>b.fi) || (a.fi==b.fi && a.se<b.se); 94 } 95 pli ans[101];unsigned alen; 96 inline void push(const pli &x) 97 { 98 ans[++alen]=x; 99 push_heap(ans+1,ans+alen+1,c2); 100 } 101 inline void pop() 102 { 103 pop_heap(ans+1,ans+alen+1,c2); 104 --alen; 105 } 106 #define OK(x) ( (alen<K) || (ans[1].fi<=(x)) ) 107 void querymx(int u) 108 { 109 pli tmp=(pli){dis2(d[u],x),d[u].n}; 110 if(alen<K) 111 push(tmp); 112 else if(ans[1]<tmp) 113 pop(),push(tmp); 114 ll dl=-inf,dr=-inf; 115 if(LC) dl=getmx(LC); 116 if(RC) dr=getmx(RC); 117 if(dl>dr) 118 { 119 if(LC && OK(dl)) querymx(LC); 120 if(RC && OK(dr)) querymx(RC); 121 } 122 else 123 { 124 if(RC && OK(dr)) querymx(RC); 125 if(LC && OK(dl)) querymx(LC); 126 } 127 } 128 /* 129 void out() 130 { 131 printf("size:%d ",int(ans.size())); 132 auto qqq=ans; 133 while(!qqq.empty()) 134 { 135 auto t=qqq.top();qqq.pop(); 136 printf("%lld ",t); 137 } 138 { 139 char s[233]; 140 scanf("%s",s); 141 } 142 } 143 */ 144 } 145 int main() 146 { 147 int i; 148 scanf("%d",&n); 149 { 150 using namespace KK; 151 for(i=1;i<=n;++i) 152 { 153 scanf("%lld%lld",&p[i].d[0],&p[i].d[1]); 154 p[i].n=i; 155 } 156 rt=build(1,n,0); 157 scanf("%d",&m); 158 while(m--) 159 { 160 scanf("%lld%lld%u",&x.d[0],&x.d[1],&K); 161 alen=0; 162 querymx(rt); 163 //out(); 164 printf("%d ",ans[1].se); 165 } 166 } 167 return 0; 168 }