问题重述

给定一个包含 $n$ 个顶点、初始无边的无向图。重复执行以下操作直到无法继续：

从所有满足 $\deg_u=0$ 或 $\deg_v=0$ 的点对 $(u,v)$ 中等概率选取一对，添加边 $(u,v)$。

给定整数 $m$，求恰好执行 $m$ 次操作后停止的概率，结果对质数 $P$ 取模。其中 $n,m$ 分别对应原题中的 $M$ 和 $M-K$。

小规模情况：$n \leq 5$

对于极小规模，直接采用深度优先搜索枚举所有可能的加边序列。维护每个顶点的度数，每次从合法边集中等概率选择，递归计算概率。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

int n, target, mod;
vector<int> degree;
ll result = 0;

ll power(ll base, ll exp) {
    ll ans = 1;
    while (exp) {
        if (exp & 1) ans = ans * base % mod;
        base = base * base % mod;
        exp >>= 1;
    }
    return ans;
}

void search(int step, vector<int>& deg, ll prob) {
    bool has_isolated = false;
    for (int i = 1; i <= n; i++) {
        if (deg[i] == 0) { has_isolated = true; break; }
    }
    if (!has_isolated) {
        if (step == target) result = (result + prob) % mod;
        return;
    }
    
    vector<pair<int,int>> candidates;
    for (int u = 1; u < n; u++) {
        for (int v = u + 1; v <= n; v++) {
            if (deg[u] == 0 || deg[v] == 0) {
                candidates.emplace_back(u, v);
            }
        }
    }
    if (candidates.empty()) return;
    
    ll inv = power(candidates.size(), mod - 2);
    for (auto& [u, v] : candidates) {
        deg[u]++; deg[v]++;
        search(step + 1, deg, prob * inv % mod);
        deg[u]--; deg[v]--;
    }
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> n >> target >> mod;
    degree.resize(n + 1, 0);
    search(0, degree, 1);
    cout << result % mod;
    return 0;
}

中等规模：$n \leq 4000$ 的动态规划解法

设 $dp[i][j]$ 表示执行 $i$ 次操作后，图中仍有 $j$ 个孤立点的概率。状态转移考虑两种加边方式：

• 孤立点-孤立点：减少 2 个孤立点，方案数 $\binom{j}{2}$
• 孤立点-非孤立点：减少 1 个孤立点，方案数 $j \times (n-j)$

设当前总边数为 $t = \frac{n(n-1)}{2} - \frac{(n-j)(n-j-1)}{2}$，则转移概率需乘以 $t$ 的逆元。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

const int MAXN = 4005;
ll dp[MAXN][MAXN];
int n, m;
ll mod;

ll power(ll a, ll b) {
    ll res = 1;
    while (b) {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> n >> m >> mod;
    
    dp[0][n] = 1;
    for (int i = 0; i < m; i++) {
        for (int j = n; j >= 1; j--) {
            if (!dp[i][j]) continue;
            ll edges = (1LL * n * (n - 1) / 2 - 1LL * (n - j) * (n - j - 1) / 2 + mod) % mod;
            ll inv = power(edges, mod - 2);
            if (j >= 2) {
                ll ways = 1LL * j * (j - 1) / 2 % mod;
                dp[i + 1][j - 2] = (dp[i + 1][j - 2] + ways * inv % mod * dp[i][j]) % mod;
            }
            if (j >= 1) {
                ll ways = 1LL * j * (n - j) % mod;
                dp[i + 1][j - 1] = (dp[i + 1][j - 1] + ways * inv % mod * dp[i][j]) % mod;
            }
        }
    }
    cout << (dp[m][0] + mod) % mod;
    return 0;
}

特殊情况：$m = n - 1$

当操作次数恰好为 $n-1$ 时，最终形成一棵树。通过打表或组合推导可得答案为：

$$\frac{2^{n-2}}{C_{n-1}}$$

其中 $C_{n-1} = \frac{1}{n}\binom{2n-2}{n-1}$ 为第 $(n-1)$ 个卡特兰数。

线性时间复杂度正解

分析操作类型：设进行 $c_0$ 次"孤立点-非孤立点"操作，$c_1$ 次"孤立点-孤立点"操作，则：

$$\begin{cases} c_0 + c_1 = m \\ c_0 + 2c_1 = n \end{cases}$$

解得 $c_1 = n - m$，记 $k = n - m$ 为"孤立点-孤立点"操作的次数。

将问题转化为：在完全图 $K_n$ 的所有 $\binom{n}{2}$ 条边的随机排列中，要求恰好有 $k$ 条边满足"两端点在该边之前均未被连接"（即该边出现时两端仍为孤立点）。

容斥原理的应用：先计算至少 $x$ 条边满足条件的方案数，再容斥得到恰好 $k$ 条的答案。容斥系数为 $(-1)^{x-k}\binom{x}{k}$。

对于固定的 $x$ 条边，将其视为优先约束，构建有向无环图。通过树形结构分析，拓扑序计数结合阶乘与双阶乘化简，最终得到：

$$\sum_{x=k}^{\lfloor n/2 \rfloor} (-1)^{x-k}\binom{x}{k}\binom{n}{2x}(2x-1)!!\frac{(2n-3-2x)!!}{(2n-3)!!}$$

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

const int MAXN = 5000005;
ll fact[MAXN], inv_fact[MAXN];
ll prod[MAXN], inv_prod[MAXN];
ll prefix[MAXN];
int n, k, mod;

ll power(ll a, ll b) {
    ll res = 1;
    while (b) {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}

ll comb(int n, int m) {
    if (n < m || m < 0) return 0;
    return fact[n] * inv_fact[m] % mod * inv_fact[n - m] % mod;
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    cin >> n >> k >> mod;
    k = n - k;
    if (2 * k > n) { cout << 0 << '\n'; return 0; }
    
    int lim = n >> 1;
    for (int i = 1; i <= lim; i++) {
        prefix[i] = (prefix[i - 1] + 2LL * (n - 2 * i) + 1) % mod;
    }
    
    fact[0] = prod[0] = 1;
    for (int i = 1; i <= lim; i++) {
        fact[i] = i * fact[i - 1] % mod;
        prod[i] = prefix[i] * prod[i - 1] % mod;
    }
    
    inv_fact[lim] = power(fact[lim], mod - 2);
    inv_prod[lim] = power(prod[lim], mod - 2);
    for (int i = lim; i >= 1; i--) {
        inv_fact[i - 1] = i * inv_fact[i] % mod;
        inv_prod[i - 1] = prefix[i] * inv_prod[i] % mod;
    }
    
    for (int i = 1; i <= lim; i++) {
        ll pairs = 1LL * (n - 2 * i + 1) * (n - 2 * i + 2) / 2;
        prod[i] = (pairs % mod) * prod[i - 1] % mod;
    }
    
    ll ans = 0;
    for (int i = k; i <= lim; i++) {
        ll coef = ((i - k) & 1) ? (mod - comb(i, k)) : comb(i, k);
        ll term = coef * inv_prod[i] % mod * prod[i] % mod;
        ans = (ans + term) % mod;
    }
    cout << ans % mod;
    return 0;
}