容易发现,当某一次游戏成功后,一定是一直选择 \(p_ib_i\) 最大的游戏玩。设 \(s=\max\limits_{i=1}^n p_ib_i\)
定义 \(dp_i\) 为还有 \(i\) 次操作时,最大的期望。
那么 \(dp_i=\max\limits_{j=1}^n(1-p_j)dp_{i-1}+p_j(a_j+(i-1)s)=dp_{i-1}+p_j((i-1)s-dp_{i-1})+a_j\)
求 \(dp_{i}\) 的答案时,我们发现要想办法找到这个 \(j\)。把第 \(j\) 个游戏看为平面直角坐标系上 \(p_j,p_j\times b_j\) 的一个点的话,你会发现所有可能的决策一定在这些点的凸包上。容易证明 \((i-1)s-dp_{i-1}\) 一定是单调递增的,所以有决策单调性,可以扫过去判断哪个是最优决策。复杂度 \(O(n+T)\)
然后考虑把 \(T\) 给省掉。由于很长的一段 \(dp\) 都是一样的操作,这个操作可以用矩阵快速幂加速,所以考虑上矩阵乘法。倍增结束位置(也可以推数学),然后矩阵乘法计算判断有没有超过斜率的范围,复杂度降到 \(O(nlogT)\)
注意卡精度,把斜率的除法要全部换成乘法。
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const int N=1e5+5;
const double eps=1e-12;
struct matrix{
double a[3][3];
matrix operator*(const matrix&m)const{
matrix s;
memset(s.a,0,sizeof(s.a));
for(int i=0;i<3;i++)
for(int j=0;j<3;j++)
for(int k=0;k<3;k++)
s.a[i][k]+=a[i][j]*m.a[j][k];
return s;
}
void operator=(const matrix&m){
memcpy(a,m.a,sizeof(a));
}
}pw[50],h,k;
struct line{
double k,b;
int id;
}g[N];
int n,st[N],tp=1,a[N];
LL t,nw;
double s,p[N],dp;
int cmp1(line x,line y)
{
if(fabs(x.k-y.k)>eps)
return x.k<y.k;
return x.b>y.b;
}
double sl(int x,int y)
{
return (g[y].b-g[x].b)/(g[x].k-g[y].k);
}
int main()
{
scanf("%d%lld",&n,&t);
for(int i=1,b;i<=n;i++)
scanf("%d%d%lf",a+i,&b,p+i),s=max(s,p[i]*b);
for(int i=1;i<=n;i++)
g[i]=(line){p[i],p[i]*a[i],i};
sort(g+1,g+n+1,cmp1);
st[1]=1;
for(int i=2;i<=n;i++)
{
if(fabs(g[i].k-g[st[tp]].k)<=eps)
continue;
while(tp^1&&(g[i].b-g[st[tp]].b)*(g[st[tp-1]].k-g[st[tp]].k)<(g[st[tp]].b-g[st[tp-1]].b)*(g[st[tp]].k-g[i].k))
--tp;
st[++tp]=i;
}
pw[0].a[1][1]=1;
pw[0].a[2][1]=1;
pw[0].a[2][2]=1;
h.a[0][2]=1;
for(int i=1;i<=tp;i++)
{
while(i^tp&&(nw*s-h.a[0][0])*(g[st[i+1]].k-g[st[i]].k)>g[st[i]].b-g[st[i+1]].b)
++i;
pw[0].a[0][0]=1-p[g[st[i]].id];
pw[0].a[1][0]=p[g[st[i]].id]*s;
pw[0].a[2][0]=p[g[st[i]].id]*a[g[st[i]].id];
for(int i=1;i<=40;i++)
pw[i]=pw[i-1]*pw[i-1];
for(int j=40;~j;--j)
{
if(nw+(1LL<<j)<t&&(i==tp||((nw+(1LL<<j))*s-(h*pw[j]).a[0][0])*(g[st[i+1]].k-g[st[i]].k)<=g[st[i]].b-g[st[i+1]].b))
{
nw+=1LL<<j;
h=h*pw[j];
}
}
h=h*pw[0],++nw;
if(nw==t)
break;
}
printf("%.6lf",h.a[0][0]);
}