01分数规划
Z_Healer
·
·
个人记录
前言
本文基础知识主要来源于OI Wiki ,所有题目来源于luogu
分数规划
分数规划用来求一个分式的极值。
形象一点就是,给出a_i和b_i,求一组w_i \in \left \{0,1 \right \} ,最大/小化
\dfrac{ \sum_{i=1}^n a_i*w_i}{\sum_{i=1}^n b_i*w_i}
另外一种描述:每种物品有两个权值a和b,选出若干个物品使得\dfrac{\sum{}{}a}{\sum{}{}b}最大/小。
二分法
分数规划问题的通用方法是二分。
假设我们求最大值,二分一个答案mid,当存在一组w_i,满足\dfrac{\sum{}{}a_i*w_i}{\sum{}{}b_i*w_i}>mid时,则说明当前的mid可以。对原不等式进行推导:
\frac{\sum{}{}a_i*w_i}{\sum{}{}b_i*w_i}>mid
\Rightarrow \sum{}{}a_i*w_i-mid*\sum{}{}b_i*w_i>0
\Rightarrow \sum{}{}w_i*(a_i-mid*b_i)>0
那么只要求出不等号左边的式子的最大值就行了。如果最大值比0要大,说明mid是可行的,否则不可行。求最小值只要比0要小就行了。
注意:二分答案mid一般都要比题目要求的精度要多两位。
Dinkelbach算法
$L$来输入,不断地迭代,直至答案满足精度。
分数规划的主要难点就在于**如何求**$\sum{}{}w_i*(a_i-mid*b_i)$**的最大值/最小值**。下面通过一些例题来讲解该式子的最大/小值的求法,对于其最大/小值的表达式的推导,还请读者自行完成,代码也着重展示$check$部分。
## [P1570 KC喝咖啡](https://www.luogu.com.cn/problem/P1570)
>有$n$个物品,每个物品有两个权值$a$和$b$,选$m$个物品使得$\frac{\sum{}{}a_i}{\sum{}{}b_i}$最大。
模板题,二分一个答案$mid$,把$a_i-mid*b_i$作为第$i$个物体的权值,排序后贪心地选出前$m$大的,如果大于$0$就说明当前$mid$可以。
```cpp
#include<bits/stdc++.h>
#define N 205
using namespace std;
int n,m;
double a[N],b[N],c[N];
inline int read(){
int w=0;char ch=getchar();
while(ch<'0'||ch>'9') ch=getchar();
while(ch>='0'&&ch<='9'){
w=(w<<3)+(w<<1)+(ch^48);
ch=getchar();
}
return w;
}
inline bool check(double mid){
for(int i=1;i<=n;i++) c[i]=a[i]-b[i]*mid;
sort(c+1,c+n+1);
double res=0;
for(int i=1;i<=m;i++) res+=c[n-i+1];
return res>=0;
}
inline double work(){
double l=0,r=1000,eps=1e-5,mid;
while(r-l>=eps){
mid=(l+r)/2.0;
if(check(mid)) l=mid;
else r=mid;
}
return l;
}
int main(){
n=read();m=read();
for(int i=1;i<=n;i++) a[i]=read();
for(int i=1;i<=n;i++) b[i]=read();
printf("%.3lf",work());
return 0;
}
```
## [P4377 [USACO18OPEN]Talent Show G](https://www.luogu.com.cn/problem/P4377)
>有$n$个物品,每个物品有两个权值$a$和$b$,选出一些物品使得$\frac{\sum{}{}a_i}{\sum{}{}b_i}$最大,要求$\sum b_i\ge W$。
加上分子比$W$大的限制后,就不能像上一道题一样直接贪心求最大值了,先将条件像上一道题转换,则题目转化为有一个容量为$W$的背包,有$n$个价值为$a_i-mid*b_i$,体积为$b_i$的物品,求最大价值,就变成了一个经典的[$01$背包问题](https://www.luogu.com.cn/problem/P1048)。
```cpp
inline bool check(double x){
for(int i=1;i<=W;i++) f[i]=-1e9;//double类型不能memset
for(int i=1;i<=n;i++){//滚动数组优化空间
for(int j=W;j>=0;j--){
int k=min(W,(int)(j+b[i]));
f[k]=max(f[k],f[j]+a[i]-x*b[i]);
}
}
return f[W]>0;
}
```
## 最优比率生成树
>图$G$中每条边有两个权值$a$和$b$,求一棵生成树$T$使得 $\frac{ \sum_{e\in T}^{} a_e}{\sum_{e\in T}{}b_e}$最小。
把$a_i-mid*b_i$作为每条边的权值,那么[最小生成树](https://www.luogu.com.cn/problem/P3366)就是最小值,最小值小于$0$就说明$mid$可以。
```cpp
struct QAQ{
int to,next;double a,b;
}e[M<<1];
int h[N],vis[N],tot;
double dis[N];
inline void add(int x,int y,double za,double zb){
e[++tot]=(QAQ){y,h[x],za,zb};h[x]=tot;
e[++tot]=(QAQ){x,h[y],za,zb};h[y]=tot;
}
inline bool check(double mid){
for(int i=1;i<=n;i++) dis[i]=1e9,vis[i]=0;
dis[1]=0;double res=0;
for(int i=1;i<=n;i++){
int k=0;
for(int j=1;j<=n;j++)
if(!vis[j]&&dis[k]>dis[j]) k=j;
vis[k]=1;res+=dis[k];
for(int j=h[k];j;j=e[j].next)
dis[e[j].to]=min(dis[e[j].to],e[j].a-mid*e[i].b);
}
return res<0;
}
```
## [P3199 [HNOI2009]最小圈](https://www.luogu.com.cn/problem/P3199)
>每条边的边权为$w$,求一个环$C$使得$\frac{\sum_{e\in C}{}w_e}{|C|}$最小。
把$a_i-mid$作为边权,那么权值最小的环就是最小值。因为我们只需要判最小值是否小于$0$,所以只需要判断图中是否存在[负环](https://www.luogu.com.cn/problem/P3385)(要用$DFS$判负环,否则会$TLE$)即可。另外本题存在一种复杂度$O(nm)$的算法,如果有兴趣可以阅读[这篇文章](https://www.cnblogs.com/y-clever/p/7043553.html)
```cpp
bool spfa(int x,LD mid){
vis[x]=1;
for(int i=h[x];i;i=e[i].next){
int y=e[i].to;
if(dis[y]>dis[x]+e[i].w-mid){
dis[y]=dis[x]+e[i].w-mid;
if(vis[y]||spfa(y,mid)) return 1;
}
}
vis[x]=0;
return 0;
}
inline bool check(LD mid){
for(int i=1;i<=n;i++) vis[i]=dis[i]=0;
for(int i=1;i<=n;i++) if(spfa(i,mid)) return 1;
return 0;
}
```
类似的还有:
[P1768 天路](https://www.luogu.com.cn/problem/P1768)
[P3288 [SCOI2014]方伯伯运椰子](https://www.luogu.com.cn/problem/P3288)
[P2868 [USACO07DEC]Sightseeing Cows G](https://www.luogu.com.cn/problem/P2868)
## [P1730 最小密度路径](https://www.luogu.com.cn/problem/P1730)
>有一个$N$个点$M$条边的有向无环图$G$,$Q$次询问,每次询问$x$到$y$的路径上边权和除以边的数量的最小值。
把$w_{u,v}-mid$作为边权,以$x$为起点跑[$SPFA$](https://www.luogu.com.cn/problem/P3371)(因为有负权边),如果$dis[y]\le 0$则说明$mid$可以,要先预处理好所有询问,否则会超时。
```cpp
inline bool check(int s,int t,double mid){
queue <int> q;
for(int i=0;i<=n;i++) dis[i]=1e9,vis[i]=0;
q.push(s);dis[s]=0;
while(!q.empty()){
int x=q.front();q.pop();
vis[x]=0;
for(int i=h[x];i;i=e[i].next){
int y=e[i].to;
if(dis[y]>dis[x]+e[i].w-mid){
dis[y]=dis[x]+e[i].w-mid;
if(!vis[y]){
vis[y]=1;
q.push(y);
}
}
}
}
return dis[t]<=0;
}
inline double work(int x,int y){
check(x,y,0);
if(dis[y]==1e9) return -1;//无解
LD l=0,r=5e6,eps=1e-5;
while(r-l>eps){
LD mid=(l+r)/2.0;
if(check(x,y,mid)) r=mid;
else l=mid;
}
return r;
}
```
## [P3705 [SDOI2017]新生舞会](https://www.luogu.com.cn/problem/P3705)
>有$n$个男生和$n$对女生,男生和女生可以任意组成舞伴,第$i$位男生和第$j$为女生组成舞伴有两个值$a_{i,j}$和$b_{i,j}$,求$\frac{\sum{}{}a}{\sum{}{}b}$的最大值
把$a_{i,j}-mid*b_{i,j}$作为第$i$位男生和第$j$为女生组成舞伴的权值,要求男女两两匹配,求最大权值。这个问题就转换为网络流,男女间连流量为$1$,费用为$a_{i,j}-mid*b_{i,j}$的边,起点向男生,女生向终点连流量为$1$,费用为$0$的边,跑[最大费用](https://www.luogu.com.cn/problem/P3381)。
```cpp
inline int spfa(){
for(int i=0;i<=t;i++) dis[i]=-1e9;
memset(pre,0,sizeof(pre));
memset(vis,0,sizeof(vis));
dis[s]=0;
queue <int> q;
q.push(s);
while(!q.empty()){
int x=q.front();q.pop();vis[x]=0;
for(int i=h[x];i;i=e[i].next){
int y=e[i].to;
if(e[i].w>0&&dis[y]<dis[x]+e[i].v){
dis[y]=dis[x]+e[i].v;
pre[y][0]=x;
pre[y][1]=i;
if(!vis[y]){
q.push(y);
vis[y]=1;
}
}
}
}
return dis[t]!=-1e9;
}
inline int EK(){
int ans=0;cost=0;
while(spfa()){
int mi=1e9;
for(int i=t;i!=s;i=pre[i][0])
mi=min(mi,e[pre[i][1]].w);
for(int i=t;i!=s;i=pre[i][0]){
e[pre[i][1]].w-=mi;
e[pre[i][1]^1].w+=mi;
}
ans+=mi;
cost+=(double)mi*dis[t];
}
return ans;
}
inline bool check(double x){
memset(h,0,sizeof(h));tot=1;
for(int i=1;i<=n;i++) add(s,i,1,0),add(i+n,t,1,0);
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
add(i,j+n,1,a[i][j]-x*b[i][j]);
EK();
return cost>=0;
}
```
## [P4322 [JSOI2016]最佳团体](https://www.luogu.com.cn/problem/P4322)
>在一颗树内选$K$个节点(不包括$0$节点),要求选该节点其父亲节点一定要被选(一开始$0$节点被选),每个节点有两个值$a_i$和$b_i$,求$\frac{\sum_{i=1}^{k}a_i}{\sum_{i=1}^{k}b_i}$最大值。
把$a_{i,j}-mid*b_{i,j}$作为第i个节点的权值,那么题目就变为了[树上背包](https://www.luogu.com.cn/problem/P2014)
```cpp
double f[N][N],a[N],b[N];//i的子树中选j个
vector<int> e[N];
void dfs(int x,int t,double c){
if(t<=0) return;
for(int i=0;i<(int)e[x].size();i++){
int y=e[x][i];
for(int j=0;j<t;j++) f[y][j]=f[x][j]+a[y]-b[y]*c;
dfs(y,t-1,c);
for(int j=1;j<=t;j++) f[x][j]=max(f[x][j],f[y][j-1]);
}
}
inline bool check(double c){
for(int i=0;i<=n;i++)
for(int j=1;j<=k;j++) f[i][j]=-1e9;
dfs(0,k,c);
return f[0][k]>=0;
}
```
## [P6087 [JSOI2015]送礼物](https://www.luogu.com.cn/problem/P6087)
>给你一个序列,求区间长度在$[L,R]$之间的美观度最大,区间$[l,r]$的美观度为$\dfrac{\max_{i=l}^{r}a_i-\min_{i=l}^{r}a_i}{r-l+k}$。
当看到求一个**分式的最大值**时,第一时间想到**分数规划**,但这道题与普通的分数规划不一样,不能直接定义权值,因为区间最大/小值(分子)与区间边界(分母)**没有直接关系**。
对于最大/小值在区间内部的,缩小区间长度只会使答案更优。那么我们尝试将最大/小值放到区间的边界上,则美观度可以变为$\dfrac{\left|a[r]-a[l]\right| }{r-l+k}$,但不过这样可能满足不了$r-l+1\ge L$的限制。我们可以先将小于$L$处理掉,再处理大于的。
对于一组最大最小值距离小于$L$的,我们将长度直接扩展到$L$是最优的,那么分母就确定为$L-1+k$,问题就变为求$L$个连续的数中的最大值和最小值,是[单调队列的模板题](https://www.luogu.com.cn/problem/P1886)。
对于长度大于$L$的区间,我们通过$01$分数规划,二分一个答案$mid$对式子$\dfrac{\left|a[r]-a[l]\right| }{r-l+k}\ge mid$进行推导。
当$a[r]>a[l]$,$(a[r]-r*mid)-(a[l]-l*mid)\ge k*mid
当a[r]<a[l],(a[l]+l*mid)-(a[r]+r*mid)\ge k*mid
由区间长度小于L的方法,我们可以想到用单调队列来处理。
当a[r]>a[l],l从1开始,r从L开始,维护(a[l]-l*mid)单调递减,当队头距离i大于R时弹出,每次从队尾加入a[i-L+1],当(a[i]-i*mid)-(a[q[h]]-q[h]*mid)\ge k*mid成立时,当前mid可以。
当a[r]<a[l],维护(a[l]-l*mid)单调递增,当(a[q[h]]-q[h]*mid)-(a[i]-q[i]*mid)\ge k*mid成立时,当前mid可以。
inline LD pre(){//处理区间长度恰好为L的
LD res=0;
int he1=1,he2=1,ta1=0,ta2=0,q1[N],q2[N];
for(int i=1;i<L;i++){
while(he1<=ta1&&a[i]>=a[q1[ta1]]) ta1--;//维护递增
while(he2<=ta2&&a[i]<=a[q2[ta2]]) ta2--;//维护递减
q1[++ta1]=i;q2[++ta2]=i;
}
for(int i=L;i<=n;i++){
while(he1<=ta1&&a[i]>=a[q1[ta1]]) ta1--;
while(he2<=ta2&&a[i]<=a[q2[ta2]]) ta2--;
q1[++ta1]=i;q2[++ta2]=i;
while(he1<=ta1&&i-q1[he1]>=L) he1++;
while(he2<=ta2&&i-q2[he2]>=L) he2++;
res=max(res,(a[q1[he1]]-a[q2[he2]])/(LD)(L-1+k));
}
return res;
}
inline bool check(LD x){
LD c[N];
int he=1,ta=0,q[N];
memset(q,0,sizeof(q));
for(int i=1;i<=n;i++) c[i]=a[i]-x*(LD)i;
for(int i=L;i<=n;i++){//维护区间中的最小值最小
while(he<=ta&&i-q[he]>=R) he++;//保证区间长度小于R
while(he<=ta&&c[i-L+1]<=c[q[ta]]) ta--;
q[++ta]=i-L+1;//保证区间长度大于L
if(c[i]-c[q[he]]>=x*k) return 1;
}
he=1;ta=0;
memset(q,0,sizeof(q));
for(int i=1;i<=n;i++) c[i]=a[i]+x*(LD)i;
for(int i=L;i<=n;i++){//维护区间中的最大值最大
while(he<=ta&&i-q[he]>=R) he++;
while(he<=ta&&c[i-L+1]>=c[q[ta]]) ta--;
q[++ta]=i-L+1;
if(c[q[he]]-c[i]>=x*k) return 1;
}
return 0;
}
inline LD work(){
LD l=0,r=1e3,mid,eps=1e-6;
while(r-l>eps){
mid=(l+r)/2.0;
if(check(mid)) l=mid;
else r=mid;
}
return max(l,pre());
}
UVA1389 Hard Life
定义一个无向图G<V,E>的密度D=\dfrac{|E|}{|V|},给你一个无向图,求出它密度最大的子图。
最大密度子图模板题。
看到求一个分式的最大值,先进行分数规划,得到:|E|-mid*|V|\ge 0。要求左边式子最大值,我们由之前做过的题出发,都是考虑\sum内单个的性质,所以我们对一条边与其连接 的两个点出发。发现选一条边就要将它连接的两个点选上,想到最大权闭合子图,把边看成建模中点,将原图中边连的两个点当做后继,则选一个边的权值为1,两个点的权值为-mid。根据最大权闭合子图建图:(s,v,1), (v_e,x,+∞),(v_e,y,+∞),(x,t,-mid),(y,t,-mid)。
在二分mid时精度过小,可能会跳过某些解,所以我们要分析两个子图的密度差最小为多少。二分答案mid,下界为\dfrac{1}{n},上界为m。设D_1=\dfrac{m_1}{n_1},D_2=\dfrac{m_2}{n_2},则有
\dfrac{m_1}{n_1}-\dfrac{m_2}{n_2}=\dfrac{m_1n_2-m_2n_1}{n_1n_2}\geq \dfrac{1}{n_1n_2}\geq \dfrac{1}{n^2}
我们可以得到一个性质:任意两个不同密度的子图 G_1,G_2的密度差 D_1-D_2 \geq \dfrac{1}{n^2}。
若对于无向图G中有一个密度为D的子图G',且不存在一个密度超过D+\dfrac{1}{n^2}的子图,则G'为最大密度子图。二分答案时只需要将精度eps调成\dfrac{1}{100n^2}(多两位保证正确性)。
最后输出方案时,先check(ans),则网络中流经的点就是要被选中的。
double eps=1e-6;
struct QAQ{
int to,next;double w;
}e[M<<1];
inline bool bfs(){
queue <int> q;
for(int i=0;i<=t;i++) dep[i]=0,cur[i]=h[i];
q.push(s);
dep[s]=1;
while(!q.empty()){
int x=q.front();q.pop();
for(int i=h[x];i;i=e[i].next){
int y=e[i].to;
if(!dep[y]&&e[i].w>eps){
dep[y]=dep[x]+1;
if(y==t) return 1;
q.push(y);
}
}
}
return 0;
}
double dfs(int x,double mi){
if(x==t) return mi;
double used=0;
for(int i=cur[x];i;i=e[i].next){
cur[x]=i;
int y=e[i].to;
if(dep[y]==dep[x]+1&&e[i].w>eps){
double tmp=dfs(y,min(e[i].w,mi-used));
used+=tmp;
e[i].w-=tmp;
e[i^1].w+=tmp;
if(used==mi) break;
}
}
return used;
}
inline bool check(double mid){
memset(h,0,sizeof(h));tot=1;
for(int i=1;i<=m;i++)
add(i,px[i]+m,INF),add(i,py[i]+m,INF);
for(int i=1;i<=n;i++) add(i+m,t,mid);
for(int i=1;i<=m;i++) add(s,i,1);
double res=m;
while(bfs())res-=dfs(s,INF);
return res>eps;
}
inline double work(){
double l=0,r=1000,mid;
while(r-l>=eps){
mid=(l+r)/2.0;
if(check(mid)) l=mid;
else r=mid;
}
return l;
}
void get(int x){
vis[x]=1;
for(int i=h[x];i;i=e[i].next)
if(e[i].w&&!vis[e[i].to]) get(e[i].to);
}
inline void out(double x){
int res=0;
memset(vis,0,sizeof(vis));
check(x);get(s);
for(int i=1;i<=n;i++) if(vis[i+m]) res++;
printf("%d\n",res);
for(int i=1;i<=n;i++) if(vis[i+m]) printf("%d\n",i);
}
总结