“21 天好习惯”第一期－2

分配颜色

Description

好耶！今天不是星期天，但是又是活动课！

今天小A带同学们玩游戏，需要分为两组。每名同学为红队或蓝队。

为了方便分配颜色，小A让同学们站成n排m列。初始时，所有人都为红队成员。

小A可以进行p次操作1和q次操作2，操作解释如下：

• 操作1:把某一行的同学进行取反操作：即红队变为蓝队，蓝队变成红队。
• 操作2:把某一列的同学进行取反操作：即红队变为蓝队，蓝队变成红队。

小A想知道，有多少种方案满足：执行完所有的操作1和操作2之后，蓝队同学恰好有t名。

由于答案可能会很大，所以只需要输出结果对555555555取模的结果即可。

PS：若两个方案中的某一行或列被小A进行操作的次数不同时，则视为两种不同的方案。

Input

输入一行五个整数：n、m、p、q、t。

$1\le n,m,p,q\le 1555$

$0\le t\le 1555^2$

Output

输出一个整数，表示问题的解。

Sample Input 1

2 2 1 1 2

Sample Output 1

4

solution

当选择一次操作1，蓝队的数量将会增加 $cn{t}_b=m-o{p}_2$

当选择一次操作2，蓝队的数量将会增加 $cn{t}_b=n-o{p}_1$

$o{p}_1,o{p}_2$ 分别表示在此次操作前执行了多少次 1 操作或 2 操作。

当有一行或一列被连续操作了两次那么相当于是无效操作，那么最终假设有效操作了 $o{p}_1,o{p}_2$ 蓝队的数量为 $sum=o{p}_1\ast m+o{p}_2\ast n-2\ast o{p}_1\ast o{p}_2$

之后枚举有效操作次数 $(i,j)$ 即可，对于有效的操作次数可以选的方案数为 $tmp=C_n^i\ast C_m^j$。

此时还有剩余的操作次数，通过连续操作两次同一行列用掉剩余的操作次数。

对于剩余的操作次数 $r$ 即为有 $r$ 相同的物体放入 $n$ 个不同的位置的方案数，为 $C_{n+r-1}^{n-1}$

综上可以得到所有可能的方案数为
$$ans=\sum _{(i,j)}^{(p,q)}{C}_n^i\ast C_m^j\ast C_{n+\frac{p-i}{2}-1}^{n-1}\ast C_{m+\frac{q-j}{2}-1}^{m-1}$$
注意组合数为非质数，且考虑到数据范围，可以通过杨辉三角预处理出前5000组 组合数，也可采用 exlucas 求解。

#include <bits/stdc++.h>
using namespace std;
#define endl '\n'
#define int ll
#define rep(i, j, k) for (ll(i) = (j); (i) <= (k); (++i))
#define rrep(i, j, k) for (ll i = (ll)(j); i >= (ll)(k); i--)
typedef long double ld;
typedef long long ll;
const int mod = 555555555;
ll C[5005][5005];
void init() {
    C[1][1] = C[1][0] = 1;
    rep(i, 2, 5000) {
        C[i][0] = C[i][i] = 1;
        rep(j, 1, i - 1) C[i][j] = (C[i - 1][j - 1] + C[i - 1][j]) % mod;
    }
}
int n, m, p, q, t;
int doit(int n, int m) { return C[n + m - 1][m - 1]; }
int get(int i, int j) {
    if (((p - i) & 1) || ((q - j) & 1))
        return 0;
    return C[n][i] * doit((p - i) / 2, n) % mod * C[m][j] % mod *
           doit((q - j) / 2, m) % mod;
}
signed main() {
    init();
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    cin >> n >> m >> p >> q >> t;
    int ans = 0;
    rep(i, 0, min(p, n)) rep(j, 0, min(q, m)) {
        if (i * m + j * n - 2 * i * j != t)
            continue;
        int tmp = get(i, j);
        ans = (ans + tmp + mod) % mod;
    }
    cout << ans << endl;
    return 0;
}

// exlucas 写的，板子有点乱
// 可以预处理模数质因子加速
#include <bits/stdc++.h>
using namespace std;
#define endl '\n'
#define int ll
#define rep(i, j, k) for (ll(i) = (j); (i) <= (k); (++i))
typedef long long ll;
const int mod = 555555555;
// bool prime[N];
// int p[N];
// int cnt;
ll quick_mod(ll a, ll b, ll m) {
    ll ans = 1;
    a %= m;
    while (b) {
        if (b & 1) {
            ans = ans * a % m;
            b--;
        }
        b >>= 1;
        a = a * a % m;
    }
    return ans;
}
// ll Work(ll n, ll p) {
//     ll ans = 0;
//     while (n) {
//         ans += n / p;
//         n /= p;
//     }
//     return ans;
// }
// ll C(ll n, ll m, ll P) {
//     ll ans = 1;
//     for (int i = 0; i < cnt && p[i] <= n; i++) {
//         ll x = Work(n, p[i]);
//         ll y = Work(n - m, p[i]);
//         ll z = Work(m, p[i]);
//         x -= (y + z);
//         ans *= quick_mod(p[i], x, P);
//         ans %= P;
//     }
//     return ans;
// }
ll p = mod, w[105];
const double eps = 1e-8;
ll mul(ll a, ll b, ll t) {
    ll res = a * b - (ll)((long double)a / t * b + eps) * t;
    return (res % t + t) % t;
}
ll exgcd(ll a, ll b, ll &x, ll &y) {
    if (!b) {
        x = 1, y = 0;
        return a;
    }
    ll d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}
ll inv(ll a, ll b) {
    ll x, y;
    exgcd(a, b, x, y);
    return (x % b + b) % b;
}
ll CRT(ll b, ll t) { return mul(mul(b, inv(p / t, t), p), p / t, p); }
ll qpow(ll a, ll b, ll t) {
    ll p = 1;
    while (b) {
        if (b & 1)
            p = mul(p, a, t);
        a = mul(a, a, t), b >>= 1;
    }
    return p;
}

ll fac(ll n, ll pi, ll pk) {
    if (!n)
        return 1;
    ll res = 1;
    rep(i, 2, pk) if (i % pi) res *= i, res %= pk;
    res = qpow(res, n / pk, pk);
    rep(i, 2, n % pk) if (i % pi) res *= i, res %= pk;
    return res * fac(n / pi, pi, pk) % pk;
}

ll C(ll n, ll m, ll pi, ll pk) {
    ll d = fac(n, pi, pk), d1 = fac(m, pi, pk), d2 = fac(n - m, pi, pk);
    ll k = 0;
    for (ll i = n; i; i /= pi)
        k += i / pi;
    for (ll i = m; i; i /= pi)
        k -= i / pi;
    for (ll i = n - m; i; i /= pi)
        k -= i / pi;
    return mul(mul(d, inv(d1, pk), pk), mul(qpow(pi, k, pk), inv(d2, pk), pk),
               pk);
}

ll exlucas(ll n, ll m) {
    ll res = 0, tmp = p, pk;
    ll lim = sqrt(p);
    rep(i, 2, lim) if (!(tmp % i)) {
        pk = 1;
        while (!(tmp % i))
            pk *= i, tmp /= i;
        res += CRT(C(n, m, i, pk), pk), res %= p;
    }
    if (tmp > 1)
        res += CRT(C(n, m, tmp, tmp), tmp), res %= p;
    return res;
}

signed main() {
    int n, m, p, q, t;
    cin >> n >> m >> p >> q >> t;
    int ans = 0;
    rep(i, 0, min(p, n)) {
        rep(j, 0, min(q, m)) {
            if (i * m + j * n - 2 * i * j != t)
                continue;
            int r = p - i, c = q - j;
            if (r % 2 == 0 and c % 2 == 0)
                r /= 2, c /= 2;
            else
                continue;
            // cout << ans << endl;
            // int tmp = quick_mod(n, r, mod) * quick_mod(m, c, mod) % mod *
            //   exlucas(n, i) % mod * exlucas(m, j) % mod;
            int tmp = exlucas(n + r - 1, n - 1) * exlucas(m + c - 1, m - 1) %
                      mod * exlucas(n, i) % mod * exlucas(m, j) % mod;
            ans = (ans + tmp) % mod;
        }
    }
    cout << ans << endl;
    return 0;
}