【CF553E】Kyoya and Train（期望dp，SPFA，FFT）

最新推荐文章于 2021-11-18 14:26:17 发布

原创最新推荐文章于 2021-11-18 14:26:17 发布 · 1k 阅读

1 ·

本内容遵循CC 4.0 BY-SA版权协议

多项式/生成函数同时被 2 个专栏收录

38 篇文章

订阅专栏

概率dp/期望dp

18 篇文章

订阅专栏

本文介绍了一种解决特定路径费用最小化问题的方法，通过倒序动态规划结合快速傅里叶变换（FFT）来降低时间复杂度。文章详细阐述了状态定义、边界条件及状态转移方程，并给出了具体实现代码。

考虑 dp。

发现正着 dp 好像不太好做，毕竟初值不太好设，而且时间一大于 $t$ 费用就要加上 $x$ ，所以考虑倒着 dp。

设 $f_{u,j}$ 表示现在已经到达 $u$ 点，耗时 $j$ ，问到达 $n$ 点的最小期望费用。

容易得到边界条件： $f_{i,j}=dis_i+x(j>t)$ （其中 $dis_i$ 表示 $i$ 到 $n$ 的最短路，可以用dijkstra 预处理），这个很好解释，一旦时间大于 $t$ 了，那么肯定要罚款，所以我们不用再关注时间这一状态了。

同时也容易得到状态转移方程：
$f_{u,j}=\min _{(u,v)}\left(w_{(u,v)}+\sum_{k=1}^{t}p_{(u,v),k}\times f_{v,j+k}\right)$
用 SPFA 转移的话，时间复杂度是 $O(nmt^2)$ ，不能接受

令 $i = (u, v)$ ，设 $gi,j=∑k=1tpi,k×fv,j+kg_{i,j}=\sum\limits_{k=1}^{t}p_{i,k}\times f_{v,j+k}$ 。

发现 $p$ 和 $f$ 的第二维相减是个和 $k$ 无关的值 $j$ ，所以考虑把它化成卷积形式。

设 $p^i=pt−k+1\hat{p}_i=p_{t-k+1}$ ，那么 $gi,j=∑k=1tp^i,t−k+1×fv,j+kg_{i,j}=\sum\limits_{k=1}^{t}\hat{p}_{i,t-k+1}\times f_{v,j+k}$ ， $p$ 和 $f$ 的第二维相加是定值 $t + j + 1$ ，是个卷积的形式。

不妨设 $A=p^iA=\hat{p}_i$ ， $B=f_v$ ， $C=g_i$ ，那么 $C = A B$ 。

所以我们不必枚举 $j$ 和 $k$ ，而是对于每一个 $i$ ，一遍 FFT 就能求出所有的 $g_{i,j}$ 。

那么时间复杂度就成功降到了 $O(nmtlog⁡t)O(nmt\log t)$ 。

~~至于更优秀的分治 FFT 做法蒟蒻表示并不会。~~

代码如下：

#include<bits/stdc++.h>

#define N 55
#define M 110
#define T 20010
#define INF 0x7fffffff

using namespace std;

const double pi=acos(-1);

struct Complex
{
	double x,y;
	Complex(){x=y=0;};
	Complex(double a,double b){x=a,y=b;}
}A[T<<3],B[T<<3],C[T<<3];

Complex operator + (Complex a,Complex b){return Complex(a.x+b.x,a.y+b.y);}
Complex operator - (Complex a,Complex b){return Complex(a.x-b.x,a.y-b.y);}
Complex operator * (Complex a,Complex b){return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}

struct data
{
	int u,s;
	data(){};
	data(int a,int b){u=a,s=b;}
	bool operator < (const data &a) const
	{
		return s>a.s;
	}
};

int n,m,t,x;
int cnt,head[N],nxt[M],to[M],w[M];
int dis[N];
bool vis[N],inq[N];
double p[M][T<<3],f[N][T<<3],now[T<<3];
int bit=0,limit=1;
int rev[T<<3];

void FFT(Complex *a,int opt)
{
	for(int i=0;i<limit;i++)
		if(i<rev[i]) swap(a[i],a[rev[i]]);
	for(int mid=1;mid<limit;mid<<=1)
	{
		int len=mid<<1;
		Complex wn=Complex(cos(pi/mid),opt*sin(pi/mid));
		for(int i=0;i<limit;i+=len)
		{
			Complex omega=Complex(1,0);
			for(int j=0;j<mid;j++,omega=omega*wn)
			{
				Complex x=a[i+j],y=omega*a[i+mid+j];
				a[i+j]=x+y,a[i+mid+j]=x-y;
			}
		}
	}
	if(opt==-1)
		for(int i=0;i<limit;i++)
			a[i].x/=limit;
}

void mul(double *a,double *b,double *c)
{
	for(int i=0;i<=t*2;i++) A[i]=Complex(a[i],0);
	for(int i=0;i<=t*2;i++) B[i]=Complex(b[i],0);
	for(int i=t*2+1;i<limit;i++) A[i]=B[i]=Complex(0,0);
	FFT(A,1),FFT(B,1);
	for(int i=0;i<limit;i++) C[i]=A[i]*B[i];
	FFT(C,-1);
	for(int i=0;i<limit;i++) c[i]=C[i].x;
}

void adde(int u,int v,int wi)
{
	to[++cnt]=v;
	w[cnt]=wi;
	nxt[cnt]=head[u];
	head[u]=cnt;
}

void dijkstra()
{
	priority_queue<data>q;
	memset(dis,127,sizeof(dis));
	dis[n]=0;
	q.push(data(n,dis[n]));
	while(!q.empty())
	{
		data now=q.top();
		q.pop();
		if(vis[now.u]) continue;
		vis[now.u]=1;
		for(int i=head[now.u];i;i=nxt[i])
		{
			int v=to[i];
			if(dis[now.u]+w[i]<dis[v])
			{
				dis[v]=dis[now.u]+w[i];
				q.push(data(v,dis[v]));
			}
		}
	}
}

void SPFA()
{
	queue<int>q;
	for(int i=1;i<=n;i++)
	{
		for(int j=0;j<=t;j++) f[i][j]=INF;
		for(int j=t+1;j<=2*t;j++) f[i][j]=dis[i]+x;
	}
	for(int j=0;j<=t;j++) f[n][j]=0;
	q.push(n);
	inq[n]=1;
	while(!q.empty())
	{
		int u=q.front();
		q.pop();
		inq[u]=0;
		for(int i=head[u];i;i=nxt[i])
		{
			int v=to[i];
			mul(f[u],p[i],now);
			bool flag=false;
			for(int j=0;j<=t;j++)
			{
				if(now[t+j+1]+w[i]<f[v][j])
				{
					f[v][j]=now[t+j+1]+w[i];
					flag=1;
				}
			}
			if(flag)
			{
				if(!inq[v])
				{
					inq[v]=1;
					q.push(v);
				}
			}
		}
	}
}

int main()
{
	scanf("%d%d%d%d",&n,&m,&t,&x);
	for(int i=1;i<=m;i++)
	{
		int u,v,w;
		scanf("%d%d%d",&u,&v,&w);
		adde(v,u,w);
		for(int j=t;j>=1;j--)
		{
			scanf("%lf",&p[i][j]);
			p[i][j]/=1e5;
		}
	}
	while(limit<=(t<<2)) limit<<=1,bit++;
	for(int i=0;i<limit;i++)
		rev[i]=(rev[i>>1]>>1)|((i&1)<<(bit-1));
	dijkstra();
	SPFA();
	printf("%.10lf\n",f[1][0]);
	return 0;
}