Loading [MathJax]/extensions/TeX/newcommand.js
\newcommand{\ord}[1]{\mathcal{O}\left(#1\right)} \newcommand{\abs}[1]{\lvert #1 \rvert} \newcommand{\floor}[1]{\lfloor #1 \rfloor} \newcommand{\ceil}[1]{\lceil #1 \rceil} \newcommand{\opord}{\operatorname{\mathcal{O}}} \newcommand{\argmax}{\operatorname{arg\,max}} \newcommand{\str}[1]{\texttt{"#1"}}

2017年4月30日 星期日

[ Steiner tree problem in graphs ] 斯坦納樹

斯坦納樹問題是一個世界知名的NP-hard問題。在圖論上的斯坦納樹是對於一張無向圖G=(V,E)以及一個點集合P \subseteq V,通常會稱P集合為terminal \; set,對於每條邊e=(u,v) \in E,令w(e)表示它的權重。我們的目標是要從G的所有子圖中找出一棵生成樹T=(V',E'),使得P \subseteq V'\sum_{e \in E'} w(e)最小。

簡單來說就是在圖G上找一棵子樹,可以把P中的點連通起來,且邊權總和最小

如果我們枚舉所有子圖,對每個子圖做最小生成樹演算法,就一定可以找到斯坦納樹,但是複雜度是\ord{(\abs E + \abs V log \abs V ) \times 2^{\abs V}},非常糟糕。

如果w(e)>0 ,e \in E,且\abs P \ll \abs V,我們可以找到一個動態規劃的方法:
dp[S][i]表示以點i為根,以S \subseteq Pterminal \; set構造出來的斯坦納樹,這樣我們最後的答案就會是dp[P][u \in P]

狀態轉移式可以寫成這樣
  • dp[S][i]=min(dp[T][j]+dp[S-T][j]+dis(i,j)\; : \; j \in V,T \subset S)
    dis(i,j)表示i \sim j的最短路徑
任兩點間的最短路徑可以用floyd在\ord {\abs V^3}預先算出來,狀態有2^{\abs P}\times \abs V個,狀態轉移為\ord{\abs V \times 枚舉子集合的時間},因此總複雜度為\ord{\abs V^3+2^{\abs P} \times \abs V^2 \times 枚舉子集合的時間 }

其中 2^{\abs P} \times 枚舉子集合的時間 只是粗略的計算,實際上它是
\sum_{i=1}^{\abs P} \binom{\abs P}{i} \times (2^i -1) \simeq \sum_{i=0}^{\abs P} \binom{\abs P}{i} \times 2^i = (1+2)^{\abs P} = 3^{\abs P}因此總複雜度可以表示為\ord{V^3+3^{\abs P} \times \abs V^2},但是這其實還可以優化,令H[j] = min(dp[T][j]+dp[S-T][j] \; : \;T \subset S)
dp[S][i]=min(H[j]+dis(i,j)\; : \;j \in \abs V)
H是可以被預先算出來的,因此總複雜度就降為\ord{\abs V^3 + \abs V 3^{\abs P}+\abs V^2 2^{\abs P}}
以下附上程式碼:
//n個點,其中r個要構成斯坦納樹
//p表示要構成斯坦納樹的點集
#define REP(i,n) for(int i=0;i<(int)n;++i)
const int MAXN=30,MAXM=8;// 0-base
const int INF=0x3f3f3f3f;
int dp[1<<MAXM][MAXN];
int g[MAXN][MAXN];//圖
void init(){memset(g,0x3f,sizeof(g));}
void add_edge(int u,int v,int w){
g[u][v]=g[v][u]=min(g[v][u],w);
}
void steiner(int n,int r,int *p){
REP(k,n)REP(i,n)REP(j,n)
g[i][j]=min(g[i][j],g[i][k]+g[k][j]);
REP(i,n)g[i][i]=0;
REP(i,r)REP(j,n)dp[1<<i][j]=g[p[i]][j];
for(int i=1;i<(1<<r);++i){
if(!(i&(i-1)))continue;
REP(j,n)dp[i][j]=INF;
REP(j,n){
int tmp=INF;
for(int s=i&(i-1);s;s=i&(s-1))
tmp=min(tmp,dp[s][j]+dp[i^s][j]);
REP(k,n)dp[i][k]=min(dp[i][k],g[j][k]+tmp);
}
}
}

有的時候圖是稀疏圖,也就是\ord V=\ord E,這種時候用這種算法效率其實不是很好,我們可以在dp的過程才用一些單源最短路徑算法算出最短路徑,這樣複雜度可以變成\ord{\abs V 3^{\abs P}+ShortestPath(G) 2^{\abs P}}其中ShortestPath(G)是在圖G中計算最短路徑的時間,用dijkstra的話是\ord{\abs E+\abs V log \abs V},這裡我用SPFA實作:
//n個點,其中r個要構成斯坦納樹
//p表示要構成斯坦納樹的點集
#define x first
#define y second
#define REP(i,n) for(int i=0;i<n;++i)
const int MAXN=30,MAXM=8;// 0-base
const int INF=0x3f3f3f3f;
int dp[1<<MAXM][MAXN];
vector<pair<int,int> > g[MAXN];//圖
void init(int n){REP(i,n)g[i].clear();}
void add_edge(int u,int v,int w){
g[u].push_back({v,w});
g[v].push_back({u,w});
}
queue<int> Q;
bool inq[MAXN];
void SPFA(int *d){
while(Q.size()){
int u=Q.front();Q.pop(),inq[u]=0;
for(auto &i:g[u]){
if(d[i.x]>d[u]+i.y){
d[i.x]=d[u]+i.y;
if(inq[i.x])continue;
Q.push(i.x),inq[i.x]=1;
}
}
}
}
void steiner(int n,int r,int *p){
memset(inq,0,sizeof(inq));
memset(dp,0x3f,sizeof(dp));
REP(i,r)dp[1<<i][p[i]]=0;
for(int i=1;i<(1<<r);++i){
REP(j,n){
for(int s=i&(i-1);s;s=i&(s-1))
dp[i][j]=min(dp[i][j],dp[s][j]+dp[i^s][j]);
if(dp[i][j]<INF)Q.push(j),inq[j]=1;
}
SPFA(dp[i]);
}
}

沒有留言:

張貼留言