C++14(g++5.4) 解法, 执行用时: 1360ms, 内存消耗: 376K, 提交时间: 2020-01-26 15:55:58

#include <cstdio>
#include <cstring>
const int N=25;
struct Matrix{double a[N][N];};
using namespace std;
inline Matrix Mul(Matrix m1,Matrix m2){
    Matrix ret;
    memset(ret.a,0,sizeof(ret.a));
    for(int i=1;i<N;i++)for(int k=1;k<N;k++){
        if(m1.a[i][k]==0) continue;
        for(int j=1;j<N;j++) ret.a[i][j]+=m1.a[i][k]*m2.a[k][j];
    }return ret;
}inline Matrix ksm(Matrix m,int p){
    Matrix ret;
    memset(ret.a,0,sizeof(ret.a));
    for(int i=1;i<N;i++) ret.a[i][i]=1;
    while(p){if(p&1)ret=Mul(ret,m);m=Mul(m,m),p>>=1;}
    return ret;
}int main(){
    int n,m,t,k,ni;
    scanf("%d%d%d%d",&n,&m,&t,&k);
    Matrix m1;
    memset(m1.a,0,sizeof(m1.a));
    m1.a[1][1]=1;
    for(int i=1;i<=t;i=ni+1){
        ni=t/(t/i);
        Matrix m2;
        double val=1.0*(t/i)/t;
        for(int j=1;j<=m-n;j++) m2.a[j][j]=1-val,m2.a[j+1][j]=val;
        m2.a[m-n+1][m-n+1]=1;
        for(int j=1;j<=m-n+1;j++) m2.a[m-n+2][j]=n+j-1;
        m2.a[m-n+2][m-n+2]=1;
        m1=Mul(ksm(m2,ni-i+1),m1);
    }printf("%.8lf\n",m1.a[m-n+2][1]*k);
}

C++11(clang++ 3.9) 解法, 执行用时: 970ms, 内存消耗: 608K, 提交时间: 2019-10-04 14:17:46

#include<bits/stdc++.h>
using namespace std;
bool s1;
int n,m,T,k;
double p;
struct ll{
	double A[22][22];
	void Init(){
		memset(A,0,sizeof(A));
		for(int i=n;i<=m;i++){
			if(i==m)A[i][i]=1;
			else {
				A[i][i]=1-p;
				A[i+1][i]=p;
			}
			A[m+1][i]=1.0*i;
		}
		A[m+1][m+1]=1;
	}
	ll operator *(const ll &a){
		ll res;
		memset(res.A,0,sizeof(res.A));
		for(int i=n;i<=m+1;i++)for(int j=n;j<=m+1;j++)for(int k=n;k<=m+1;k++)res.A[i][j]+=A[i][k]*a.A[k][j];
		return res;
	}
}Ans,num;
void fast(int x){
	num.Init();
	while(x){
		if(x&1)Ans=Ans*num;
		num=num*num;
		x>>=1;
	}
}
bool s2;
int main(){
	int sum=0;
	scanf("%d%d%d%d",&n,&m,&T,&k);
	Ans.A[n][m+1]=1;
	for(int r=T;r;){
		int l=(T/(T/r+1))+1;//从l 到 r 概率相等 
		p=(T/r)*1.0/T;
		fast(r-l+1);
		r=l-1;
	}
	printf("%.6lf\n",Ans.A[n][n]*k);
	return 0;
}