树形DP+概率DP
动态规划报告
树形dp
树形 DP，即在树上进行的 DP。由于树固有的递归性质，树形 DP 一般都是递归进行的。一般需要在遍历树的同时维护所需的信息
以一道题目为例
2022CCPC桂林站G Group Homework
No, we don't want group homework. It's the place where 1+1<1 can be true. However, you are currently the leader of a group with three members. Luckily, your group members will write everything down for you since you are the prestigious leader. Unluckily, you have to assign the new algorithm homework to your two team members, Mr. Dian and Mr. Beng, who can't understand the ideas of each other correctly and mess up all the details in the cooperation.
The new homework is about a tree: there are n vertices on the tree with n−1 bidirectional edges connecting them. Each vertex i is a problem with a score of ai. You can assign a tree path to each member, then Mr. Dian and Mr. Beng will solve the problems on their path independently. The final score of the homework is decided by the sum of the scores of the vertices visited by exactly one member — as we mentioned before, if both of them solved one problem independently, they will argue a lot about who is right, and all the effort will be in vain.
Now, you — Mr. Ji — want to maximize the final score (as well as the chance not to drop out of school due to too low GPA). Make a wise choice!
题意为求树上两条链的不相交部分的最大权值和
题解
考虑两条链的关系，我们可以证明至多交于一点。若交于两点以上，我们可以修改为更大的答案。
所以可以确定要么是两条不相交的链，要么交于一点
- 对于两条不相交链，可以直接分类讨论
- 对于交于一点，可以枚举交点，用换根DP维护四条该点为根时的权值最大链
$mx[x][i]$:x为起点的权值最大的四条链
$dp[x][0]$:x为根节点子树中两条不相交链的最大和
$dp[x][1]$:x为根节点子树中最大的一条链 (不一定要经过x)
$dp[x][2]$:x为根节点子树中一条由x到叶子节点的链 加上一条与之不相交的链的最大和
$dp[x][3]$:x为根节点子树中的一条 不经过x 的最大链
$dp[x][4]$:x为根节点子树中一条由x到叶子节点的最大链
在具体实现时有很多细节，需要结合代码来理解
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=2e5+3;
int n;
int a[maxn];
vector<int> t[maxn];
ll dp[maxn][5],mx[maxn][5];
ll mx4;//最长的四条链之和
/*
mx[x][i]:x为起点的权值最大的四条链
dp[x][0]:x为根节点子树中两条不相交链的最大和
dp[x][1]:x为根节点子树中最大的一条链 __(不一定要经过x)__
dp[x][2]:x为根节点子树中一条由x到叶子节点的链 加上一条与之不相交的链的最大和
dp[x][3]:x为根节点子树中的一条 __不经过x__ 的最大链
dp[x][4]:x为根节点子树中一条由x到叶子节点的最大链
*/
void dfs(int x,int fa){
dp[x][0]=dp[x][1]=dp[x][2]=dp[x][4]=a[x];
dp[x][3]=0;
for(auto y:t[x]){
if(y==fa)continue;
dfs(y,x);
//由于有些状态转移不能在自己的一些状态转移之后转移 故需要注意转移状态的顺序
dp[x][0]=max(dp[x][0],dp[y][0]);
dp[x][0]=max(dp[x][0],dp[x][4]+dp[y][2]);//连x-y
dp[x][0]=max(dp[x][0],dp[x][2]+dp[y][4]);//连x-y
dp[x][0]=max(dp[x][0],dp[x][1]+dp[y][1]);//不连x-y
dp[x][1]=max(dp[x][1],dp[y][1]);
dp[x][1]=max(dp[x][1],dp[y][4]+dp[x][4]);
dp[x][2]=max(dp[x][2],dp[y][2]+a[x]);
dp[x][2]=max(dp[x][2],dp[x][4]+dp[y][1]);
dp[x][2]=max(dp[x][2],dp[x][3]+dp[y][4]+a[x]);
dp[x][3]=max(dp[x][3],dp[y][1]);
dp[x][4]=max(dp[x][4],dp[y][4]+a[x]);
}
}
void dfs4(int x,int fa,ll fmx){//fmx维护父节点除x节点外最长链加上a[fa]的值
mx[x][1]=mx[x][2]=mx[x][3]=mx[x][4]=0;
for(auto y:t[x]){
if(y==fa)continue;
if(dp[y][4]>=mx[x][1]){
mx[x][4]=mx[x][3];
mx[x][3]=mx[x][2];
mx[x][2]=mx[x][1];
mx[x][1]=dp[y][4];
}
else if(dp[y][4]>=mx[x][2]){
mx[x][4]=mx[x][3];
mx[x][3]=mx[x][2];
mx[x][2]=dp[y][4];
}
else if(dp[y][4]>=mx[x][3]){
mx[x][4]=mx[x][3];
mx[x][3]=dp[y][4];
}
else if(dp[y][4]>=mx[x][4]){
mx[x][4]=dp[y][4];
}
}
if(fmx>=mx[x][1]){
mx[x][4]=mx[x][3];
mx[x][3]=mx[x][2];
mx[x][2]=mx[x][1];
mx[x][1]=fmx;
}
else if(fmx>=mx[x][2]){
mx[x][4]=mx[x][3];
mx[x][3]=mx[x][2];
mx[x][2]=fmx;
}
else if(fmx>=mx[x][3]){
mx[x][4]=mx[x][3];
mx[x][3]=fmx;
}
else if(fmx>=mx[x][4]){
mx[x][4]=fmx;
}
mx4=max(mx4,mx[x][1]+mx[x][2]+mx[x][3]+mx[x][4]);
for(auto y:t[x]){
if(y==fa) continue;
if(dp[y][4]==mx[x][1]) dfs4(y,x,max(mx[x][2],fmx)+a[x]);
else dfs4(y,x,max(mx[x][1],fmx)+a[x]);
}
}
void solve(){
cin>>n;
for(int i=1;i<=n;i++) cin>>a[i];
if(n==1){
cout<<0;
return ;
}//特判1个节点
for(int i=1;i<n;++i){
int u,v;
cin>>u>>v;
t[u].push_back(v);
t[v].push_back(u);
}
dfs(1,0);
dfs4(1,0,0);
cout<<max(dp[1][0],mx4);
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
在桂林G题目中使用到了换根DP的方法，接下来介绍一下换根DP的一些用法
换根DP
树形 DP 中的换根 DP 问题又被称为二次扫描，通常不会指定根结点，并且根结点的变化会对一些值，例如子结点深度和、点权和等产生影响。
通常需要两次 DFS，第一次 DFS 预处理诸如深度，点权和之类的信息，在第二次 DFS 开始运行换根动态规划。
还是以一道题目为例 [POI2008] STA-Station
题目描述
给定一个 $n$ 个点的树，请求出一个结点，使得以这个结点为根时，所有结点的深度之和最大。
一个结点的深度之定义为该节点到根的简单路径上边的数量。
如果有多个结点符合要求，输出任意一个即可。
对于全部的测试点，保证 $1 \leq n \leq 10^6,1\leq u, v\leq n$，给出的是一棵树。
题解
以 $siz[x]$ 记录以1为根时x子树大小，先通过一次dfs预处理出 $siz$ 数组
以 $dp[x]$ 记录以 x为根 时所有节点的总深度
$dp[x]\rightarrow dp[y]$ 的转移即为换根的过程。显然在换根的转移过程中，以 $y$ 为根或以 $x$ 为根会导致其子树中的结点的深度产生改变。
其中y的子树结点深度都减少一，其余节点深度增加一，即转移方程为
$dp[y]=dp[x]-siz[y]+n-siz[y]=dp[x]+n-2 \times siz[y]$
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e6+3;
int n;
ll siz[maxn],dp[maxn],dep[maxn],ans=0;
vector<int> t[maxn];
void pre(int x,int fa){
dep[x]=dep[fa]+1;
siz[x]=1;
for(auto y:t[x]){
if(y==fa) continue ;
pre(y,x);
siz[x]+=siz[y];
}
}
void dfs(int x,int fa){
ans=max(ans,dp[x]);
for(auto y:t[x]){
if(y==fa) continue ;
dp[y]=dp[x]+n-2*siz[y];
dfs(y,x);
}
}
void solve(){
cin>>n;
for(int i=1;i<n;i++){
int u,v;
cin>>u>>v;
t[u].push_back(v);
t[v].push_back(u);
}
dep[0]=0;
pre(1,0);
for(int i=1;i<=n;i++) dp[1]+=dep[i]-1;
dfs(1,0);
for(int i=1;i<=n;i++) if(dp[i]==ans){cout<<i;return ;}
return ;
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
树上背包
除了换根DP，树上背包也是树形DP中的一种常见问题
树上的背包问题，简单来说就是背包问题与树形 DP 的结合。
以 洛谷 P2014 CTSC1997 选课 为例
题目描述
在大学里每个学生，为了达到一定的学分，必须从很多课程里选择一些课程来学习，在课程里有些课程必须在某些课程之前学习，如高等数学总是在其它课程之前学习。现在有 $N$ 门功课，每门课有个学分，每门课有一门或没有直接先修课（若课程 $a$ 是课程 $b$ 的先修课即只有学完了课程 $a$，才能学习课程 $b$）。一个学生要从这些课程里选择 $M$ 门课程学习，问他能获得的最大学分是多少？
输入格式
第一行有两个整数 $N , M$ 用空格隔开。($ 1 \leq N \leq 300, 1 \leq M \leq 300$)
接下来的 $N$ 行,第 $I+1$ 行包含两个整数 $k_i$和 $s_i$, $k_i$ 表示第$I$门课的直接先修课，$s_i$ 表示第$I$门课的学分。若 $k_i=0$ 表示没有直接先修课（$1\leq k_i \leq N, 1 \leq s_i \leq 20$）。
题解
根据题意每门课至多只有一门先修课，那么不妨添加第 $0$ 门课程作为没有先修课的先修课，这样就可以以先修关系为边，建立一棵0为根的树。需要注意最后统计结果时要统计 $m+1$ 个节点(包含 $0$ 号根节点)
设 $dp[x][i][j]$ 表示 $x$ 号点的子树中前 $i$ 棵子树选择 $j$ 门课程的最大学分
转移的过程就是在遍历树的过程中进行背包DP的转移。枚举 $x$ 的每个子节点 $y$ ，同时枚举以 $y$ 为根的子树选择了几门课程，将子树的结果合并到 $x$ 上
设 $siz[x]$ 表示 $x$ 的子树大小，可以写出转移方程
$dp[x][i][j]=\max {k<j,k<siz[y]} dp[x][i-1][j-k] + dp[y][siz][k]$
其中表示选择前 $i$ 棵子树的这一维可以用01背包相同的滚动数组方式省略掉
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e3+10;
int n,m;
int s[maxn],f[maxn];
int dp[maxn][maxn],siz[maxn];
vector<int> t[maxn];
void dfs(int x){
dp[x][1]=s[x];
siz[x]=1;
for(auto y:t[x]){
dfs(y);
siz[x]+=siz[y];
for(int i=m+1;i;i--){
for(int j=1;j<=siz[y];j++){
dp[x][i+j]=max(dp[x][i+j],dp[x][i]+dp[y][j]);
}
}
}
}
void solve(){
cin>>n>>m;
for(int i=1;i<=n;i++){
cin>>f[i]>>s[i];
t[f[i]].push_back(i);
}
dfs(0);
cout<<dp[0][m+1];
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
再来看一道树上背包的题目
2022CCPC广州站I. Infection
A highly propagating bacterium infects a tree of n nodes (with $n−1$ edges, no cycles). These nodes are indexed from $1$ to $n$.
Exactly one node will be infected at the beginning. Each node on the tree has an initial susceptibility weight ai, which represents that node i has a probability of $\frac{ a_i }{\Sigma ^n_{i=1}ai}$ to become the initial infected node of the tree.
In addition, each node has an infection probability $p_i$, which represents its probability of being infected by adjacent nodes.
Specifically, starting from the initial infected node, if a node is infected, the uninfected node $x$ that is adjacent to it will have a probability of $p_x$ to become a new infected node, and then $x$ can continue to infect its adjacent nodes.
Now, your task is to calculate the probability that exactly k nodes are eventually infected. You need to output an answer for each $k=1,…,n.$
You need to output the answer modulo $10^9+7,$ which means if your answer is $\frac a b (gcd(a,b)=1)$ , you need to output $a⋅b{−1}mod109+7$, where $b⋅b^{-1}\equiv 1 (mod 10^9+7)$.
Input
The first line contains an integer $n (2\leq n\leq 2000)$, denoting the number of nodes of the tree.
The next $n−1$ lines, each line contains two positive integers $u$ and $v (1\leq u,v\leq n)$, denoting that there is an edge $(u,v)$ on the tree.
Next n lines, the i-th line contains three positive integers $a_i,b_i,c_i$ , where ai means as above and $p_i=\frac{b_i}{c_i}$. It is guaranteed that $1\leq a_i,b_i,c_i\leq 10^9,\sum _{i=1}^n a_i\leq 10^9,b_i\leq c_i,gcd(b_i,c_i)=1$.
Output
Output $n$ lines, each line contains single integer. The integer on the $i$-th line should be the answer modulo $10^9+7$ when $k=i$.
题解
$dp[x][i][1/0]$ 表示 $x$ 的子树感染了 $i$ 个点，是否包含初始感染点的概率之和, $siz[x]$ 表示 $x$ 子树大小
由于感染是相邻的，所以 $j>0$ 时 $x$ 点一定是被感染了的点
$\frac{a[x]}{\sum a}$ 为 $x$ 点初始感染概率，$p[x]=\frac{b[x]}{c[x]}$ 为被传染概率
所以在初始化背包时$dp[x][1][0]=\frac{a[x]}{\sum a},dp[x][1][1]=p[x]$
$dp[x][0][1]$ 是不合法状态，概率为 $0$
$dp[x][0][0]=1-p[x]$ ,该子树不被感染只需要该子树根节点 $x$ 不被传染
$dp[y]\rightarrow dp[x]$ 转移时根据初始感染点的位置分三种情况
- 初始感染点不在 $x$ 子树时
$dp[x][i+j][0]=\sum _{i\leq siz[x],j\leq siz[y]}dp[x][i][0]*dp[y][j][0]$ - 初始感染点在 $x$ 子树但不在 $y$ 子树时
$dp[x][i+j][1]=\sum _{i\leq siz[x],j\leq siz[y]}dp[x][i][1]*dp[y][j][0]$ - 初始感染点在 $y$ 子树时
$dp[x][i+j][1]=\sum _{i\leq siz[x],j\leq siz[y]}dp[x][i][0]*dp[y][j][0]$
在统计答案时分别统计只感染 $x$ 子树的贡献，即$dp[x]*(1-p[fa])$
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=2e3+10;
const ll mod=1e9+7;
int n;
int siz[maxn];
ll a[maxn],b[maxn],c[maxn],p[maxn],suma;
ll dp[maxn][maxn][2],ans[maxn];
vector<int> t[maxn];
ll inv(ll x){
ll k=mod-2,base=x;
ll ans=1;
while(k){
if(k&1) (ans*=base)%=mod;
(base*=base)%=mod;
k>>=1;
}
return ans;
}
ll f[maxn][2];
void dfs(int x,int fa){
dp[x][1][0]=p[x];
dp[x][1][1]=a[x]*suma%mod;
siz[x]=1;
for(auto y:t[x]){
if(y==fa) continue ;
dfs(y,x);
for(int i=1;i<=siz[x]+siz[y];i++) f[i][1]=f[i][0]=0;
for(int i=1;i<=siz[x];i++){
for(int j=0;j<=siz[y];j++){
(f[i+j][0]+=dp[x][i][0]*dp[y][j][0])%=mod;
(f[i+j][1]+=dp[x][i][1]*dp[y][j][0])%=mod;
if(j) (f[i+j][1]+=dp[x][i][0]*dp[y][j][1])%=mod;
}
}
siz[x]+=siz[y];
for(int i=1;i<=siz[x];i++){
dp[x][i][0]=f[i][0];
dp[x][i][1]=f[i][1];
}
}
for(int i=1;i<=siz[x];i++){
(ans[i]+=(1-p[fa]+mod)%mod*dp[x][i][1])%=mod;
}
dp[x][0][0]=(1-p[x]+mod)%mod;
}
void solve(){
cin>>n;
for(int i=1;i<n;i++){
int u,v;
cin>>u>>v;
t[u].push_back(v);
t[v].push_back(u);
}
suma=0;
for(int i=1;i<=n;i++){
cin>>a[i]>>b[i]>>c[i];
suma+=a[i];
p[i]=b[i]*inv(c[i])%mod;
}
suma=inv(suma%mod);
dfs(1,0);
for(int i=1;i<=n;i++) cout<<ans[i]<<endl;
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
概率DP
DP求概率
这类题目采用顺推，也就是从初始状态推向结果。同一般的 DP 类似的，难点依然是对状态转移方程的刻画，只是这类题目经过了概率论知识的包装。
以 Codeforces 148 D Bag of mice 为例
The dragon and the princess are arguing about what to do on the New Year's Eve. The dragon suggests flying to the mountains to watch fairies dancing in the moonlight, while the princess thinks they should just go to bed early. They are desperate to come to an amicable agreement, so they decide to leave this up to chance.
They take turns drawing a mouse from a bag which initially contains w white and b black mice. The person who is the first to draw a white mouse wins. After each mouse drawn by the dragon the rest of mice in the bag panic, and one of them jumps out of the bag itself (the princess draws her mice carefully and doesn't scare other mice). Princess draws first. What is the probability of the princess winning?
If there are no more mice in the bag and nobody has drawn a white mouse, the dragon wins. Mice which jump out of the bag themselves are not considered to be drawn (do not define the winner). Once a mouse has left the bag, it never returns to it. Every mouse is drawn from the bag with the same probability as every other one, and every mouse jumps out of the bag with the same probability as every other one.
Input
The only line of input data contains two integers $w$ and $b (0 \leq w, b \leq 1000)$.
Output
Output the probability of the princess winning. The answer is considered to be correct if its absolute or relative error does not exceed $10^{-9}$.
题目大意：袋子里有 $w$ 只白鼠和 $b$ 只黑鼠，公主和龙轮流从袋子里抓老鼠。谁先抓到白色老鼠谁就赢，如果袋子里没有老鼠了并且没有谁抓到白色老鼠，那么算龙赢。公主每次抓一只老鼠，龙每次抓完一只老鼠之后会有一只老鼠跑出来。每次抓的老鼠和跑出来的老鼠都是随机的。公主先抓。问公主赢的概率。
题解
设 $dp[i][j]$ 表示 轮到公主 时袋子里有 $i$ 只白鼠，$j$ 只黑鼠，公主赢的概率。对抓到的老鼠进行分类讨论，考虑转移方程。
- 公主抓到白鼠，公主获胜。$dp[i][j]+=\frac{i}{i+j}$
- 公主抓到黑鼠，龙抓到白鼠，龙赢了。无贡献
- 都抓到黑鼠，跑了一只黑鼠，$dp[i][j] \leftarrow dp[i][j-3] , dp[i][j]+=\frac{j}{i+j}\times \frac{j-1}{i+j-1}\times \frac{j-2}{i+j-2}\times dp[i][j-3]$
- 都抓到黑鼠，跑了一只白鼠，$dp[i][j] \leftarrow dp[i-1][j-2] , dp[i][j]+=\frac{j}{i+j}\times \frac{j-1}{i+j-1}\times \frac{i}{i+j-2}\times dp[i-1][j-2]$
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e3+10;
double dp[maxn][maxn];
void solve(){
int w,b;
cin>>w>>b;
for(int i=1;i<=w;i++) dp[i][0]=1;
for(int i=1;i<=b;i++) dp[0][i]=0;
for(int i=1;i<=w;i++){
for(int j=1;j<=b;j++){
dp[i][j]+=(double)i/(i+j);
if(j>=3) dp[i][j]+=(double)j/(i+j)*(j-1)/(i+j-1)*(j-2)/(i+j-2)*dp[i][j-3];
if(i>=1&&j>=2) dp[i][j]+=(double)j/(i+j)*(j-1)/(i+j-1)*i/(i+j-2)*dp[i-1][j-2];
}
}
cout<<fixed<<setprecision(10)<<dp[w][b];
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
DP求期望
POJ2096 Collecting Bugs
Description
Ivan is fond of collecting. Unlike other people who collect post stamps, coins or other material stuff, he collects software bugs. When Ivan gets a new program, he classifies all possible bugs into n categories. Each day he discovers exactly one bug in the program and adds information about it and its category into a spreadsheet. When he finds bugs in all bug categories, he calls the program disgusting, publishes this spreadsheet on his home page, and forgets completely about the program.
Two companies, Macrosoft and Microhard are in tight competition. Microhard wants to decrease sales of one Macrosoft program. They hire Ivan to prove that the program in question is disgusting. However, Ivan has a complicated problem. This new program has s subcomponents, and finding bugs of all types in each subcomponent would take too long before the target could be reached. So Ivan and Microhard agreed to use a simpler criteria --- Ivan should find at least one bug in each subsystem and at least one bug of each category.
Macrosoft knows about these plans and it wants to estimate the time that is required for Ivan to call its program disgusting. It's important because the company releases a new version soon, so it can correct its plans and release it quicker. Nobody would be interested in Ivan's opinion about the reliability of the obsolete version.
A bug found in the program can be of any category with equal probability. Similarly, the bug can be found in any given subsystem with equal probability. Any particular bug cannot belong to two different categories or happen simultaneously in two different subsystems. The number of bugs in the program is almost infinite, so the probability of finding a new bug of some category in some subsystem does not reduce after finding any number of bugs of that category in that subsystem.
Find an average time (in days of Ivan's work) required to name the program disgusting.
Input
Input file contains two integer numbers, n and s (0 < n, s <= 1 000).
Output
Output the expectation of the Ivan's working days needed to call the program disgusting, accurate to 4 digits after the decimal point.
题目大意：一个软件有 $s$ 个子系统，会产生 $n$ 种 bug。某人一天发现一个 bug，这个 bug 属于某种 bug 分类，也属于某个子系统。每个 bug 属于某个子系统的概率 $\frac{1}{s}$，属于某种 bug 分类的概率是 $\frac{1}{n}$。求发现 $n$ 种 bug，且 $s$ 个子系统都找到 bug 的期望天数。
题解
设$dp[i][j]$ 为已经找到了 $i$ 种 bug，属于$j$ 种分类时的期望天数。显然 $dp[n][s]=0$,题目所求答案为 $dp[0][0]$
对新发现的 bug 进行分类讨论，考虑状态转移
- 若该 bug 种类已知，分类已知,$dp[i][j]\rightarrow dp[i][j]$：概率为$\frac{i}{n}\times\frac{j}{s}$
- 若该 bug 种类已知，分类未知,$dp[i][j]\rightarrow dp[i][j+1]$：概率为$\frac{i}{n}\times\frac{s-j}{s}$
- 若该 bug 种类未知，分类已知,$dp[i][j]\rightarrow dp[i+1][j]$：概率为$\frac{n-i}{n}\times\frac{j}{s}$
- 若该 bug 种类未知，分类未知,$dp[i][j]\rightarrow dp[i+1][j+1]$：概率为$\frac{n-i}{n}\times\frac{s-j}{s}$
由期望的线性性质，可以得到方程
$dp[i][j]=\frac{i}{n}\times\frac{j}{s}\times (dp[i][j]+1)+\frac{i}{n}\times\frac{s-j}{s}\times (dp[i][j+1]+1)+\frac{n-i}{n}\times\frac{j}{s}\times (dp[i+1][j]+1)+\frac{n-i}{n}\times\frac{s-j}{s}\times (dp[i+1][j+1]+1)$
化简后可得$dp[i][j]=\frac{\frac{i}{n}\times\frac{s-j}{s}\times dp[i][j+1]+\frac{n-i}{n}\times\frac{j}{s}\times dp[i+1][j]+\frac{n-i}{n}\times\frac{s-j}{s}\times dp[i+1][j+1]+1}{1-\frac{i}{n}\times\frac{j}{s}}$
$=\frac{i\times(s-j)\times dp[i][j+1]+(n-i)\times i\times dp[i+1][j]+(n-i)\times(s-j)\times dp[i+1][j+1] + n*s}{n\times s-i\times j}$
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e3+10;
double dp[maxn][maxn];
void solve(){
int n,s;
cin>>n>>s;
for(int i=n;i>=0;i--){
for(int j=s;j>=0;j--){
if(i==n&&j==s) continue ;
dp[i][j]=(double)(i*(s-j)*dp[i][j+1]+(n-i)*j*dp[i+1][j]+(n-i)*(s-j)*dp[i+1][j+1]+n*s)/(n*s-i*j);
}
}
cout<<fixed<<setprecision(10)<<dp[0][0];
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
// cin>>T;
while(T--){
solve();
}
return 0;
}
可以看出，与DP概率相比，DP求期望在得到期望的方程后往往还需要进行一些变化才能变成最终的转移方程
Codeforces Round #848 (Div. 2) D. Flexible String Revisit
You are given two binary strings $a$ and $b$ of length $n$. In each move, the string $a$ is modified in the following way.
An index $i (1\leq i\leq n)$ is chosen uniformly at random. The character ai will be flipped. That is, if $a_i$ is $0$, it becomes $1$, and if $a_i$ is $1$, it becomes $0$.
What is the expected number of moves required to make both strings equal for the first time?
A binary string is a string, in which the character is either $0$ or $1$.
Input
The first line contains a single integer $t (1\leq t\leq 10^5)$ — the number of test cases. The description of the test cases follows.
The first line of each test case contains a single integer $n (1\leq n\leq 10^6)$ — the length of the strings.
The second line of each test case contains the binary string $a$ of length $n$.
The third line of each test case contains the binary string $b$ of length $n$.
It is guaranteed that the sum of n over all test cases does not exceed $10^6$.
Output
For each test case, output a single line containing the expected number of moves modulo $998244353$.
Formally, let $M=998244353$. It can be shown that the answer can be expressed as an irreducible fraction $\frac{p}{q}$, where $p$ and $q$ are integers and $q \equiv 0(modM)$. Output the integer equal to $p⋅q^{−1}modM$. In other words, output such an integer $x$ that $0\leq x\lneq M$ and $x⋅q≡p(modM)$.
题解
设 $dp[i]$ 表示有 $i$ 个位置不同时的期望
显然 $dp[0]=0$ 已经到达目标状态, $dp[n]=1+dp[n-1]$ 所有位置都不同
考虑对于 $dp[i]$
- 若选择不同的位置，$dp[i]\rightarrow dp[i-1]$ ,概率为 $\frac{i}{n}$
- 若选择相同的位置，$dp[i]\rightarrow dp[i+1]$ ,概率为 $\frac{n-i}{n}$
可以得到 $dp[i]=\frac{i}{n}\times (dp[i-1]+1)+\frac{n-i}{n}\times (dp[i+1]+1)$
令 $x=i-1$ 化简可得转移方程 $dp[x]=\frac{n\times dp[x-1]-x\times dp[x-2]-n}{n-x-1}$
但是这个方程得转移需要两个初始状态 $dp[0]$ 和 $dp[1]$，而计算出 $dp[1]=2^n-1$ 并不容易
CF官方题解给出的办法是设 $dp[i]=a[i]+b[i]\times dp[i+1]$
将 $dp[i]=\frac{i}{n}\times (dp[i-1]+1)+\frac{n-i}{n}\times (dp[i+1]+1)$ 式中的 $dp[i-1]$ 替换为 $a[i-1]+b[i-1]\times dp[i]$ ,从而得到 $dp[i]\rightarrow dp[i+1]$ 的方程
$dp[i]=\frac{i}{n}\times (a[i-1]+b[i-1]\times dp[i])+\frac{n-i}{n}\times dp[i+1]+1$
化简后得到 $dp[i]=\frac{a[i-1]\times i+n }{n-i\times b[i-1]}+\frac{n-i}{n-i\times b[i-1]}\times dp[i+1]$
已知 $dp[0]=0$ 所以有 $dp[1]=1+\frac{n-1}{n}\times dp[2]$
可以得到 $a[1]=1,b[1]=\frac{n-1}{n}$
再通过$a[i]=\frac{a[i-1]\times i+n }{n-i\times b[i-1]},b[i]=\frac{n-i}{n-i\times b[i-1]}$，就可以求得$a[i],b[i](2\leq i\leq n)$
由于无法求得 $a[0]$ 和 $b[0]$，就仍然无法得到$dp[0]\rightarrow dp[1]$
观察 $dp[n]=1+dp[n-1]$，与之前类似地可以设 $dp[i]=c[i]+d[i]\times dp[i-1]$
将 $dp[i]=\frac{i}{n}\times (dp[i-1]+1)+\frac{n-i}{n}\times (dp[i+1]+1)$ 式中的 $dp[i+1]$ 替换为 $c[i+1]+d[i+1]\times dp[i]$，从而得到 $dp[i]\rightarrow dp[i-1]$ 的方程
$dp[i]=\frac{i}{n}\times dp[i-1]+\frac{n-i}{n}\times(c[i+1]+d[i+1]\times dp[i])+1$
化简得 $dp[i]=\frac{n+(n-i)\times c[i+1]}{n-(n-i)\times d[i+1]}+\frac{i}{n-(n-i)\times d[i+1]}\times dp[i-1]$
类似得，通过 $c[n]=1,d[n]=1$ 和 $c[i]=\frac{n+(n-i)\times c[i+1]}{n-(n-i)\times d[i+1]},d[i]=\frac{i}{n-(n-i)\times d[i+1]}$，可以求得 $c[i],d[i]$
$$
\begin{cases}
dp[i]=\frac{a[i-1]\times i+n }{n-i\times b[i-1]}+\frac{n-i}{n-i\times b[i-1]}\times dp[i+1]\
dp[i]=\frac{n+(n-i)\times c[i+1]}{n-(n-i)\times d[i+1]}+\frac{i}{n-(n-i)\times d[i+1]}\times dp[i-1]
\end{cases}
$$
解方程组可得 $dp[i]=\frac{c[i]+d[i]\times a[i-1]}{1-d[i]\times b[i-1]}$
由于 $dp[1]=c[1]+d[1]\times dp[0]=c[1]$，代入 $dp[i]=\frac{c[i]+d[i]\times a[i-1]}{1-d[i]\times b[i-1]}$ 可以得到 $a[0]+b[0]=0$
又因为 $dp[0]=a[0]+b[0]\times dp[1]$，即$a[0]+b[0]*c[1]=0$
联立得
$$
\begin{cases}
当 n=1 时，c[1]=1 , a[0],b[0]为任意实数\
当 n\neq 1 时，a[0]=0,b[0]=0
\end{cases}
$$
所以可以把 $a[0]=0,b[0]=0$作为初始化
code
#include <bits/stdc++.h>
using namespace std;
#define mod 998244353
#define N 1000005
template<int MOD>
struct ModInt {
unsigned x;
ModInt() : x(0) { }
ModInt(signed sig) : x(sig) { }
ModInt(signed long long sig) : x(sig%MOD) { }
int get() const { return (int)x; }
ModInt pow(ll p) { ModInt res = 1, a = *this; while (p) { if (p & 1) res *= a; a *= a; p >>= 1; } return res; }
ModInt &operator+=(ModInt that) { if ((x += that.x) >= MOD) x -= MOD; return *this; }
ModInt &operator-=(ModInt that) { if ((x += MOD - that.x) >= MOD) x -= MOD; return *this; }
ModInt &operator*=(ModInt that) { x = (unsigned long long)x * that.x % MOD; return *this; }
ModInt &operator/=(ModInt that) { return (*this) *= that.pow(MOD - 2); }
ModInt operator+(ModInt that) const { return ModInt(*this) += that; }
ModInt operator-(ModInt that) const { return ModInt(*this) -= that; }
ModInt operator*(ModInt that) const { return ModInt(*this) *= that; }
ModInt operator/(ModInt that) const { return ModInt(*this) /= that; }
bool operator<(ModInt that) const { return x < that.x; }
friend ostream& operator<<(ostream &os, ModInt a) { os << a.x; return os; }
};
typedef ModInt<998244353> mint;
mint a[N], b[N], c[N], d[N];
int main(){
ios::sync_with_stdio(false);
cin.tie(0);
int t;
cin >> t;
while(t--){
int n;
cin >> n;
string s1, s2;
cin >> s1 >> s2;
int diff = 0;
for(int i = 0; i < n; i++){
diff += s1[i]!=s2[i];
}
c[n] = 1, d[n] = 1, a[1] = 1, b[1] = ((mint) n-1)/n;
for(int i = 2; i <= n; i++){
a[i] = ((mint) n + a[i-1]*i)/((mint) n - b[i-1]*i);
b[i] = ((mint) n-i)/((mint) n - b[i-1]*i);
}
for(int i = n-1; i >= 1; i--){
c[i] = ((mint) n + c[i+1]*(n-i))/((mint) n - d[i+1]*(n-i));
d[i] = (mint) i / ((mint) n - d[i+1]*(n-i));
}
mint ans = (c[diff]+d[diff]*a[diff-1])/((mint) 1 - d[diff]*b[diff-1]);
cout << ans << endl;
}
return 0;
}
换一种倒推的DP方式将更加简单
设 $dp[i]$ 为从 $i$ 个不同位置到 $i-1$ 个不同位置的期望
- 若选择不同位置，期望为 $\frac{i}{n}$
- 若选择相同位置，期望为 $\frac{n-i}{n}\times(1+dp[i+1]+dp[i])$
意为这次选择到达了 $dp[i+1]$ 的状态，需要再选 $dp[i+1]$ 次回到 $dp[i]$
可以直接得到转移方程
$dp[i]=\frac{i}{n}+\frac{n-i}{n}\times(dp[i]+dp[i+1]+1)$
$dp[i]=\frac{(n-i)\times dp[i+1]+n}{i}$
初始状态 $dp[n]=1$ ，因为所有位置都不同
最后的结果即为 $\sum^k_1 dp[i]$ ，$k$ 为不同位置数
code
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int maxn=1e6+10;
const ll mod=998244353;
ll dp[maxn];
ll getinv(ll x){
ll k=mod-2;
ll inv=1;
while(k){
if(k&1) (inv*=x)%=mod;
(x*=x)%=mod;
k>>=1;
}
return inv;
}
void solve(){
int n;
string a,b;
cin>>n>>a>>b;
int dif=0;
for(int i=0;i<n;i++) if(a[i]!=b[i]) dif++;
dp[n]=1;
for(int i=n-1;i>=0;i--){
dp[i]=(n+(n-i)*dp[i+1])%mod*getinv(i)%mod;
}
ll ans=0;
for(int i=1;i<=dif;i++) (ans+=dp[i])%=mod;
cout<<ans<<endl;
}
int main(){
ios::sync_with_stdio(false);
cin.tie(nullptr);
int T=1;
cin>>T;
while(T--){
solve();
}
return 0;
}
