Range Minimum Query in O(1) con O(n) precalcolo
Prerequisiti:
- Conoscere le sparse table
- Sapere cos’e’ l’LCA (Lowest Common Ancestor) di due nodi in un albero
Problema:
Viene dato un array composto da N elementi ordinabili, una query e’ composta da due indici l e r, la risposta alla query e’ l’indice dell’elemento minimo in [l,r]
Come molti sanno con le sparse table si possono fare RMQ (Range Minimum Query) su un array in O(1) con O(N\log(N)) precalcolo.
Viene presentato un algoritmo per poter fare RMQ in tempo costante con precalcolo lineare.
I principali 3 step dell’algoritmo sono:
Ridurre RMQ a LCA
Per ridurre RMQ a LCA, basta costruire un cartesian tree tale che la x dei nodi del cartesian tree sia l’indice i dell’elemento e la y sia il valore a_i dell’elemento.
Un cartesian tree e’ un albero binario in cui ogni nodo ha una x e una y.
Per ogni nodo con coordinate x e y, se x_l, y_l e x_r, y_r sono le coordinate di figlo sinistro e destro, allora x_l < x < x_r (proprieta’ di bst) e y < y_l\; e\; y < y_r (proprieta’ di heap)
Per costruire il cartesian tree in tempo lineare, possiamo scorrere l’array iniziale, tenendo uno stack dei nodi raggiungibili dalla radice senza andare mai a sinistra, e ogni volta che aggiungiamo un nodo, poppiamo dallo stack tutti i nodi con y maggiore, l’ultimo che abbiamo poppato diventa figlio sinistro del nuovo nodo e il nuovo nodo diventa figlio destro del primo che non abbiamo poppato.
Dimostro che l’LCA tra due nodi a e b in questo cartesian tree equivale al minimo in [a,b] nell’array originale.
E’ facile notare che l’LCA di a e b deve essere in [a,b], questo segue direttamente dalle proprieta’ di bst del cartesian tree.
L’LCA tra a e b e’ l’elemento con y minima, questo perche’, data l’ultima proprieta’, questo LCA deve essere ancestor di tutti i nodi in [a,b], e data la proprieta’ di heap del cartesian tree, questo puo’ solo essere il nodo con y minima.
Ridurre LCA a una particolare istanza di RMQ
Questo e’ un metodo molto conosciuto per risolvere LCA su un albero, consiste nel linearizzare l’albero con un euler-tour (dfs aggiungendo un nodo ogni volta che ci passi, anche quando “torni in su” dai figli).
Una volta linearizzato l’albero, l’LCA tra due nodi corrisponde alla RMQ sulle altezze fra la prima occorrenza dei due nodi nell’albero linearizzato.
Dimostrazione: chiamiamo a e b i due nodi e l il loro LCA.
Tra la prima e l’ultima occorrenza di un nodo nell’albero linearizzato e’ presente tutto e solo il suo sottoalbero, quindi se l e’ presente tra le prime occorrenze di a e b, e’ sicuramente il nodo con altezza minima.
l e’ presente tra le prime occorrenze di a e b, questo perche’ essendo il loro LCA, a e b sono in sottoalberi di figli diversi di l.
Risolvere il caso particolare di RMQ
Molti algoritmi per l’LCA si fermano qui e usano un algoritmo generico per RMQ.
Invece l’algoritmo di Farach-Colton e Bender sfrutta il fatto che nell’albero linearizzato delle altezze la differenza di altezze fra due elementi successivi e’ o +1 o -1.
Per risolvere questa istanza dividiamo l’array in \frac{2N}{\log(N)} blocchi grandi \frac{\log(N)}{2}.
Sull’array dei minimi dei blocchi costruiamo una sparse table in \frac{2N}{\log(N)} \log\left(\frac{2N}{\log(N)}\right) = \frac{2N}{\log(N)} \left(1 + \log\left(\frac{N}{\log(N)}\right)\right) \leq \frac{2N}{\log(N)} + 2N = O(N)
Per avere invece il minimo in range che non prendono blocchi interi, ci serve poter avere il minimo all’interno di un blocco.
Per fare questo, dato che la differenza e’ sempre +1 o -1, notiamo che esistono 2^{\frac{\log(N)}{2}-1} = \frac{\sqrt{N}}{2} blocchi diversi.
Quindi per ognuno, per ogni l e r al suo interno, possiamo precalcolarci la RMQ, questo con una complessita’ di O(\sqrt{N} \log^2(N)) = O(N)
Cosi’ per rispondere a una query l,r, possiamo usare la sparse table per avere il minimo nei blocchi interi in [l,r] e il secondo precalcolo per avere il minimo negli elementi che la sparse table non ha coperto.
Implementazione:
Questo e’ come implemento io le cose quando faccio competitive programming, ma ho messo un po’ di commenti, quindi con dell’impegno dovreste riuscire a capire
#include <bits/stdc++.h>
using namespace std;
template<class T>
class rmq {
private:
int bs; // block size
int nb; // number of blocks
int * tour; // euler tour
int * fo; // first occurence of i in eulet tour
int * bm; // block's mask
int * prem; // precomputed values for each mask
int * st; // sparse table
int * h; // heights of nodes (distance from root)
inline int& pat(int mask, int l, int r){ // to linearize prem
return prem[(mask*bs+l)*bs+r];
}
inline int& stat(int j, int i){ // to linearize sparse table
return st[j*nb+i];
}
void dfs(int p, pair<int,int>* g, int& cnt){ // euler tour
fo[p]=cnt;
tour[cnt++]=p;
if(g[p].first>=0){
h[g[p].first]=h[p]+1;
dfs(g[p].first,g,cnt);
tour[cnt++]=p;
}
if(g[p].second>=0){
h[g[p].second]=h[p]+1;
dfs(g[p].second,g,cnt);
tour[cnt++]=p;
}
}
public:
// precalcolo O(N)
rmq(T* a, int n){
int tmp;
// build cartesian tree
pair<int,int> * g = (pair<int,int>*)malloc(n*sizeof(pair<int,int>)); // cartesian tree (first is left child, second is right)
stack<int> s;
for(int i=0; i<n; ++i){
tmp=-1;
while(!s.empty() && a[s.top()]>=a[i]){
tmp = s.top();
s.pop();
}
g[i].first=tmp;
g[i].second=-1;
if(!s.empty())g[s.top()].second = i;
s.push(i);
}
while(s.size()>1)
s.pop();
// get cartesian tree's euler tour
fo = (int*)malloc(n*sizeof(int));
h = (int*)malloc(n*sizeof(int));
n=(n+n-1); // changing n to be the size of euler tour
tour = (int*)malloc(n*sizeof(int));
tmp=h[s.top()]=0;
dfs(s.top(),g,tmp);
free(g);
// set some constants
bs = ((31-__builtin_clz(n))>>1)+1; // set block size to O(log(n)/2)
nb = (n+bs-1)/bs; // set number of blocks to ceil(n/bs)
// precompute mask for each block
bm = (int*)calloc(nb,sizeof(int));
for(int i=0; i<nb; ++i){
for(int j=i*bs+1; j<(i+1)*bs; ++j){
bm[i]<<=1;
if(j<n&&h[tour[j]]>h[tour[j-1]])bm[i]|=1;
}
}
// precompute RMQ for each mask, l and r
vector<bool> computed(1<<(bs-1)); // true if mask i was already computed
prem = (int*)malloc(sizeof(int)*(1<<(bs-1))*bs*bs);
for (int i=0; i<nb; ++i) {
if(!computed[bm[i]]){
computed[bm[i]]=1;
for(int l=0; l<bs; ++l){
pat(bm[i],l,l)=l;
for(int r=l+1; r<bs; ++r){
if(i*bs+r<n&&h[tour[i*bs+pat(bm[i],l,r-1)]]>h[tour[i*bs+r]])
pat(bm[i],l,r)=r;
else pat(bm[i],l,r)=pat(bm[i],l,r-1);
}
}
}
}
// build sparse table
int sth = 33-__builtin_clz(nb); // "height" of sparse table
st = (int*)malloc(sizeof(int)*nb*sth);
for(int i=0; i<nb; ++i)
stat(0,i)=i*bs+pat(bm[i],0,bs-1);
for(int j=1; j<sth; ++j)
for(int i=0; i+(1<<j)-1<nb; ++i)
stat(j,i) = h[tour[stat(j-1,i)]] > h[tour[stat(j-1,i+(1<<(j-1)))]] ? stat(j-1,i+(1<<(j-1))):stat(j-1,i);
}
// RMQ O(1)
int operator()(int l, int r){
l=fo[l];
r=fo[r];
if(l>r)swap(l,r);
int ll=l/bs,rr=r/bs;
if(ll==rr){ // if you are fully inside a block
return tour[pat(bm[ll],l%bs,r%bs)+ll*bs];
} else {
int x = tour[pat(bm[ll],l%bs,bs-1)+ll*bs]; // inside l's block
int t = tour[pat(bm[rr],0,r%bs)+rr*bs]; // inside r's block
if(h[t]<h[x])x=t;
if(rr-ll>1){ // at least 1 full block not included with l and r
int stl=ll+1,str=rr-1;
int j = 31-__builtin_clz(str-stl+1);
t=tour[stat(j,stl)];
if(h[t]<h[x])x=t;
t=tour[stat(j,str-(1<<j)+1)];
if(h[t]<h[x])x=t;
}
return x;
}
}
};
Problemi
- preoii_biliardino
- ois_desk (si puo’ risolvere molto piu’ facilmente con la stessa complessita’)
- qualsiasi problema con RMQ o LCA