概率dp

首先是正着递推的计算概率的dp问题

https://ac.nowcoder.com/acm/contest/28263/A

纯数学题

对随机的数字大小分类讨论,计算概率的时候利用高中几何概型的线性规划手法进行计算。

double g=0.5;
void solve(){
	double k,x;cin>>k>>x;
	double ans=0;
	if(k==x)ans=g;
	else if(k>=2*x)ans=0;
	else if(k<x){
		ans=g*(x*x-k*k);
		ans/=x*x;
		ans+=g;
	}
	else {
		ans=g*(2*x-k)*(2*x-k);
		ans/=x*x;;
	}
	baoliu(ans,2);
	cout<<endl;
}

CF148D Bag of mice

img

状态:$f[i][j]$表示袋中有i只白鼠j只黑鼠时,A获胜的概率

起点:$f[0][i]=0,f[i][0]=1$
终点:$f[w][b]$
转移:
1.先手拿到白鼠:$f[i][j]+=i/(i+j)$
2.先手黑鼠,后手白鼠:$f[i][j]+=0,$ 这种情况不用处理
3.先手黑鼠,后手黑鼠,跑白鼠:$f[i][j]+=j/(i+j) *(j-1)/(i+j-1) *i/(i+j-2) *f[i-1][j-2]$
4.先手黑鼠,后手黑鼠,跑黑鼠:$f[i][j]+=j/(i+j) *(j-1)/(i+j-1) *(j-2)/(i+j-2) *f[i][j-3]$

int n, m;
int a[N];
//两个人各操作一次为1轮
double f[N][N];
void solve(){
	cin>>n>>m;
	for(int i=1;i<=n;i++){
		f[i][0]=1;
		
	}
	for(int i=1;i<=m;i++){
		f[0][i]=0;
	}
	for(int i=1;i<=n;i+=1.0){
		for(int j=1;j<=m;j+=1.0){
			f[i][j]+=(double)i/(i+j);
  if(i>=1 && j>=2)
  f[i][j]+=f[i-1][j-2]*(double)j*(j-1)*(i)/(double)((i+j)*(i+j-1)*(i+j-2));
  if(j>=3)
  f[i][j]+=f[i][j-3]*(double)j*(j-1)*(j-2)/(double)((i+j)*(i+j-1)*(i+j-2));
  //cerr<<f[i][j]<<" ";
		}
	}
	baoliu(f[n][m],9);
}

POJ3071 Football

题意:有2^n支球队比赛,每次和相邻的球队踢,两两淘汰,给定任意两支球队相互踢赢的概率,求最后哪只球队最可能夺冠

img

img

Solution:考虑$f[i][j]$是第i次大混战编号为j的队伍获胜的概率,考虑转移需要找到本轮的对手,并调用其获胜的概率。由于n范围很小,直接暴力枚举即可。核心问题是如何找到合法可能的对手,我们观察后利用位运算得到。

#include<iostream>
#include<cstring>
#include<algorithm>
using namespace std;
#define baoliu(x,y) cout<<fixed<<setprecision(y)<<x;
const int N=130;
double w[N][N];
double f[8][N];
//f[i][j]:第i场,第j队win的概率
int main(){
	int n;
	//从0下标开始才能在后续进行高效位运算
while(cin>>n,n!=-1){
	//对于这种多测也考虑写成solve函数,忘清空了
	int m=1<<n;
	memset(f,0,sizeof f);
	for(int i=0;i<m;i++){
		for(int j=0;j<m;j++){
			cin>>w[i][j];
		}
	}
	for(int i=0;i<m;i++)f[0][i]=1;
	for(int i=1;i<=n;i++){
		for(int j=0;j<m;j++){
			for(int k=0;k<m;k++){
				if(((j>>(i-1))^1)==(k>>(i-1))){
					//所有情况的获胜概率是累加的
					f[i][j]+=f[i-1][j]*f[i-1][k]*w[j][k];
				}
			}
		}
	}
	double mx=0;int pos=0;
	for(int i=0;i<m;i++){
		if(mx<f[n][i]){
			mx=f[n][i];
			pos=i;
		}
	}
	//仔细读题发现,队伍编号从1开始,而我从0开始是为利用位运算
	cout<<pos+1<<endl;
}
}

[CF768D Jon and Orbs](https://www.luogu.com.cn/problem/CF768D

题面翻译

有k件物品,每天随机出现一件。设每种球都出现至少一次以上为事件U.

求最小天数使得U发生的概率不小于$\frac{p}{2000}$,有 $q$ 组询问,每组询问给定这个 $p$。

Solution:calendar:

考虑$dp[i][j]=dp[i-1][j]*p_{1}+dp[i-1][j-1]*P_{2} $

其中$ p_{1}为第i-1天已经有j种物品的概率,p_{1}=\frac{j}{i} ,$

$ p_{2}$为第i-1天已经有j-1种物品的概率,第一天从剩下的j-i+1种物品的随机到了,则$p_{2}=\frac{k-(j-1)}{i} $

double f[M][N];
void solve(){
	int k;
	cin>>k>>m;
	f[0][0]=1;
	for(int i=1;i<=20000;i++){
		for(int j=1;j<=k;j++){
		f[i][j]=f[i-1][j-1]*(double)(k-j+1)/k+f[i-1][j]*(double)j/k;
		}
	}
	int pos=0;
	while(pos<=M&&f[pos][k]<0.99){
		pos++;
		
	}
	cerr<<pos<<endl;
	while(m--){
		double p;cin>>p;
		int pos=0;
		//根据题目数据范围,概率设置到0.5,7274收敛
		//概率趋于1的时候,11503收敛
		while(f[pos][k]<(double)p/(double)2000){
			pos++;
			
		}
		cout<<pos<<endl;
	}
}

下面开始比较常见的期望dp

期望dp往往是最终态的答案比较好求,所以我们设置dp数组的状态含义是:到达终点的期望步数

恶意竞争

题意:一个球,n种材质,m种花纹,每次随机从无限多的球堆里拿一个球,问期望多少次可以拿到这所有花纹和材质的球,每次等概率随机拿。

Solution:发现计算概率只需要知道当前手里已经拿了多少种足球和篮球,所以我们容易定义dp状态为$dp[i][j]$,转移分别是四种情况也就是拿到的这个球的花纹是不是新的,这个球的材质是不是新的。

debug:双重循环下标错误导致RE 。dp数组下标写错导致wa

double f[N][N];
void solve(){
	cin>>n>>m;
	f[n][m]=0;
	for(int i=n;i>=0;i--){
		for(int j=m;j>=0;j--){
			if(i==n&&j==m)continue;
			//分母为0
			double p1=1.0*i*j/(double)(n*m);
			double p2=1.0*i*(m-j)/(double)(n*m);
			double p3=1.0*j*(n-i)/(double)(n*m);
			double p4=1.0*(n-i)*(m-j)/(double)(n*m);
			//cerr<<p1<<" "<<p2<<" "<<p3<<" "<<p4<<endl;
			f[i][j]+=p2*f[i][j+1]+p3*f[i+1][j]+p4*f[i+1][j+1]+1;
			f[i][j]/=(double)(1.0-p1);
		//	cerr<<i<<" "<<j<<" "<<f[i][j]<<endl;
			
					}
	}
	baoliu(f[0][0],6);
}

[Luogu P4492 HAOI2018]苹果树

[Luogu P2473 SCOI2008] 奖励关

[Luogu P2221 HAOI2012]高速公路

[Luogu P4065 JXOI2017]颜色

拓扑排序dp递推————记忆化搜索的期望dp

最朴素经典的图上等概率随机游走问题,先从dag入手,考虑一个有向无环图,从1走到n,1是该图的起点,n是该图的终点,每次从当前点等概率选择一个出边走。

先考虑每次走到下一个点的概率是$\frac{1}{out} $,out为该点的出度,我们直接从1开始dfs搜索,考虑到dag的性质,到n的终点只有一条,对于同一个点只需要搜索一次,我们采用记忆化。

#include<bits/stdc++.h>
//dag_>dfs
using namespace std;
const int N=1e5+10;
#define int long long
struct edge{int v;double w;};
vector<edge>e[N];
int n,m;
double  ans[N];
//ans[n]=0;
void dfs(int u,int fa){
	//cerr<<u<<endl;
	//长期不写记忆化搜索导致的
	if(ans[u])return ;
    int in=e[u].size();
    for(auto [v,w]:e[u]){
        if(v==fa)continue;
        dfs(v,u);
        ans[u]+=1.0/in*(ans[v]+w);
        cerr<<u<<" "<<v<<endl;
    }
    
}

signed  main(){
    cin>>n>>m;
    ans[n]=0;
    for(int i=1;i<=m;i++){
        int u,v;double w;cin>>u>>v>>w;
        e[u].push_back({v,w});
     //   e[v].push_back({u,w});
    }
    dfs(1,0);
    for(int i=1;i<=n;i++)cerr<<ans[i]<<" ";
    cout<<fixed<<setprecision(2)<<ans[1]<<endl;
    return 0;
}

下面给出这题的第二种写法:建立反图,拓扑排序后递推更新答案

debug:累加答案的时候要用的是正图上的出度去计算概率

struct edge{int v;double w;};
vector<edge>e[N];
int din[N];
vector<int>tp;
double ans[N];
bool topsort(){
queue<int>q;
for(int i=1;i<=n;i++){
	if(din[i]==0)q.push(i);
}
  while(q.size()){
  	auto u=q.front();q.pop();
  	tp.push_back(u);
	for(auto [v,w]:e[u]){
		din[v]--;
		
		if(din[v]==0){
			q.push(v);
			
		}
	}
  }
  int cnt=tp.size();
  return cnt==n;
}
int tmp[N];
void solve(){
	cin>>n>>m;
	for(int i=1;i<=m;i++){
		int u,v;double w;cin>>u>>v>>w;
		tmp[u]++;din[u]++;
		e[v].push_back({u,w});
	}
	topsort();
	for(auto u:tp){
		//cerr<<u<<" "<<din[u]<<endl;
		for(auto [v,w]:e[u]){
			//cerr<<v<<" "<<w<<endl;
			ans[v]+=1.0/tmp[v]*(ans[u]+w);
		}
	}
	baoliu(ans[1],2);
}

HDU3853 LOOPS

HDU4035 Maze

高斯消元解决后效性dp,引入就是在一般图上随机游走

img

Solution:上面的图片已经很清楚了,思路。首先考虑排序不等式的贪心,想要编号就要知道每个边的期望次数,从而又想要知道经过每个点的期望经过次数,再考虑dp发现需相互依赖,于是高斯消元dp。

debug:下载测试数据本地测,发现存在精度问题,将eps调到1e-10才过,本来是1e-8,但dx的std是1e-6,考虑高斯消元有什么问题吗,但板子都是测过的。疑似他没有进行每次找最大值交换操作?

此外还有的bug是存边数组开小导致RE..思路总体清晰,但是不熟,此外重构了高斯消元模板

const double eps = 1e-10;


double b[N][N]; //增广矩阵
vector<int>e[N];
int gauss(double a[][N],int n){
  int c,r;//当前列,当前行 
  for(c=r=1; c<=n; c++){
    //1.找到c列的最大行t
    int t=r; 
    for(int i=r; i<=n; i++)
      if(fabs(a[i][c])>fabs(a[t][c])) t=i;
    if(fabs(a[t][c])<eps) continue; //c列已0化
    
    //2.把最大行换到上面 
    for(int i=c; i<=n+1; i++)swap(a[t][i],a[r][i]);
    
    //3.把当前行r的第一个数,变成1   
    for(int i=n+1; i>=c; i--)a[r][i]/=a[r][c];
    
    //4.把当前列c下面的所有数,全部消成0    
    for(int i=r+1; i<=n; i++)
      if(fabs(a[i][c])>eps)
        for(int j=n+1; j>=c; j--)
          a[i][j]-=a[i][c]*a[r][j];
    r++; //这一行的工作做完,换下一行 
  }
  if(r<n+1){ //说明已经提前变成梯形矩阵
    for(int i=r; i<=n; i++) 
      if(fabs(a[i][n+1])>eps) //a[i][n]==bi 
        return 2; //左边=0,右边≠0,无解
    return 1; //0==0,无穷多解 
  }
 // 5.唯一解,从下往上回代,得到方程的解
  for(int i=n; i>=0; i--) 
    for(int j=i+1; j<n+1; j++)
      a[i][n+1]-=a[i][j]*a[j][n+1];
  return 0;
}
int in[N];
//f为经过点的期望次数
//ans:边的期望次数
//从1-n的期望得分=这个边的次数*编号,编号小的次数大

double res[M];
void solve(){
	int n,m;
	cin>>n>>m;
	vector<pii>eg(m+1);
	for(int i=1;i<=m;i++){
		int u,v;cin>>u>>v;
		e[u].push_back(v);
		e[v].push_back(u);
		in[u]++;in[v]++;
		eg[i]={u,v};
	}
	//高斯消元采取1base
	for(int i=1;i<=n-1;i++){
		for(auto v:e[i]){
			if(v!=n)b[i][v]=-1.0/in[v];
		}
		b[i][i]=1;
	}
	b[1][n]=1;
	// for(int i=1;i<=n;i++)cerr<<in[i]<<" ";
  //cerr<<endl;
	gauss(b,n-1);
	//for(int i=1;i<=n;i++)cerr<<b[i][n]<<" ";
	for(int i=1;i<=m;i++){
		auto [u,v]=eg[i];
		res[i]+=b[u][n]/(double)in[u]+b[v][n]/(double)in[v];
	}
	
	//cerr<<endl;
	//for(int i=1;i<=m;i++)cerr<<res[i]<<" ";
	sort(res+1,res+m+1);
	double ans=0;
	for(int i=1;i<=m;i++){
		ans+=res[i]*(double)(m+1-i);
	}
	baoliu(ans,3);
}

CF24D Broken robot

[Luogu P4457 BJOI2018]治疗之雨

HDU4418 Time Travel