题意
给一张有向图,每条边的通过时间是一个1~t的给定离散概率分布。每条边还有一个通过费用。你现在要从1号点到n号点,假如你到达n号点的时间比t大,那么你需要额外付出X的费用。请求出最优策略下的最小期望费用。
n
≤
50
,
m
≤
100
,
t
≤
20000
n\leq 50,m\leq 100,t\leq 20000
n≤50,m≤100,t≤20000
分析
- 先搞清楚最优策略是什么意思:
- 设 f [ t ] [ x ] f[t][x] f[t][x]表示时刻t在x时,到达n的最小期望费用。
- 根据定义,则有 f [ t ] [ x ] = m i n ( c e + ∑ a f [ t + a ] [ y ] × P a ) f[t][x]=min(c_e+\sum _af[t+a][y]\times P_a) f[t][x]=min(ce+∑af[t+a][y]×Pa)
- f [ t ′ ] [ x ] , t ′ > t f[t'][x],t'>t f[t′][x],t′>t的值:此时必定付出额外费用,因此是走花费最短路,可以通过预处理得出。
- 这个dp可以通过分治fft来优化。 对于t’>t的情况我是先通过预处理算好贡献的。
- 具体来说,对时间分治:现在要计算 t ∈ [ m i d + 1 , R ] t\in[mid+1,R] t∈[mid+1,R]对 t ∈ [ L , m i d ] t\in[L,mid] t∈[L,mid]的状态的贡献。枚举每一条边,可以发现这个转移就是一个减法卷积。做一个FFT就可以了。
- 实现上,要注意将次数界控制在32768,否则有被卡常的风险。
#include <bits/stdc++.h>
using namespace std;
typedef double db;
const int N = 3.5e4;
const db E = 1e5, inf = 1e9, pi = acos(-1);
const int SZ = sizeof E;
int n,m,t,W;
struct edge{
int from, to, c;
db p[N];
} e[110];
int d[55][55];
db upd[20010][110];
db f[20010][55];
db A[N], B[N], C[N];
complex<db> w[N];
int M, h[N], MM = 32768;
void fft(complex<db> *a) {
for(int i = 0; i < M; i++)
if (h[i] > i) swap(a[i], a[h[i]]);
int st = MM;
for(int m = 1; m < M; m <<= 1) {
st >>= 1;
for(int i = 0; i < M; i += (m << 1)) {
for(int j = 0; j < m; j++) {
complex<db> z = a[i + j + m] * w[j * st];
a[i + j + m] = a[i + j] - z;
a[i + j] = a[i + j] + z;
}
}
}
}
void ifft(complex<db> *a) {
reverse(a + 1, a + M);
fft(a);
db z = 1.0 / M;
for(int i = 0; i < M; i++) a[i].real(a[i].real() * z);
}
int mx;
void calc(int sz,int cc) {
int sa = sz - cc;
memset(C, 0, (sz + 5) * SZ);
for(M = 1; M <= sa + sz; M <<= 1);
static complex<db> a[N], b[N];
for(int i = 1; i < M; i++)
h[i] = (h[i >> 1] >> 1) + (i & 1) * (M >> 1);
memset(a, 0, M * sizeof a[0]);
memset(b, 0, M * sizeof b[0]);
for(int i = 0; i <= sz; i++) {
a[i].real(A[i]);
b[i].real(B[sz - i]);
}
fft(a);
fft(b);
for(int i = 0; i < M; i++) a[i] = a[i] * b[i];
ifft(a);
for(int i = 0; i < sz; i++) {
C[i] = a[i + sz - cc].real();
}
}
void divide(int L, int R) {
if (L > R) return;
if (L == R) {
for(int i = 1; i <= n; i++) f[L][i] = inf;
for(int i = 1; i <= m; i++) {
f[L][e[i].from] = min(f[L][e[i].from], upd[L][i] + e[i].c);
}
f[L][n] = 0;
return;
}
int mid = L + R >> 1, sz = R - L + 1;
divide(mid + 1, R);
for(int i = 1; i <= m; i++) {
int x = e[i].from, y = e[i].to;
memset(A, 0, (sz + 5) * SZ);
memset(B, 0, (sz + 5) * SZ);
int cc = mid - L + 1;
for(int j = mid + 1; j <= R; j++) {
A[j - L - cc] = f[j][y];
}
for(int j = 1; j <= R - L + 1; j++) {
B[j] = e[i].p[j];
}
calc(sz, cc);
for(int j = L; j <= mid; j++) {
upd[j][i] += C[j - L];
}
}
divide(L, mid);
}
void read(int &a) {
char c; while((c=getchar()) > '9' || c < '0');
a = c - '0'; while ((c=getchar()) >= '0' && c <= '9') a = a * 10 + c - '0';
}
int main() {
int g = clock();
freopen("e.in","r",stdin);
cin >> n >> m >> t >> W;
memset(d,127/2,sizeof d);
db WW = 1.0 / E;
for(int i = 1; i <= m; i++) {
scanf("%d %d %d", &e[i].from, &e[i].to, &e[i].c);
d[e[i].from][e[i].to] = min(d[e[i].from][e[i].to], e[i].c);
for(int j = 1; j <= t; j++) {
int z; read(z);
e[i].p[j] = z * WW;
}
}
for(int i = 1; i <= n; i++) d[i][i] = 0;
for(int k = 1; k <= n; k++) {
for(int i = 1; i <= n; i++) if (i != k) {
for(int j = 1; j <= n; j++) if (j != k && j != i) {
d[i][j] = min(d[i][j], d[i][k] + d[k][j]);
}
}
}
for(int i = 1; i <= m; i++) {
db po = 0;
for(int j = 1; j <= t; j++) {
po += e[i].p[t - j + 1];
upd[j][i] = (d[e[i].to][n] + W) * po;
if (upd[j][i] > inf) upd[j][i] = inf;
}
}
w[0].real(1);
for(int i = 1; i < MM; i++) {
w[i] = complex<db>(cos(i*2*pi/MM),sin(i*2*pi/MM));
}
divide(0, t);
printf("%.10lf\n",f[0][1]);
cerr << clock() - g << endl;
}