Link
设(f_{u,i})表示(i)时刻到(u)的最小答案,那么我们有:(f_{u,i}=minlimits_{(u,v,w,id)in E}(w+sumlimits_{j=0}^tf_{v,,i+j}p_{id,j}))。
令(g_{e,i})表示(i)时刻走上(e)这条边的最小答案,那么我们有:(g_{e,i}=w_e+sumlimits_{j=1}^tf_{v_e,i+j}p_{e,j},f_{u,i}=minlimits_{u_e=u}g_{e,i})。
那么我们可以把(p)翻转一下,右边的就变成了一个卷积的形式。
注意到只有后面的时间会对前面的时间有影响,所以我们对时间分治,先计算右边,然后计算右边对左边的影响,再计算左边即可。
边界的计算需要预处理一些最短路、前缀和之类的东西。
#include<cmath>
#include<cstdio>
#include<cctype>
#include<cstring>
#include<algorithm>
using ld=double;
const int N=20007,M=107,S=57,T=65537;const ld pi=acos(-1);
namespace IO
{
char ibuf[(1<<21)+1],*iS,*iT;
char Get(){return (iS==iT? (iT=(iS=ibuf)+fread(ibuf,1,(1<<21)+1,stdin),(iS==iT? EOF:*iS++)):*iS++);}
int read(){int x=0,c=Get();while(!isdigit(c))c=Get();while(isdigit(c))x=x*10+c-48,c=Get();return x;}
}
using IO::read;
struct complex{ld x,y;}A[T],B[T],w[T];
complex operator+(const complex&a,const complex&b){return {a.x+b.x,a.y+b.y};}
complex operator-(const complex&a,const complex&b){return {a.x-b.x,a.y-b.y};}
complex operator*(const complex&a,const complex&b){return {a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};}
struct edge{int u,v,w;ld p[N],g[N],s[N];}e[M];
int n,m,t,x,lim,vis[S],dis[S],rev[T];ld f[S][N];
void init(int n)
{
int len=32-__builtin_clz(n);lim=1<<len;
for(int i=1;i<lim;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<(len-1));
for(int i=0;i<lim;++i) {ld t=pi*i/lim;w[i]={cos(t),sin(t)};}
}
void FFT(complex*a,int f)
{
static complex x;
if(!~f) std::reverse(a+1,a+lim);
for(int i=0;i<lim;++i) if(i<rev[i]) std::swap(a[i],a[rev[i]]);
for(int i=1;i<lim;i<<=1) for(int l=i<<1,j=0;j<lim;j+=l) for(int k=0;k<i;++k) x=a[i+j+k]*w[lim/i*k],a[i+j+k]=a[j+k]-x,a[j+k]=a[j+k]+x;
if(!~f) for(int i=0;i<lim;++i) a[i].x/=lim,a[i].y/=lim;
}
void solve(int l,int r)
{
if(l==r)
{
for(int i=1;i<n;++i) f[i][l]=1e9;
for(int i=1;i<=m;++i) e[i].g[l]+=e[i].s[t-l+1]*(dis[e[i].v]+x),f[e[i].u][l]=std::min(f[e[i].u][l],e[i].g[l]+e[i].w);
return ;
}
int mid=(l+r)>>1;
solve(mid+1,r),init(r-l+r-mid);
for(int j=1;j<=m;++j)
{
memset(A,0,lim<<5),memset(B,0,lim<<5);
for(int i=1;i<=r-l;++i) A[r-l-i].x=e[j].p[i];
for(int i=mid+1;i<=r;++i) B[i-mid-1].x=f[e[j].v][i];
FFT(A,1),FFT(B,1);
for(int i=0;i<lim;++i) A[i]=A[i]*B[i];
FFT(A,-1);
for(int i=0;i<mid-l+1;++i) e[j].g[l+i]+=A[r-mid-1+i].x;
}
solve(l,mid);
}
int main()
{
n=read(),m=read(),t=read(),x=read(),memset(dis,0x3f,sizeof dis),dis[n]=0;
for(int i=1;i<=m;++i)
{
e[i].u=read(),e[i].v=read(),e[i].w=read();
for(int j=1;j<=t;++j) e[i].p[j]=read()*0.00001;
for(int j=t;j;--j) e[i].s[j]=e[i].s[j+1]+e[i].p[j];
}
for(int i=1,u;i<=n;++i)
{
u=0;
for(int j=1;j<=n;++j) if(!vis[j]&&dis[j]<dis[u]) u=j;
vis[u]=1;
for(int j=1;j<=m;++j) if(e[j].v==u) dis[e[j].u]=std::min(dis[e[j].u],dis[u]+e[j].w);
}
solve(0,t),printf("%.10lf",f[1][0]);
}