考虑 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)\),设 \(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;
}
标签:Kyoya,int,double,FFT,CF553E,now,dp,dis
From: https://www.cnblogs.com/ez-lcw/p/16837309.html