最小费用最大流，但能过 HLPP 板子题，还能处理负环
namespace flow{ // }{{{
typedef long long ll;
constexpr int V = 5e3, E = 5e4;
constexpr int EDGE_NIL = 0;
struct Edge{
int to;
ll lf, cost;
int nxt;
} es[E*2+4];
ll sumcost = 0, sumflow = 0;
int is, it, iv;
ll minc, maxf;
int head[V], cnt = (EDGE_NIL|1)+1;
ll pi[V];
int fe[V], mark[V], time = 2; // fe: father edge
int fa[V]; // 200ms improve
void init(int v, int s, int t){
is = s, it = t, iv = v;
sd fill(head, head+iv, EDGE_NIL);
sd fill(mark, mark+iv, 0);
sd fill(fe, fe+iv, EDGE_NIL);
time = 2;
cnt = (EDGE_NIL|1)+1;
sumcost = sumflow = 0;
minc = maxf = 0;
}
void addflow(int s, int t, ll f, ll c){
es[cnt] = (Edge){t, f, c, head[s]}, head[s] = cnt++;
es[cnt] = (Edge){s, 0, -c, head[t]}, head[t] = cnt++;
sumflow += f, sumcost += sd abs(c);
}
void mktree(int x, int from_e){
fe[x] = from_e;
fa[x] = es[from_e^1].to;
mark[x] = 1;
for(int i=head[x]; i!=EDGE_NIL; i=es[i].nxt){
if(mark[es[i].to] == 1 || es[i].lf == 0) continue;
mktree(es[i].to, i);
}
}
ll getpi(int x){
if(mark[x] == time) return pi[x];
mark[x] = time;
pi[x] = getpi(fa[x]) - es[fe[x]].cost;
return pi[x];
}
ll pushflow(int e){ // return delta-cost
int rt = es[e].to, lca = es[e^1].to;
time++;
while(rt){
mark[rt] = time;
rt = fa[rt];
}
while(mark[lca] != time){
mark[lca] = time;
lca = fa[lca];
}
ll df = es[e].lf;
int todel = e, dir = -1; // dir: direction, 0->es[e].to
for(int i=es[e^1].to; i!=lca; i=fa[i]){
if(es[fe[i]].lf < df){
df = es[fe[i]].lf;
todel = fe[i];
dir = 1;
}
}
for(int i=es[e].to; i!=lca; i=fa[i]){
if(es[fe[i]^1].lf < df){
df = es[fe[i]^1].lf;
todel = fe[i];
dir = 0;
}
}
ll dcst = 0; // delta cost
if(df) {
for(int i=es[e].to; i!=lca; i=fa[i]){
es[fe[i]].lf += df;
es[fe[i]^1].lf -= df;
dcst += es[fe[i]^1].cost * df;
}
for(int i=es[e^1].to; i!=lca; i=fa[i]){
es[fe[i]].lf -= df;
es[fe[i]^1].lf += df;
dcst += es[fe[i]].cost * df;
}
es[e].lf -= df;
es[e^1].lf += df;
dcst += es[e].cost * df;
}
if(todel == e) return dcst;
int last = e^dir, lastu = es[e^dir^1].to;
for(int i=es[e^dir].to; i!=es[todel^1].to; ){
mark[i]=time-1;
int i_ = fa[i];
fa[i] = lastu;
lastu = i;
sd swap(fe[i], last);
last ^= 1;
i=i_;
}
return dcst;
}
void mcmf(){
ll sfl_ = sumflow, scs_ = sumcost;
addflow(it, is, sumflow+1, -(sumcost+1));
sumflow = sfl_, sumcost = scs_;
mktree(it, EDGE_NIL);
mark[it] = ++time;
fa[it] = 0;
bool run = true;
while(run){
run = false;
UP(i, (EDGE_NIL|1)+1, cnt){
int s = es[i^1].to, t = es[i].to;
if(es[i].lf && es[i].cost + getpi(t) - getpi(s) < 0){
run = true;
minc += pushflow(i);
}
}
}
maxf = es[cnt-1].lf;
minc += maxf * (scs_+1);
}
} // {}}}