I finally figured out a nice way to implement dynamic convex hull with a reasonable constant factor! I expect this will be "well known in China" though.↵
↵
Consider the problem of finding [nested convex hulls](https://en.wikipedia.org/wiki/Convex_layers) in $O(n\log^2{n})$. The simplest way I know of is to make a convex hull data structure that supports point deletions, which is what I do here. It should be possible to extend this implementation to handle insertions as well.↵
↵
This implementation is based on this [paper](https://www.sciencedirect.com/science/article/pii/002200008190012X) by Overmars and van Leeuwen that describes a data structure for online fully dynamic convex hull in $O(\log^2{n})$ per operation. Directly implementing that is complicated and has a bad constant factor, so we'll do it a bit differently using ideas from [user:FizzyDavid,2020-04-13]'s [solution](https://codeforces.me/blog/entry/72611?#comment-569105) to [1270H](https://codeforces.me/contest/1270/problem/H).↵
↵
The data structure↵
==================↵
↵
We need a data structure for points in the plane that can support↵
↵
- deleting a point in $O(\log^2{n})$↵
↵
- returning the convex hull of $k$ points in $O(k\log^2{n})$ time↵
↵
The left and right hulls are handled separately. From now on, we'll only consider the left hull.↵
↵
We can build the left hull using a divide and conquer strategy:↵
↵
1. Split the points in half by y-coordinate.↵
2. Recursively build the hull for both halves↵
3. Find a "bridge" that connects the two hulls. This can be done in $O(\log{n})$ with a binary search.↵
4. Merge the two hulls with the bridge. Edges covered by the bridge are discarded.↵
↵
To make this dynamic, we can turn this divide and conquer into a segment tree.↵
↵
The paper explicitly stores the convex hull in a balanced binary search tree that supports split and concatenate in $O(\log^2{n})$. We'll instead implicitly represent the convex hull of each node in a segment tree.↵
↵
Each node of the segment tree must store↵
↵
- Pointers to its left and right children↵
↵
- The bridge between the convex hulls of its children↵
↵
- The bottommost point of its leaves↵
↵
In a leaf node, we just pretend the bridge is between the point and itself.↵
↵
Interestingly, computing the bridge can be done by querying its children in $O(\log{n})$.↵
↵
Finding the bridge efficiently↵
==================↵
↵
A node of the segment tree must be able to compute the bridge of its children's convex hulls efficiently.↵
↵
The paper has $10$ cases, but I think we only need $6$.↵
↵
We have two segment tree nodes $p$ and $q$ that we want to find the bridge between. Let $a$ and $b$ be the points of the bridge of the first node and $c$ and $d$ be the points of the bridge of the second node.↵
↵
Repeat the following until both $p$ and $q$ are leaves.↵
↵
- If $a\ne b$ ($p$ is not a leaf) and $a,b,c$ is in counterclockwise order, we can discard $b$ and all points above it. Set $p$ to its first child.↵
↵
- If $c\ne d$ and $b,c,d$ is in counterclockwise order, we can discard $c$ and all points below it. Set $q$ to its second child.↵
↵
- Otherwise, if $a=b$, we can discard $d$ and all points above it. Set $q$ to its first child.↵
↵
- Otherwise, if $c=d$, we can discard $a$ and all points below it. Set $p$ to its second child.↵
↵
- If we get here, $a$, $b$, $c$, $d$ are clockwise.Consider anyThe horizontal line through the firbottommost point of $q$ separates the two hulls. Depending on which side of the line the intersection of lines $ab$ and $cd$ are, we can either discard $a$ and all points below it, or $d$ and all points above it. We either set $p$ to its second child or set $q$ to its first child.↵
↵
After every iteration, we move either $p$ or $q$ to its children, so this process must terminate in $O(\log{n})$ iterations.↵
↵
Extracting the convex hull↵
==================↵
↵
To compute the convex hull, we define a recursive function that does the following:↵
↵
- Given a node and two points $l$ and $r$ on the convex hull of the node, output the points on the convex hull between them, inclusive.↵
↵
We just need to handle a few cases.↵
↵
- If the bridge is contained between $l$ and $r$, we recurse on both children↵
- If the bridge is fully covered, we recurse on the child that isn't fully covered.↵
↵
It is impossible for the bridge to be partially covered.↵
↵
Deletions↵
==================↵
Deletions are kind of annoying because there is nothing that naturally acts as identity. I just removed all nodes with one child from the tree to avoid adding extra cases to the bridge-finding.↵
↵
Implementation↵
==============↵
↵
<spoiler summary="Code">↵
~~~~~↵
//"segment tree"-style implementation of online decremental dynamic convex hull↵
//Based on paper "Maintenance of configurations in the plane" by Overmars and van Leeuwen, with some modifications↵
//This implementation only supports efficient (O(log^2n)) deletion of points, which is enough to compute nested convex hulls in O(nlog^2n).↵
//It should be possible to modify this to support insertions as well.↵
↵
//Assumes all points are distinct↵
//Assumes coordinates are at most 10^6↵
//Can handle 200000 points in a few seconds.↵
↵
#include <cstdio>↵
#include <vector>↵
#include <map>↵
#include <stdint.h>↵
#include <algorithm>↵
#include <cassert>↵
#include <set>↵
↵
struct Point{↵
int64_t x,y;↵
Point operator-(Point p)const{↵
return {x-p.x,y-p.y};↵
}↵
int64_t cross(Point p)const{↵
return x*p.y-y*p.x;↵
}↵
int64_t dot(Point p)const{↵
return x*p.x+y*p.y;↵
}↵
bool operator<(Point p)const{↵
if(y!=p.y) return y<p.y;↵
return x<p.x;↵
}↵
bool operator==(Point p)const{↵
return x==p.x&&y==p.y;↵
}↵
Point operator-()const{↵
return {-x,-y};↵
}↵
};↵
↵
int64_t cross(Point a,Point b,Point c){↵
return (b-a).cross(c-a);↵
}↵
↵
↵
class LeftHull{↵
std::vector<Point> ps;↵
struct Node{↵
int bl,br;↵
int L,R;↵
int lchd,rchd;↵
};↵
std::vector<Node> nodes;↵
int root;↵
bool isleaf(int w){↵
return nodes[w].lchd==-1&&nodes[w].rchd==-1;↵
}↵
void pull(int w){↵
assert(!isleaf(w));↵
int l=nodes[w].lchd,r=nodes[w].rchd;↵
int64_t split_y=ps[nodes[r].L].y;↵
while(!isleaf(l)||!isleaf(r)){↵
int a=nodes[l].bl,b=nodes[l].br,↵
c=nodes[r].bl,d=nodes[r].br;↵
if(a!=b && cross(ps[a],ps[b],ps[c])>0){↵
l=nodes[l].lchd;↵
}else if(c!=d && cross(ps[b],ps[c],ps[d])>0){↵
r=nodes[r].rchd;↵
}else if(a==b){↵
r=nodes[r].lchd;↵
}else if(c==d){↵
l=nodes[l].rchd;↵
}else{↵
int64_t s1=cross(ps[a],ps[b],ps[c]);↵
int64_t s2=cross(ps[b],ps[a],ps[d]);↵
assert(s1+s2>=0);↵
if(s1+s2==0||s1*ps[d].y+s2*ps[c].y<split_y*(s1+s2)){↵
l=nodes[l].rchd;↵
}else{↵
r=nodes[r].lchd;↵
}↵
}↵
}↵
nodes[w].bl=nodes[l].L;↵
nodes[w].br=nodes[r].L;↵
}↵
void build(int w,int L,int R){↵
nodes[w].L=L;↵
nodes[w].R=R;↵
if(R-L==1){↵
nodes[w].lchd=nodes[w].rchd=-1;↵
nodes[w].bl=nodes[w].br=L;↵
}else{↵
int M=(L+R)/2;↵
nodes[w].lchd=w+1;↵
nodes[w].rchd=w+2*(M-L);↵
build(nodes[w].lchd,L,M);↵
build(nodes[w].rchd,M,R);↵
pull(w);↵
}↵
}↵
int erase(int w,int L,int R){↵
if(R<=nodes[w].L||L>=nodes[w].R) return w;↵
if(L<=nodes[w].L&&R>=nodes[w].R) return -1;↵
nodes[w].lchd=erase(nodes[w].lchd,L,R);↵
nodes[w].rchd=erase(nodes[w].rchd,L,R);↵
if(nodes[w].lchd==-1) return nodes[w].rchd;↵
if(nodes[w].rchd==-1) return nodes[w].lchd;↵
pull(w);↵
return w;↵
}↵
//only works for whole hull↵
void get_hull(int w,int l,int r,std::vector<int>& res){↵
if(isleaf(w)){↵
res.push_back(nodes[w].L);↵
}else if(r<=nodes[w].bl){↵
get_hull(nodes[w].lchd,l,r,res);↵
}else if(l>=nodes[w].br){↵
get_hull(nodes[w].rchd,l,r,res);↵
}else{↵
assert(l<=nodes[w].bl&&nodes[w].br<=r);↵
get_hull(nodes[w].lchd,l,nodes[w].bl,res);↵
get_hull(nodes[w].rchd,nodes[w].br,r,res);↵
}↵
}↵
public:↵
LeftHull(const std::vector<Point>& ps):ps(ps),nodes(ps.size()*2),root(0){↵
build(0,0,ps.size());↵
}↵
std::vector<int> get_hull(){↵
if(root==-1) return {};↵
std::vector<int> res;↵
get_hull(root,0,ps.size()-1,res);↵
return res;↵
}↵
void erase(int L){↵
root=erase(root,L,L+1);↵
}↵
};↵
↵
std::vector<Point> ps;↵
std::map<Point,int> id;↵
int layer[1000005];↵
int ans[1000005];↵
↵
int main(){↵
int N;↵
scanf("%d",&N);↵
for(int i=0;i<N;i++){↵
int X,Y;↵
scanf("%d %d",&X,&Y);↵
ps.push_back({X,Y});↵
id[{X,Y}]=i;↵
}↵
↵
std::sort(ps.begin(),ps.end());↵
struct LeftHull left(ps);↵
std::reverse(ps.begin(),ps.end());↵
for(auto& p:ps){↵
p=-p;↵
}↵
struct LeftHull right(ps);↵
for(auto& p:ps){↵
p=-p;↵
}↵
std::reverse(ps.begin(),ps.end());↵
for(int l=1,cnt=0;cnt<N;l++){↵
std::set<int> hull;↵
for(int i:left.get_hull()){↵
hull.insert(i);↵
}↵
for(int i:right.get_hull()){↵
hull.insert(N-1-i);↵
}↵
for(int i:hull){↵
assert(!layer[i]);↵
cnt++;↵
layer[i]=l;↵
left.erase(i);↵
right.erase(N-1-i);↵
}↵
}↵
for(int i=0;i<N;i++){↵
ans[id[ps[i]]]=layer[i];↵
}↵
for(int i=0;i<N;i++){↵
printf("%d\n",ans[i]);↵
}↵
}↵
~~~~~↵
</spoiler>↵
↵
↵
Questions↵
==================↵
Another problem mentioned in the paper is maintaining points in the plane such that no point has larger $x$ and $y$ coordinates.↵
↵
[1270H](https://codeforces.me/contest/1270/problem/H) mentioned above is similar.↵
↵
What are other problems solvable this way?↵
↵
Can all problems solvable by Overmars and van Leeuwen's technique be solved by this type of segment tree?↵
↵
Is it possible to optimize all these to $O(n\log{n})$? Apparently it is [possible](https://arxiv.org/abs/1902.11169) for dynamic convex hull at least.↵
↵
↵
↵
Consider the problem of finding [nested convex hulls](https://en.wikipedia.org/wiki/Convex_layers) in $O(n\log^2{n})$. The simplest way I know of is to make a convex hull data structure that supports point deletions, which is what I do here. It should be possible to extend this implementation to handle insertions as well.↵
↵
This implementation is based on this [paper](https://www.sciencedirect.com/science/article/pii/002200008190012X) by Overmars and van Leeuwen that describes a data structure for online fully dynamic convex hull in $O(\log^2{n})$ per operation. Directly implementing that is complicated and has a bad constant factor, so we'll do it a bit differently using ideas from [user:FizzyDavid,2020-04-13]'s [solution](https://codeforces.me/blog/entry/72611?#comment-569105) to [1270H](https://codeforces.me/contest/1270/problem/H).↵
↵
The data structure↵
==================↵
↵
We need a data structure for points in the plane that can support↵
↵
- deleting a point in $O(\log^2{n})$↵
↵
- returning the convex hull of $k$ points in $O(k\log^2{n})$ time↵
↵
The left and right hulls are handled separately. From now on, we'll only consider the left hull.↵
↵
We can build the left hull using a divide and conquer strategy:↵
↵
1. Split the points in half by y-coordinate.↵
2. Recursively build the hull for both halves↵
3. Find a "bridge" that connects the two hulls. This can be done in $O(\log{n})$ with a binary search.↵
4. Merge the two hulls with the bridge. Edges covered by the bridge are discarded.↵
↵
To make this dynamic, we can turn this divide and conquer into a segment tree.↵
↵
The paper explicitly stores the convex hull in a balanced binary search tree that supports split and concatenate in $O(\log^2{n})$. We'll instead implicitly represent the convex hull of each node in a segment tree.↵
↵
Each node of the segment tree must store↵
↵
- Pointers to its left and right children↵
↵
- The bridge between the convex hulls of its children↵
↵
- The bottommost point of its leaves↵
↵
In a leaf node, we just pretend the bridge is between the point and itself.↵
↵
Interestingly, computing the bridge can be done by querying its children in $O(\log{n})$.↵
↵
Finding the bridge efficiently↵
==================↵
↵
A node of the segment tree must be able to compute the bridge of its children's convex hulls efficiently.↵
↵
The paper has $10$ cases, but I think we only need $6$.↵
↵
We have two segment tree nodes $p$ and $q$ that we want to find the bridge between. Let $a$ and $b$ be the points of the bridge of the first node and $c$ and $d$ be the points of the bridge of the second node.↵
↵
Repeat the following until both $p$ and $q$ are leaves.↵
↵
- If $a\ne b$ ($p$ is not a leaf) and $a,b,c$ is in counterclockwise order, we can discard $b$ and all points above it. Set $p$ to its first child.↵
↵
- If $c\ne d$ and $b,c,d$ is in counterclockwise order, we can discard $c$ and all points below it. Set $q$ to its second child.↵
↵
- Otherwise, if $a=b$, we can discard $d$ and all points above it. Set $q$ to its first child.↵
↵
- Otherwise, if $c=d$, we can discard $a$ and all points below it. Set $p$ to its second child.↵
↵
- If we get here, $a$, $b$, $c$, $d$ are clockwise.
↵
After every iteration, we move either $p$ or $q$ to its children, so this process must terminate in $O(\log{n})$ iterations.↵
↵
Extracting the convex hull↵
==================↵
↵
To compute the convex hull, we define a recursive function that does the following:↵
↵
- Given a node and two points $l$ and $r$ on the convex hull of the node, output the points on the convex hull between them, inclusive.↵
↵
We just need to handle a few cases.↵
↵
- If the bridge is contained between $l$ and $r$, we recurse on both children↵
- If the bridge is fully covered, we recurse on the child that isn't fully covered.↵
↵
It is impossible for the bridge to be partially covered.↵
↵
Deletions↵
==================↵
Deletions are kind of annoying because there is nothing that naturally acts as identity. I just removed all nodes with one child from the tree to avoid adding extra cases to the bridge-finding.↵
↵
Implementation↵
==============↵
↵
<spoiler summary="Code">↵
~~~~~↵
//"segment tree"-style implementation of online decremental dynamic convex hull↵
//Based on paper "Maintenance of configurations in the plane" by Overmars and van Leeuwen, with some modifications↵
//This implementation only supports efficient (O(log^2n)) deletion of points, which is enough to compute nested convex hulls in O(nlog^2n).↵
//It should be possible to modify this to support insertions as well.↵
↵
//Assumes all points are distinct↵
//Assumes coordinates are at most 10^6↵
//Can handle 200000 points in a few seconds.↵
↵
#include <cstdio>↵
#include <vector>↵
#include <map>↵
#include <stdint.h>↵
#include <algorithm>↵
#include <cassert>↵
#include <set>↵
↵
struct Point{↵
int64_t x,y;↵
Point operator-(Point p)const{↵
return {x-p.x,y-p.y};↵
}↵
int64_t cross(Point p)const{↵
return x*p.y-y*p.x;↵
}↵
int64_t dot(Point p)const{↵
return x*p.x+y*p.y;↵
}↵
bool operator<(Point p)const{↵
if(y!=p.y) return y<p.y;↵
return x<p.x;↵
}↵
bool operator==(Point p)const{↵
return x==p.x&&y==p.y;↵
}↵
Point operator-()const{↵
return {-x,-y};↵
}↵
};↵
↵
int64_t cross(Point a,Point b,Point c){↵
return (b-a).cross(c-a);↵
}↵
↵
↵
class LeftHull{↵
std::vector<Point> ps;↵
struct Node{↵
int bl,br;↵
int L,R;↵
int lchd,rchd;↵
};↵
std::vector<Node> nodes;↵
int root;↵
bool isleaf(int w){↵
return nodes[w].lchd==-1&&nodes[w].rchd==-1;↵
}↵
void pull(int w){↵
assert(!isleaf(w));↵
int l=nodes[w].lchd,r=nodes[w].rchd;↵
int64_t split_y=ps[nodes[r].L].y;↵
while(!isleaf(l)||!isleaf(r)){↵
int a=nodes[l].bl,b=nodes[l].br,↵
c=nodes[r].bl,d=nodes[r].br;↵
if(a!=b && cross(ps[a],ps[b],ps[c])>0){↵
l=nodes[l].lchd;↵
}else if(c!=d && cross(ps[b],ps[c],ps[d])>0){↵
r=nodes[r].rchd;↵
}else if(a==b){↵
r=nodes[r].lchd;↵
}else if(c==d){↵
l=nodes[l].rchd;↵
}else{↵
int64_t s1=cross(ps[a],ps[b],ps[c]);↵
int64_t s2=cross(ps[b],ps[a],ps[d]);↵
assert(s1+s2>=0);↵
if(s1+s2==0||s1*ps[d].y+s2*ps[c].y<split_y*(s1+s2)){↵
l=nodes[l].rchd;↵
}else{↵
r=nodes[r].lchd;↵
}↵
}↵
}↵
nodes[w].bl=nodes[l].L;↵
nodes[w].br=nodes[r].L;↵
}↵
void build(int w,int L,int R){↵
nodes[w].L=L;↵
nodes[w].R=R;↵
if(R-L==1){↵
nodes[w].lchd=nodes[w].rchd=-1;↵
nodes[w].bl=nodes[w].br=L;↵
}else{↵
int M=(L+R)/2;↵
nodes[w].lchd=w+1;↵
nodes[w].rchd=w+2*(M-L);↵
build(nodes[w].lchd,L,M);↵
build(nodes[w].rchd,M,R);↵
pull(w);↵
}↵
}↵
int erase(int w,int L,int R){↵
if(R<=nodes[w].L||L>=nodes[w].R) return w;↵
if(L<=nodes[w].L&&R>=nodes[w].R) return -1;↵
nodes[w].lchd=erase(nodes[w].lchd,L,R);↵
nodes[w].rchd=erase(nodes[w].rchd,L,R);↵
if(nodes[w].lchd==-1) return nodes[w].rchd;↵
if(nodes[w].rchd==-1) return nodes[w].lchd;↵
pull(w);↵
return w;↵
}↵
//only works for whole hull↵
void get_hull(int w,int l,int r,std::vector<int>& res){↵
if(isleaf(w)){↵
res.push_back(nodes[w].L);↵
}else if(r<=nodes[w].bl){↵
get_hull(nodes[w].lchd,l,r,res);↵
}else if(l>=nodes[w].br){↵
get_hull(nodes[w].rchd,l,r,res);↵
}else{↵
assert(l<=nodes[w].bl&&nodes[w].br<=r);↵
get_hull(nodes[w].lchd,l,nodes[w].bl,res);↵
get_hull(nodes[w].rchd,nodes[w].br,r,res);↵
}↵
}↵
public:↵
LeftHull(const std::vector<Point>& ps):ps(ps),nodes(ps.size()*2),root(0){↵
build(0,0,ps.size());↵
}↵
std::vector<int> get_hull(){↵
if(root==-1) return {};↵
std::vector<int> res;↵
get_hull(root,0,ps.size()-1,res);↵
return res;↵
}↵
void erase(int L){↵
root=erase(root,L,L+1);↵
}↵
};↵
↵
std::vector<Point> ps;↵
std::map<Point,int> id;↵
int layer[1000005];↵
int ans[1000005];↵
↵
int main(){↵
int N;↵
scanf("%d",&N);↵
for(int i=0;i<N;i++){↵
int X,Y;↵
scanf("%d %d",&X,&Y);↵
ps.push_back({X,Y});↵
id[{X,Y}]=i;↵
}↵
↵
std::sort(ps.begin(),ps.end());↵
struct LeftHull left(ps);↵
std::reverse(ps.begin(),ps.end());↵
for(auto& p:ps){↵
p=-p;↵
}↵
struct LeftHull right(ps);↵
for(auto& p:ps){↵
p=-p;↵
}↵
std::reverse(ps.begin(),ps.end());↵
for(int l=1,cnt=0;cnt<N;l++){↵
std::set<int> hull;↵
for(int i:left.get_hull()){↵
hull.insert(i);↵
}↵
for(int i:right.get_hull()){↵
hull.insert(N-1-i);↵
}↵
for(int i:hull){↵
assert(!layer[i]);↵
cnt++;↵
layer[i]=l;↵
left.erase(i);↵
right.erase(N-1-i);↵
}↵
}↵
for(int i=0;i<N;i++){↵
ans[id[ps[i]]]=layer[i];↵
}↵
for(int i=0;i<N;i++){↵
printf("%d\n",ans[i]);↵
}↵
}↵
~~~~~↵
</spoiler>↵
↵
↵
Questions↵
==================↵
Another problem mentioned in the paper is maintaining points in the plane such that no point has larger $x$ and $y$ coordinates.↵
↵
[1270H](https://codeforces.me/contest/1270/problem/H) mentioned above is similar.↵
↵
What are other problems solvable this way?↵
↵
Can all problems solvable by Overmars and van Leeuwen's technique be solved by this type of segment tree?↵
↵
Is it possible to optimize all these to $O(n\log{n})$? Apparently it is [possible](https://arxiv.org/abs/1902.11169) for dynamic convex hull at least.↵
↵
↵
