【luogu P7418】Counting Graphs P（DP）（思维）（容斥）

Counting Graphs P

题目链接：luogu P7418

题目大意

给你一个图，然后 fi,j 表示是否存在一个从 1 到 i 的路径经过的边数是 j。
然后问你能构造出来多少个图使得它的 f 函数跟给出的图的一样。
会有自环，没有重边。

思路

首先参考 Minimizing Edges P 的思路，搞出奇偶最短路。

首先还是一样，特判二分图。
那如果二分图，每一层的点数分别是 $a_1,a_2,...$。
那我们考虑对于每两层之间算连边的方案数，然后乘起来。
那对于第 $i$ 层和第 $i+1$ 层，然后其实要 $i+1$ 每个都连到，那其实对于每个点可以枚举令一行与它的连边状态，那就是 $(2^{a_i}-1{)}^{a_{i+1}}$。

然后我们考虑继续用上一题的思路，$x+y$ 为 $y$ 坐标，$x$ 为 $x$ 坐标。
然后考虑在这个上面 DP。
那往上一层的时候还好，但一层内转移是两边都要的，所以我们要多设一维，表示 $f_{x,y,p}$ 为当前到 $(x,y)$ 这个状态，$p$ 个点需要下一个状态连过去（这个转移要是内部的转移）。（$s1$ 为 $(x,y)$ 这个状态的个数，$s2$ 为 $(x-1,y-1)$ 这个状态的个数）

$f_{x,y,p}=\sum _{q=0}^{s1-p}{g}_{x,y,q}{C}_{s1-p}^q(2^{s2}-1{)}^{s1-p}$

解释：
我们枚举的 $q$ 是在这里两层传递 $q$ 个。
$C_{s1-p}^q$：在剩下不内部转移的点中选 $q$ 个进行一层的转移。
$(2^{s2}-1{)}^{s1-p}$：两层之间转移的方案数，跟二分图的计算是一样的。
$g_{x,y,q}$：表示 $s1-q$ 个点都是由 $(x-1,y+1)$ 转移过来的方案数。

然后接着你考虑如何转移 $g_{x,y,p}$。
考虑枚举上一个状态有多少个点需要现在这个状态把边连过去，然后枚举为 $q$ 个，然后：（$s1$ 为 $(x,y)$ 这个状态的个数，$s3$ 为 $(x+1,y-1)$ 这个状态的个数）

$g_{x,y,p}=C_{s1}^p\sum _{q=0}^{s3}{f}_{x+1,y-1,q}{h}_{s3,q,s1-p}$

解释：
$C_{s1}^p$：先要选 $p$ 个。
$f_{x+1,y-1,q}$：上一个状态的方案数。
$h_{s3,q,s1-p}$：表示大小为 $s3$ 的集合跟大小为 $s1-p$ 的集合之间连边，保证 $s2$ 中的 $q$ 个点和大小为 $s1-p$ 的集合的度数至少是 $1$，它的方案数。

不难看出又要求 $h_{i,j,k}$，考虑用容斥，枚举 $p$ 个 $j$ 中的点没有出度。

$h_{i,j,k}=\sum _{p=0}^j(-1{)}^p{C}_j^p(2^{i-p}-1{)}^k$

解释：
$(-1{)}^p$：容斥
$C_j^p$：在要的 $j$ 个点中选 $p$ 个不要。
$(2^{i-p}-1{)}^k$：匹配的方案数，跟二分图的搞法是一样的。

然后考虑如何统计答案，不难看出就是每一行的最后一个（$x+1=y$ 或者就是最后一个状态）
那明显的看到这两个肯定统计方法是不一样的（一个可以内部消化，一个不行）
那如果不能内部消化，那就不能往右边连，所以贡献就是 $f_{x,y,0}$。

那如果是 $x+1=y$ 呢？
那我们就可以同类消化，使得任意个点要往右连都是可以的。
设 $(x,y)$，即 $(x,x+1)$ 的状态个数是 $s$。

那贡献就是 $\sum _{i=0}^s{T}_{s,i}{f}_{x,y,i}$。

解释：
$f_{x,y,i}$：就是原来状态的方案数。
$T_{s,i}$：大小为 $s$ 的集合中 $i$ 个点通过同类连边消化掉的方案数。

然后考虑求 $T_{i,j}$，也是容斥，设 $k$ 个点的度数是 $0$。

$T_{i,j}=\sum _{k=0}^j(-1{)}^k{C}_j^k{2}^{\frac{(i-k)(i-k+1)}{2}}$

解释：
$(-1{)}^k$：容斥
$C_j^k$：选 $k$ 个度数为 $0$
$2^{\frac{(i-k)(i-k+1)}{2}}$：除去 $k$ 个点剩下点任意连边的方案数。
就是对于每一条可能的边都有选或者不选，点数是 $i-k=x$，可能的边就是 $\frac{x(x-1)}{2}$ 即 $\frac{(i-k)(i-k-1)}{2}$，然后每条边都可以选或不选就是 $2$ 的那么多次方了。

然后自此，你就可以搞出来了！！！

至于 $T,h$ 你可以直接预处理，也可以开个记忆化，要用再算，然后 $f,g$ 就是每次询问都会变，记得初始化。
然后你可以预处理出来阶乘（以及它的逆元），组合数，$2^x$（$x⩽n^2$），以及 $(2^x-1{)}^y$（$x,y⩽n$）

然后就可以搞啦！

代码

#include<map>
#include<queue>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define ll long long
#define mo 1000000007
#define INF 0x3f3f3f3f3f3f3f3f

using namespace std;

const int N = 205;
const int M = 80005;
struct node {
	int to, nxt;
}e[M];
int T, n, m, x, y, le[N], KK;
int cunt[N];
ll ans, f[N][N][N >> 1], g[N][N][N >> 1], t[N][N], h[N][N][N];
ll jc[N], cff[N][N], inv[N], cf[M], s1, s2, s3, re;
bool gogo;

struct st {
	int x, y;
}tp[N];
int nm[N][N];
int tpn, tpx, tpy;

struct ztzt {
	int dis, now;
};
bool operator <(ztzt x, ztzt y) {
	return x.dis > y.dis;
}
priority_queue <ztzt> q;
int dis[N << 1];
bool in[N << 1];

void add(int x, int y) {
	e[++KK] = (node){y, le[x]}; le[x] = KK;
	e[++KK] = (node){x, le[y]}; le[y] = KK;
}

void dij() {//跑出奇偶最短路
	for (int i = 0; i <= 2 * n; i++)
		in[i] = 0, dis[i] = INF;
	dis[1] = 0; q.push((ztzt){0, 1});
	while (!q.empty()) {
		int now = q.top().now;
		q.pop();
		if (in[now]) continue;
		in[now] = 1;
		for (int i = le[now]; i; i = e[i].nxt)
			if (!in[e[i].to] && dis[e[i].to] > dis[now] + 1) {
				dis[e[i].to] = dis[now] + 1;
				q.push((ztzt){dis[e[i].to], e[i].to});
			}
	}
}

bool cmp1(st x, st y) {
	if (x.x + x.y != y.x + y.y) return x.x + x.y < y.x + y.y;
	return x.x < y.x;
}

ll ksm(ll x, ll y) {
	ll re = 1;
	while (y) {
		if (y & 1) re = (re * x) % mo;
		x = (x * x) % mo;
		y >>= 1;
	}
	return re;
}

ll C(int n, int m) {
	if (n < m) return 0ll;
	return jc[n] * inv[m] % mo * inv[n - m] % mo; 
}

ll H(int i, int j, int k) {
	if (h[i][j][k] != -1) return h[i][j][k];
	h[i][j][k] = 0;
	for (int p = 0; p <= j; p++) {
		h[i][j][k] = (h[i][j][k] + ((p & 1) ? -1 : 1) * C(j, p) * cff[i - p][k] % mo) % mo;
		if (h[i][j][k] < 0) h[i][j][k] += mo;
	}
	return h[i][j][k];
}

int main() {	
	jc[0] = 1;//预处理
	for (int i = 1; i < N; i++)
		jc[i] = (jc[i - 1] * i) % mo;
	inv[N - 1] = ksm(jc[N - 1], mo - 2);
	for (int i = N - 2; i >= 0; i--)
		inv[i] = (inv[i + 1] * (i + 1)) % mo;
	cf[0] = 1;
	for (int i = 1; i < M; i++) cf[i] = (cf[i - 1] << 1) % mo;
	for (int i = 0; i < N; i++) {
		cff[i][0] = 1;
		for (int j = 1; j < N; j++)
			cff[i][j] = (cff[i][j - 1] * (cf[i] - 1)) % mo;
	}
	
	memset(h, -1, sizeof(h));
	
	for (int i = 0; i < N; i++)
		for (int j = 0; j <= i; j++) {
			for (int k = 0; k <= j; k++) {
				t[i][j] = (t[i][j] + ((k & 1) ? -1 : 1) * C(j, k) * cf[(i - k) * (i - k + 1) / 2] % mo) % mo;
				if (t[i][j] < 0) t[i][j] += mo;
			}
		}
	
	scanf("%d", &T);
	while (T--) {
		scanf("%d %d", &n, &m);
		
		for (int i = 0; i <= 2 * n; i++)//清空数组
			for (int j = 0; j <= 2 * n; j++)
				for (int k = 0; k <= n; k++)
					f[i][j][k] = g[i][j][k] = 0;
		KK = 0;
		for (int i = 0; i <= 2 * n; i++) le[i] = 0;
		
		for (int i = 1; i <= m; i++) {
			scanf("%d %d", &x, &y);
			add(x, y + n); add(x + n, y);
		}
		
		dij();
		
		gogo = 0;
		for (int i = 1; i <= n; i++) {
			if (dis[i] == dis[0] || dis[i + n] == dis[0]) {
				gogo = 1;
				break;
			}
		}
		ans = 1; tpn = 0;
		if (gogo) {//特判二分图
			for (int i = 0; i <= n * 2; i++) cunt[i] = 0;
			for (int i = 1; i <= n; i++) {
				if (dis[i] != dis[0]) cunt[dis[i]]++;
					else if (dis[i + n] != dis[0]) cunt[dis[i + n]]++;
			}
			for (int i = 1; i <= n * 2; i++) {
				ans = (ans * cff[cunt[i - 1]][cunt[i]]) % mo;
			}
			printf("%lld\n", ans);
		}
		else {
			for (int i = 1; i <= n; i++) {
				tpx = dis[i]; tpy = dis[i + n];
				if (tpx > tpy) swap(tpx, tpy);
				nm[tpx][tpy]++; tp[++tpn] = (st){tpx, tpy};
			}
			sort(tp + 1, tp + tpn + 1, cmp1);
			
			for (int i = 1; i <= tpn; i++) {//DP
				int st = i;
				while (i < tpn && tp[i].x == tp[i + 1].x && tp[i].y == tp[i + 1].y) i++;
				x = tp[i].x; y = tp[i].y;
				
				s1 = nm[x][y];
				if (!x || !y) s2 = 0;
					else s2 = nm[x - 1][y - 1];
				if (!x) s3 = 0;
					else s3 = nm[x - 1][y + 1];
				if (!s3) g[x][y][st == 1 ? 0 : s1] = 1;
				else {
					for (int p = 0; p <= s1; p++) {
						for (int q = 0; q <= s3; q++)
							g[x][y][p] = (g[x][y][p] + f[x - 1][y + 1][q] * H(s3, q, s1 - p) % mo) % mo;
						g[x][y][p] = (g[x][y][p] * C(s1, p)) % mo;
					}
				}
				for (int p = 0; p <= s1; p++) {
					for (int q = 0; q <= s1 - p; q++)
						f[x][y][p] = (f[x][y][p] + C(s1 - q, p) * g[x][y][q] % mo * cff[s2][s1 - p] % mo) % mo;
				}
				
				if (i == tpn || tp[i + 1].x != x + 1 || tp[i + 1].y != y - 1) {//把每一层的贡献算了
					if (x + 1 == y) {
						re = 0;
						for (int p = 0; p <= s1; p++)
							re = (re + f[x][y][p] * t[s1][p] % mo) % mo;
						ans = (ans * re) % mo;
					}
					else {
						ans = (ans * f[x][y][0]) % mo;
					}
				}
			}
			
			printf("%lld\n", ans);
			
			for (int i = 1; i <= n; i++) {
				tpx = dis[i]; tpy = dis[i + n];
				if (tpx > tpy) swap(tpx, tpy);
				nm[tpx][tpy] = 0;
			}
		}
	} 
	
	return 0;
}