HDU 4870 Rating 概率DP(高斯消元)

题意:初始rating为0,每次比赛有p的概率变为min(x+50,1000),(1-p)的概率变为(x-100,0).
0.3<=p<=1.0,某人用两个账号 每次用rating小的账号比赛,问有1个账号rating为1000时,比赛次数的期望?

S的下一个状态为S'=ax+by,E(ax+by)=aE(x)+bE(y)
rating不超过1000 设dp[a][b] rating为a,b时 有一个号等于1000时 参加比赛次数的期望
因为每次用较小的一个,比赛的最后一步肯定为(950,950)->(1000,950) 每个状态又都是50的倍数,获胜一次50分离散化成1分 则设dp[i] i->20的期望次数.

dp[i]=1+p*dp[i+1]+(1-p)*dp[i-2] 有环 因为只有20个方程,移项后高斯消元,最后答案为2x[0]-x[19].

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=25;
double p,a[N][N],x[N];
int equ,var;
void init()
{
	memset(a,0,sizeof(a));
	//d[i]=1+p*d[i+1]+(1-p)*d[i-2]
	//(p-1)d[i-2]+d[i]-pd[i+1]=1
	for(int i=0;i<=20;i++)
	{
		x[i]=1;
		a[i][i]=1;
		if(i<20)
			a[i][i+1]=-p;
		if(i>=2)
			a[i][i-2]=p-1.0;
		else
			a[i][0]+=p-1.0;
	}
//	cout<<a[0][0]<<' '<<a[1][0]<<endl;
	equ=var=20;
}
int Gauss()
{
	for(int k=0,col=0;k<equ&&col<var;k++,col++)
	{
		int r=k;
		//精度
		for(int i=k+1;i<equ;i++)
			if(fabs(a[i][col])>fabs(a[r][col]))
				r=i;
		//if(fabs(a[r][col]<eps) return 0;		
		if(r!=k)
		{
			for(int j=col;j<var;j++)
				swap(a[k][j],a[r][j]);
			swap(x[k],x[r]);	
		}
		//a[k][col]->1
		x[k]/=a[k][col];
		for(int j=col+1;j<var;j++)
			a[k][j]/=a[k][col];
		a[k][col]=1;
		
		for(int i=0;i<equ;i++)
		{
			if(i!=k)
			{
				//a[i][col]消为0
				x[i]-=x[k]*a[i][col];
				for(int j=col+1;j<var;j++)
					a[i][j]-=a[k][j]*a[i][col];
				a[i][col]=0;
			}
		}
	}
	return 1;
}
int main()
{
	while(cin>>p)
	{
		init();
		Gauss();
		double ans=2.0*x[0]-x[19];
		printf("%.7lf\n",ans);
	}
	return 0;
}