[CF1067D]Computer Game[凸包/斜率优化+倍增+矩阵乘法]

题意

你有 (n) 个任务，初始收益为 (a) ,共 (t) 轮游戏，每轮可以选择完成一个任务(可以做多次)，完成之后可以给任意任务升级，升级之后的任务收益为 (b) ，每个任务还有完成的概率 (p) ，问期望最大收益是多少。

(nleq 10^5,1leq a< bleq 10^8,tleq 10^9)

分析

• 一个重要而显然的结论是如果我们有了一次升级的机会，一定升级 (b*p) 最大的那一个，之后一直选择完成这个任务,记 (M) 表示最大的 (b*p) 。然后只关注没有一次完成时的最大期望收益。

• 定义状态 (f_t​) 表示还剩 (t​) 轮游戏，还没有完成一个任务的期望最大收益。

• 转移： (f_{t+1}=max{p_i*(tM+a_i)+(1-p_i)*f_t}).

• 可以变形成斜率优化的形式(注意这里的 (f_t) 要看成 (k) 的一部分) ,也可以写成这种形式：

[f_{t+1}=max{p_i*(tM-f_t)+p_ia_i}+f_t ]

• 将 (p_i) 看成 (k) , (p_ia_i) 看成 (b) ，然后每个转移都是形如直线 (y=kx+b) 的形式。

• 我们可以维护一个下凸壳，每次可以二分当前的 (x=tM-f_t) 在哪一段，可以做到 (Tlogn) 的时间。

• 发现 (tM-f_t)是单调不降的，证明如下：

设 (x_{t+1} geq x_t)
有：(tM-f_tgeq (t-1)M-f_{t-1})
(Mgeq f_t-f_{t-1})
考虑 (f​) 的实际意义可以发现上式一定成立，因为一轮的期望收益不会超过 (M) 。
可以考虑构造，对于 (t+1) 轮的最优决策， (t) 轮可以复制过来(除开第一轮)，而在 (t+1) 轮产生的每种情况，(t) 轮都可以效仿(发生的概率是相同的)，发现期望收益只有最后一轮不一样,最多差为 (M)。

• 然后在每一段倍增矩乘即可。

• 总时间复杂度为 (O(nlog T))。

代码

#include<bits/stdc++.h>
using namespace std;
#define go(u) for(int i=head[u],v=e[i].to;i;i=e[i].lst,v=e[i].to)
#define rep(i,a,b) for(int i=a;i<=b;++i)
#define pb push_back
typedef long long LL;
inline int gi(){
	int x=0,f=1;char ch=getchar();
	while(!isdigit(ch))	{if(ch=='-') f=-1;ch=getchar();}
	while(isdigit(ch)){x=(x<<3)+(x<<1)+ch-48;ch=getchar();}
	return x*f;
}
template<typename T>inline bool Max(T &a,T b){return a<b?a=b,1:0;}
template<typename T>inline bool Min(T &a,T b){return b<a?a=b,1:0;}
const int N=1e5 + 7;
typedef double db;
const db eps=1e-12;
int n,len,tp,st[N];
LL t;
db a[N],b[N],p[N],M;
int dcmp(db x){
	if(fabs(x)<eps) return 0;
	return x<0?-1:1;
}
struct line{
	db k,b;int id;
	line(){}line(db k,db b,int id):k(k),b(b),id(id){}
	bool operator <(const line &rhs)const{
		if(dcmp(k-rhs.k)!=0) return  dcmp(k-rhs.k)<0;
		return dcmp(b-rhs.b)<0;
	}
}l[N];
db X(int a,int b){
	return (l[b].b-l[a].b)/(l[a].k-l[b].k);
}
struct mat{
	db v[5][5];
	mat(){memset(v,0,sizeof v);}
	void init(){memset(v,0,sizeof v);}
	mat operator *(const mat &rhs)const{
		mat res;
		rep(i,0,4)
		rep(j,0,4)
		rep(k,0,4)
		res.v[i][j]+=v[i][k]*rhs.v[k][j];
		return res;
	}
}A,B[34],C;
mat Pow(mat a,LL b){
	mat res;
	rep(i,1,4) res.v[i][i]=1;
	for(;b;b>>=1,a=a*a) if(b&1) res=res*a;
	return res;
}
int main(){
	scanf("%d%I64d
",&n,&t);
	rep(i,1,n) {
		scanf("%lf%lf%lf",&a[i],&b[i],&p[i]);
		l[i]=line(p[i],p[i]*a[i],i);
		Max(M,b[i]*p[i]);
	}
	sort(l+1,l+1+n);
	rep(i,1,n) {
		if(dcmp(l[i].k-l[i+1].k)==0) continue;
		l[++len]=l[i];
	}
	rep(i,1,len){
		while(tp>1&&dcmp(X(st[tp-1],i)-X(st[tp],st[tp-1]))<=0) --tp;
		st[++tp]=i;
	}
	int now=1;LL cnt=0;
	for(int now=1;now<=tp&&cnt^t;++now){
		double R=cnt*M-A.v[1][1];
		while(now<tp&&X(st[now],st[now+1])<=R) ++now;
		
		int i=st[now];R=X(st[now],st[now+1]);
		B[0].init();
		B[0].v[1][1]=1-l[i].k;
		B[0].v[2][1]=l[i].k;
		B[0].v[2][2]=B[0].v[3][1]=B[0].v[3][3]=B[0].v[4][2]=B[0].v[4][4]=1;
		A.v[1][3]=l[i].b,A.v[1][4]=M;
		
		rep(i,1,33) B[i]=B[i-1]*B[i-1];
		
		for(int i=33;~i;--i)if(t-cnt>(1ll<<i)){
			C=A*B[i];
			if(now==tp||dcmp((cnt+(1ll<<i))*M-C.v[1][1]-R)<0) A=C,cnt+=(1ll<<i);
		}
		cnt++,A=A*B[0];
	}
	printf("%.13lf
",A.v[1][1]);
	return 0;
}