【CF553E】Kyoya and Train（期望dp，SPFA，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\) 了，那么肯定要罚款，所以我们不用再关注时间这一状态了。

同时也容易得到状态转移方程：

用 SPFA 转移的话，时间复杂度是 \(O(nmt^2)\)，不能接受

令 \(i=(u,v)\)，设 \(g_{i,j}=\sum\limits_{k=1}^{t}p_{i,k}\times f_{v,j+k}\)。

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

设 \(\hat{p}_i=p_{t-k+1}\)，那么 \(g_{i,j}=\sum\limits_{k=1}^{t}\hat{p}_{i,t-k+1}\times f_{v,j+k}\)，\(p\) 和 \(f\) 的第二维相加是定值 \(t+j+1\)，是个卷积的形式。

不妨设 \(A=\hat{p}_i\)，\(B=f_v\)，\(C=g_i\)，那么 \(C=AB\)。

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

那么时间复杂度就成功降到了 \(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;
}

原文地址：http://www.cnblogs.com/ez-lcw/p/16837309.html