0/1分數規劃？我不會啊！

時間  2019-11-17
標籤 1分 規劃 不會 简体版
原文   原文鏈接

最近入門了 01 分數規劃，因而寫篇博客記念html

(若是cnblogs上的格式您實在受不了，請轉至個人洛谷blog，而且那裏的講解會稍微詳細一些，另外，更新內容也在個人洛谷博客中)ios

分數規劃是一項不經常使用到的（但還蠻實用的）算法，通常來說就是求一個最優比率。git

分數規劃的簡單介紹

分數規劃顧名思義就是求一個分數表達式的最大（小）值。算法

好比說有 n 個物品，每一個物品有兩個權值 a 和 b ，而後讓你選出任意件數（但可能會有限制）的物品，使得兩個權值和間的比值最大，即求 \(\dfrac{\sum_{i=1}^{k} a[i]}{\sum_{j=1}^{k} b[j]}\) （在這裏 \(1-k\) 爲挑出的 k 件物品）的最大值，而後對選擇物品方面給出必定的限制條件，那麼一道 0/1 分數規劃的題目就完成了網絡

分數規劃的求解方法

分數規劃有兩種求解方法，一種是 二分法，而另外一種則是 Dinkelbach 算法，通常來說咱們選用第一種方法進行分數規劃求解函數

1.二分法

問題同上，求 \(\dfrac{\sum_{i=1}^{k} a[i]}{\sum_{j=1}^{k} b[j]}\) 的最大值，而後加上一個限制：\(k>=W\)post

咱們令 \(sum=\sum_{i=1}^{k} a[i] ,\ tot=\sum_{j=1}^{k} b[j]\) ，\(k>=W\)優化

而後假設問題中的最優解爲 \(ans\) ，那麼必然有：
\[\dfrac{sum}{tot}<=ans\]ui

移項得：
\[sum<=ans\times tot\]spa

繼續移就獲得：

\[sum-ans\times tot<=0\]

這樣轉化有什麼用呢？那咱們嘗試將 \(sum\) 和 \(tot\) 帶回去，就能夠獲得這麼一個式子：

\[\sum_{i=1}^{k} (a[i]-b[i]\times ans) <=0\]

這個式子不難理解，就是把總體的貢獻轉化爲了單件物品的貢獻。

那麼咱們只須要二分這個 ans， 計算出每件物品的 \(a-b\times ans\)，而後排個序，貪心取前 \(W\) 個加起來，看看最後的值是否 \(<=0\) ，而後就能夠根據結果移動左右邊界了。

2.Dinkelbach 算法

這個算法其實就是以二分函數圖像的單調性爲原理，利用了check中計算得出的結果，經過改變二分枚舉的中間點在圖像上的位置來優化二分。

而後先給一張圖：

在這裏咱們看到每次枚舉的中間點都在一條直線上（但單一的直線並非二分的函數圖像），那麼對於枚舉出的中間點其實能夠不着急扔，先計算出當前點所在的直線，而後找到這條直線在 x 軸上的截距，用這個截距做爲下一次二分的中間點，這樣可能會更快逼近正確答案

儘管理解起來不難，但仍是配套圖吧。下面是算法實現的圖組

這裏顯然紅點就是答案

這裏咱們第一次能夠直接二分中點，獲得一個初始點

而後咱們過這個獲得的點作一條直線，與函數圖像相切，

接着咱們將該直線與 \(x\) 軸交點的 \(x\) 座標做爲新點的 \(x\) 座標（好像離答案更遠了啊）

而後新點重複上述過程，作直線

這時咱們發現已經獲得答案了

然鵝這個算法的效果卻不見的有多麼優秀，由於二分的函數圖像有多是坑坑窪窪的，因此有時這個算法跑數據的時間反而比二分大，可是若是圖像是較爲光滑的弧線，或許\(Dinkelbach\) 算法就能充分展示它的優點了

至於其餘的地方（包括證實），和二分其實差很少，就不囉嗦了

分數規劃的有機結合

分數規劃通常來說不會單獨成題，通常來說有如下幾種形式：

0.不與任何算法結合，即分數規劃裸題

1.與\(0/1\)揹包結合，即最優比率揹包

2.與生成樹結合，即最優比率生成樹

3.與負環斷定結合，即最優比率環

4.與網絡流結合，即最大密度子圖

5.與費用流結合，即最優比率流（這個是我亂叫的）

6.與其餘的各類帶選擇的算法亂套，即最優比率啥啥的...

對於上面的後三個結合圖論算法，其實相似是把供選擇的物品變爲了圖中的邊，而後限制條件求解最優比率。

0.分數規劃裸題

這個在介紹二分算法的時候講過了

題目

Dropping tests

代碼

//by Judge
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define eps 1e-4
using namespace std;
const int M=1005;
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
char buf[1<<21],*p1=buf,*p2=buf;
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,k; double a[M],b[M],c[M],l,r=1e6,mid,sum;
inline bool check(){
    for(int i=1;i<=n;++i)
        c[i]=a[i]-b[i]*mid;
    sort(c+1,c+1+n),sum=0;
    for(int i=k;i<=n;++i)
        sum+=c[i];
    return sum>=0;
}
int main(){
    while(1){
        n=read(),k=read();
        if(!n&&!k) return 0;
        l=0,r=1e6,++k;
        for(int i=1;i<=n;++i)
            a[i]=read();
        for(int i=1;i<=n;++i)
            b[i]=read();
        while(r-l>eps){
            mid=(l+r)/2;
            check()?l=mid:r=mid;
        }
        printf("%.0lf\n",100.0*l);
    }
    return 0;
}

1.最優比率揹包

題目大意和裸題的差很少，也是求 n 個物品中，\(\dfrac{\sum_{i=1}^{k} a[i]}{\sum_{j=1}^{k} b[j]}\) 的最大值，可是條件有變化，此次咱們要對權值 b 進行限制，首先咱們仍是令 \(sum=\sum_{i=1}^{k} a[i] ,\ tot=\sum_{j=1}^{k} b[j]\) ，那麼限制條件就是 \(tot>=W\)（\(W\)已給出）

那麼這個時候咱們設完 \(ans\) 以後獲得的式子一樣是：
\[\sum_{i=1}^{k} (a[i]-b[i]\times ans) <=0\]

因而 \(ans\) 照常二分，對於條件成立的判斷方式變一變就行了。

其實你應該也猜到了判斷方式，原先咱們不是貪心取前面的 k 個數麼，那此次咱們把貪心改爲揹包就行了啊

那樣的話，對於每一個物品來講，體積就是 b ，而價值就是上面的 \(a[i]-b[i]\times ans\)，而後咱們跑 0/1 揹包，判斷滿揹包價值的正負性。

題目

Talent Show

代碼

//by Judge
#include<cstdio>
#include<cstring>
#include<iostream>
#define ll long long
using namespace std;
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
inline void cmax(ll& a,ll b){if(a<b)a=b;}
inline int Min(int a,int b){return a<b?a:b;}
char buf[1<<21],*p1=buf,*p2=buf;
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,l,r=1e6,mid,V,t[305]; ll w[305],f[10005];
inline bool check(int tim){
    memset(f,-0x3f,sizeof(f)),f[0]=0;
    ll tmp=f[V];
    for(int i=1,j,k;i<=n;++i)
        for(j=V;~j;--j) if(f[j]!=tmp){
            cmax(f[Min(j+w[i],V)],f[j]+t[i]-w[i]*tim);
        }
    return f[V]>=0;
}
inline int calc(){
    while(l<=r){
        mid=l+r>>1;
        check(mid)?l=mid+1:r=mid-1;
    }
    return r;
}
int main(){
    n=read(),V=read();
    for(int i=1;i<=n;++i){
        w[i]=read(),
        t[i]=read()*1000;
    }
    return !printf("%d\n",calc());
}

2.最優比率生成樹

對於構造最優比率生成樹的分數規劃，其實差很少相似是將 \(n\) 物品轉換爲 \(m\) 條邊，而後求解，只不過限制改了，選出的 \(n-1\) 條邊要構成一棵樹

具體點說，題目通常會給你一張圖（有時候是徹底圖），而後給你 n 個點 m 條邊，求出原圖包含 n 個點的一棵樹，使得所選的邊的兩個權值和比值最大

那麼求法呢，其實也很相似：

咱們考慮構造最小生成樹的時候，邊是從小到大加的

那麼對於每一條邊，咱們讓他的邊權爲 \(a[i]-b[i]\times ans\)

而後和克魯斯卡爾同樣排個序一條一條嘗試加就行了（然鵝有時候你得嘗試普里姆算法）

若是加成功了就累計入 \(sum\)

最後根據 \(sum\) 的正負性來判斷 \(check\) 的返回值

其實和前面兩個沒高明到哪裏去，基本想法同樣的

題目

Desert King

代碼

（慢着，kruskal會炸？仍是我打開方式不對？只好開 \(Prim\) ）

//by Judge
#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#define eps 1e-5
using namespace std;
const int M=1005;
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
char buf[1<<21],*p1=buf,*p2=buf;
inline bool cmin(double& a,double b){return a>b?a=b,1:0;}
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,pat,x[M],y[M],z[M],vis[M];
double f[M][M],g[M][M],minv[M];
double l,r,mid,sum,tmp;
inline bool check(){
    memset(vis,0,sizeof(vis));
    sum=0,vis[1]=1;
    for(int i=1;i<=n;++i)
        minv[i]=g[1][i]-f[1][i]*mid;
    for(int i=2,j,k;i<=n;++i){
        for(tmp=1e16,k=-1,j=2;j<=n;++j)
            if(!vis[j]&&cmin(tmp,minv[j])) k=j;
        if(k==-1) break;
        for(vis[k]=1,sum+=tmp,j=2;j<=n;++j) if(!vis[j])
            cmin(minv[j],g[k][j]-f[k][j]*mid);
    }
    return sum>=0;
}
inline double dis(int i,int j){
    static double dx,dy;
    dx=x[i]-x[j],dy=y[i]-y[j];
    return sqrt(dx*dx+dy*dy);
}
int main(){
    while(1){
        n=read();
        if(!n) return 0;
        for(int i=1;i<=n;++i){
            x[i]=read(),
            y[i]=read(),
            z[i]=read();
        }
        for(int i=1;i<=n;++i)
            for(int j=i+1;j<=n;++j)
                f[j][i]=f[i][j]=dis(i,j),
                g[j][i]=g[i][j]=abs(z[i]-z[j]);
        l=0,r=100;
        while(r-l>eps){
            mid=(l+r)/2,
            (check()?l:r)=mid;
        }
        printf("%.3lf\n",l);
    }
    return 0;
}

3.最優比率環

簡而言之，最優比率環就是給你一個圖，圖上每條邊有兩個權值（固然，第二個權值可能恆爲 1 ），而後讓你在圖中找到一個環，令環上邊的個兩個權值和的比值最大（或是最小）

而後就是要用上面生成樹相似的思路，獲得每條邊邊權爲： \(a[i]-b[i]\times ans\) ，而後在圖上找負環（有時是正環，可是正負環能夠人工轉化的嘛），找到就 \(check\) 成功

題目

[HNOI2009]最小圈

代碼

//by Judge
#include<cstdio>
#include<cstring>
#include<iostream>
#define eps 1e-10
using namespace std;
const int M=1e5+7;
inline bool cmin(double& a,double b){return a>b?a=b,1:0;}
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,m,tim,pat,head[M],vis[M]; double l,r,mid,d[M];
struct Edge{ int to,next; double val; Edge(){}
    Edge(int a,double b,int c){ to=a,val=b,next=c; }
}e[M<<1];
inline void add(int u,int v,double z){
    e[++pat]=Edge(v,z,head[u]),head[u]=pat;
}
bool dfs(int u){ vis[u]=1;
    for(int i=head[u];i;i=e[i].next)
        if(cmin(d[e[i].to],d[u]+e[i].val-mid))
            if(vis[e[i].to]||dfs(e[i].to))
                return vis[u]=0,1;
    return vis[u]=0;
}
inline bool check(){
    memset(d,0,sizeof(d));
    for(int i=1;i<=n;++i)
        if(dfs(i)) return 1;
    return 0;
}
int main(){ n=read(),m=read(); double z;
    for(int i=1,x,y;i<=m;++i){
        x=read(),y=read();
        scanf("%lf",&z);
        add(x,y,z);
        if(z<0) l+=z; else r+=z;
    }
    while(r-l>eps){
        mid=(l+r)/2,
        (check()?r:l)=mid;
    }
    return !printf("%.8lf\n",l);
}

4.最大密度子圖

題目通常就是給出 \(n\) 個點而後讓你導出一個子圖，讓這個子圖中邊數與點數的比值最大

不過有時候題目中的邊（甚至是點）可能會帶上權值，那麼具體問題具體分析就行了

題目

Hard Life

（話說這道題翻譯照搬的\(bzoj\)，\(UVA\) 是要輸出方案的啊...）

代碼

（好像就這個代碼最長了，題目講得簡單然鵝碼量巨大）

//by Judge
#include<cmath>
#include<queue>
#include<cstdio>
#include<cstring>
#include<iostream>
#define db double
using namespace std;
const db eps=1e-8;const int N=1005;
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,m,S,T,pat,ans,h[N],q[N],head[N],cur[N];
int u[N],v[N],du[N],vis[N]; db U,l,r,mid;
struct Edge { int to,nxt; db flow; } e[N<<5];
void insert(int u,int v,db w) {
    e[++pat]=(Edge){v,head[u],w},head[u]=pat;
    e[++pat]=(Edge){u,head[v],0},head[v]=pat;
}
#define v e[i].to
inline bool bfs() { int hd=0,tl=1,u;
    memset(h,0,sizeof(h)),h[q[1]=S]=1;
    while(hd<tl) { u=q[++hd];
        for(int i=head[u]; i; i=e[i].nxt)
            if(e[i].flow>eps&&!h[v])
                h[v]=h[u]+1,q[++tl]=v;
    } return h[T];
}
db dfs(int x,db F) {
    if(x==T) return F; db w,used=0;
    for(int &i=cur[x]; i; i=e[i].nxt)
        if(h[x]+1==h[v]) {
            w=dfs(v,min(e[i].flow,F-used));
            e[i].flow-=w,e[i^1].flow+=w;
            used+=w; if(used==F) return F;
        }
    if(!used) h[x]=0;
    return used;
}
db dinic() { db res=0;
    for(;bfs();res+=dfs(S,1e16))
        memcpy(cur,head,sizeof(cur));
    return res;
}
void DFS(int x) { ++ans,vis[x]=1;
    for(int i=head[x]; i; i=e[i].nxt)
        if(e[i].flow>eps&&!vis[v]) DFS(v);
}
#undef v
inline bool check() {
    pat=1,S=0,T=n+1;
    memset(head,0,sizeof(head));
    for(int i=1; i<=n; ++i){
        insert(S,i,U),
        insert(i,T,U+mid+mid-du[i]);
    }
    for(int i=1; i<=m; ++i){
        insert(u[i],v[i],1),
        insert(v[i],u[i],1);
    }
    return U*n-dinic()>eps;
}
int main() {
    while(scanf("%d%d",&n,&m)!=EOF) {
        for(int i=0; i<=n; ++i) du[i]=0;
        memset(vis,0,sizeof(vis)),ans=-1,U=m;
        if(!m) {puts("1\n1\n\n");continue;}
        for(int i=1; i<=m; ++i)
            u[i]=read(),v[i]=read(),
            ++du[u[i]],++du[v[i]];
        for(l=0,r=m;r-l>1.0/n/n;){
            mid=(l+r)/2,
            (check()?l:r)=mid;
        }
        mid=l,check(),DFS(S);
        printf("%d\n",ans);
        for(int i=1; i<=n; ++i)
            if(vis[i]) printf("%d\n",i);
    } return 0;
}

5.最優比率流

仍是那句老話，和以前同樣的，就是根據二分出來的答案賦邊權，而後通常是每次建一邊圖跑費用流根據最小費用的正負性判斷 \(ans\) 合法性（用 \(zkw\) 跑費用流有時候會掛，好比這個例題，固然多是我天生大常數那就 \(GG\) ）

題目

[SDOI2017]新生舞會

（這裏是一道二分圖最大費用最大流，可是代碼實現中能夠邊權取反而後跑最小費用）

代碼

//by Judge
#include<queue>
#include<cstdio>
#include<cstring>
#include<iostream>
#define eps 1e-8
#define ll long long
using namespace std;
const double inf=1e18+7;
const int M=205;
#ifndef Judge
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#endif
inline int Min(int a,int b){return a<b?a:b;}
inline bool cmax(double& a,double b){return a<b?a=b,1:0;}
char buf[1<<21],*p1=buf,*p2=buf;
inline int read(){ int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0'; return x*f;
} int n,S,T,pat=1,tot,head[M],cur[M]; ll a[M][M],b[M][M];
int inq[M],vis[M]; double l,r=1e4,mid,ans,dis[M]; queue<int> q;
struct Edge{ int to,flow,next; double val;
    Edge(int a,double b,int c,int d){
        to=a,val=b,flow=c,next=d;
    } Edge(){}
}e[M*M*M];
inline void add(int u,int v,double c,int w){
    e[++pat]=Edge(v,c,w,head[u]),head[u]=pat;
    e[++pat]=Edge(u,-c,0,head[v]),head[v]=pat;
}
#define v e[i].to
inline bool spfa(){
    for(int i=S;i<=T;++i)
        dis[i]=-inf,vis[i]=0;
    dis[S]=0,q.push(S);
    while(!q.empty()){
        int u=q.front();
        q.pop(),vis[u]=0;
        for(int i=head[u];i;i=e[i].next)
            if(e[i].flow&&cmax(dis[v],dis[u]+e[i].val))
                if(!vis[v]) vis[v]=1,q.push(v);
    } return dis[T]!=-inf;
}
int dfs(int u,int flow){
    if(u==T) return ans+=dis[T]*flow,flow;
    int res=0; vis[u]=1;
    for(int &i=cur[u],k;i;i=e[i].next)
        if(dis[v]==dis[u]+e[i].val)
            if(!vis[v]&&e[i].flow){
                k=dfs(v,Min(e[i].flow,flow-res));
                if(k){
                    res+=k;
                    e[i].flow-=k;
                    e[i^1].flow+=k;
                    if(res==flow)
                        break;
                }
            }
    return res;
}
inline bool check(){
    pat=1,ans=0,
    memset(head,0,sizeof(head));
    for(int i=1;i<=n;++i){
        add(S,i,0,1),
        add(n+i,T,0,1);
    }
    for(int i=1;i<=n;++i) 
        for(int j=1;j<=n;++j)
            add(i,n+j,a[i][j]-b[i][j]*mid,1);
    for(;spfa();dfs(S,1e9))
        memcpy(cur,head,sizeof(cur));
    return ans<=0;
}
int main(){
    n=read(),T=n+1<<1;
    for(int i=1;i<=n;++i)
        for(int j=1;j<=n;++j)
            a[i][j]=read();
    for(int i=1;i<=n;++i)
        for(int j=1;j<=n;++j)
            b[i][j]=read();
    for(tot=pat;r-l>eps;){
        mid=(l+r)/2,
        (check()?r:l)=mid;
    }
    return !printf("%.6lf\n",l);
}

6.各類其餘結合

洛谷 blog 已添加

分數規劃的一點感想

講道理，其實分數規劃這個東西哪兒都能套，然鵝如今關於分數規劃的題目卻並非不少，至少沒有到氾濫的地步，不過相信在不（知道多）久的未來，分數規劃會成爲人盡皆知的一種算法

\[tHe\ EnD\]

哦，對了，參考文獻？

01分數規劃問題相關算法與題目講解

poj 3155 (最大密度子圖)

另外的另外，有時間的話會把本身出的分數規劃的題目放到6這塊內容上來

相關文章
相關標籤/搜索
每日一句
    每一个你不满意的现在，都有一个你没有努力的曾经。
本站公眾號
   歡迎關注本站公眾號,獲取更多信息
聯繫我們 最近搜索 最新文章 沪ICP备13005482号-10 MyBatis教程 SQL 教程 MySQL教程 Java 教程 Thymeleaf 教程 Hibernate教程 Spring教程 Redis教程