题面：
有一个初始为空的序列，在序列末尾随机添加1或2，有p的概率添加1，1 - p的概率添加2，如果序列末尾有连着的两个相同的数k，那么他们会合并成k + 1，这个合并只要可行，可以一直持续下去
序列有个长度限制n，当序列凑到n位且无法合并，添加操作就结束，问结束时所有权值的和的期望。
solution:
假设数字i出现的概率为p(i)，则数字i + 1出现的期望可以粗略地估计为p(i)^2，那么，值大于50的数出现的概率就非常小了，对答案的贡献可以忽略不计
定义a[i][j]:在长度限制为i的情况下，最左端出现过数字j，最后让序列长度停止在i的概率
那么有a[i][j] = a[i][j - 1] * a[i - 1][j - 1]，即先出现一次j - 1，再出现j - 1的概率
类似地，定义b[i][j]，但加一个约束，序列初始时的第一个数字必须是2，有b[i][j] = b[i][j - 1] * a[i - 1][j - 1]
边界为a[i][1] = p，a[i][2] = b[i][2] = 1 - p
定义f[i][j]:构造一个长度为i的序列，构造结束时序列最左端为j，序列的和的期望
一般地，有f[i][j] = j + ∑(f[i - 1][k] * a[i - 1][k] * (1 - a[i - 2][k])) / ∑(a[i - 1][k] * (1 - a[i - 2][k])) (k = 1 ~ j - 1)
因为第一位确定为j时，第二位出现过的数字一定不能大于等于j，否则会和这个j合并
对于第三位，因为第二位的值已经确定，所以只需保证第三位不为k即可，乘上概率转移一下
特别地，当j = 1，
f[i][1] = 1 + ∑(f[i - 1][k] * b[i - 1][k] * (1 - a[i - 2][k])) / ∑(b[i - 1][k] * (1 - a[i - 2][k])) (k = 2 ~ 50)
最后统计答案为Ans = ∑f[n][i] * a[n][i] * (1 - a[n][i])
注意到转移超过五十次时，a，b数组的值就不在更新了
所以暴力前50次转移，后面的转移就可以用矩阵乘法了
#include<iostream> #include<cstdio> #include<algorithm> #include<cmath> #include<cstring> #include<vector> #include<queue> #include<set> #include<map> #include<stack> #include<bitset> #include<ext/pb_ds/priority_queue.hpp> using namespace std; const int N = 51; typedef long double DB; const DB One = 1.000; struct data{ DB a[N][N]; data(){memset(a,0,sizeof(a));} data operator * (const data &b) { data c; for (int k = 0; k < N; k++) for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) c.a[i][j] += a[i][k] * b.a[k][j]; return c; } }; int n; DB Ans,p,q,a[N][N],b[N][N],f[N][N]; void Calc_50() { for (int i = 2; i < N; i++) { for (int j = 1; j < N; j++) { if (j == 1) a[i][j] += p; if (j == 2) a[i][j] += q,b[i][j] += q; a[i][j] += a[i][j - 1] * a[i - 1][j - 1]; b[i][j] += b[i][j - 1] * a[i - 1][j - 1]; } DB s1,s2; s1 = s2 = 0; for (int k = 2; k <= 50; k++) { s2 += b[i - 1][k] * (One - a[i - 2][k]); s1 += f[i - 1][k] * b[i - 1][k] * (One - a[i - 2][k]); } f[i][1] = One + s1 / s2; for (int j = 2; j <= 50; j++) { s1 = s2 = 0; for (int k = 1; k < j; k++) { s2 += a[i - 1][k] * (One - a[i - 2][k]); s1 += f[i - 1][k] * a[i - 1][k] * (One - a[i - 2][k]); } f[i][j] = j + s1 / s2; } } } data ksm(data x) { data ret; for (int i = 0; i < N; i++) ret.a[i][i] = 1; for (; n; n >>= 1) { if (n & 1) ret = ret * x; x = x * x; } return ret; } void Calc() { data k; DB sum = 0; f[50][0] = k.a[0][0] = 1; for (int i = 2; i <= 50; i++) sum += b[50][i] * (One - a[50][i]); k.a[0][1] = 1; for (int i = 2; i <= 50; i++) k.a[i][1] = b[50][i] * (One - a[50][i]) / sum; for (int i = 2; i <= 50; i++) { sum = 0; k.a[0][i] = (DB)(i); for (int j = 1; j < i; j++) sum += a[50][j] * (One - a[50][j]); for (int j = 1; j < i; j++) k.a[j][i] = a[50][j] * (One - a[50][j]) / sum; } k = ksm(k); for (int i = 1; i <= 50; i++) for (int j = 0; j <= 50; j++) Ans += f[50][j] * k.a[j][i] * a[50][i] * (One - a[50][i]); printf("%.12lf\n",(double)(Ans)); } int main() { #ifdef DMC freopen("DMC.txt","r",stdin); #endif cin >> n >> p; p /= 1E9; q = One - p; a[1][1] = p; a[1][2] = b[1][2] = q; f[1][1] = 1; f[1][2] = 2; Calc_50(); if (n <= 50) { for (int i = 1; i <= 50; i++) Ans += f[n][i] * a[n][i] * (One - a[n - 1][i]); printf("%.12lf\n",(double)(Ans)); return 0; } n -= 50; Calc(); return 0; }