p356.ans1
======================
1/6

=================
p356.ans2
======================
0

=================
p356.cpp
======================
#include <bits/stdc++.h>

using namespace std;

template<typename T1, typename T2>
ostream& operator<<(ostream& out, const pair<T1, T2>& x) {
    return out << x.first << ' ' << x.second;
}

template<typename T1, typename T2>
istream& operator>>(istream& in, pair<T1, T2>& x) {
    return in >> x.first >> x.second;
}

template<typename T>
istream& operator>>(istream& in, vector<T>& a) {
    for(auto& x: a) {
        in >> x;
    }
    return in;
};

template<typename T>
ostream& operator<<(ostream& out, const vector<T>& a) {
    for(auto x: a) {
        out << x << ' ';
    }
    return out;
};

// base and base_digits must be consistent
const int base = 1000000000;
const int base_digits = 9;

struct bigint {
    vector<int> z;
    int sign;

    bigint() : sign(1) {}

    bigint(long long v) { *this = v; }

    bigint(const string& s) { read(s); }

    void operator=(const bigint& v) {
        sign = v.sign;
        z = v.z;
    }

    void operator=(long long v) {
        sign = 1;
        if(v < 0) {
            sign = -1, v = -v;
        }
        z.clear();
        for(; v > 0; v = v / base) {
            z.push_back(v % base);
        }
    }

    bigint operator+(const bigint& v) const {
        if(sign == v.sign) {
            bigint res = v;

            for(int i = 0, carry = 0;
                i < (int)max(z.size(), v.z.size()) || carry; ++i) {
                if(i == (int)res.z.size()) {
                    res.z.push_back(0);
                }
                res.z[i] += carry + (i < (int)z.size() ? z[i] : 0);
                carry = res.z[i] >= base;
                if(carry) {
                    res.z[i] -= base;
                }
            }
            return res;
        }
        return *this - (-v);
    }

    bigint operator-(const bigint& v) const {
        if(sign == v.sign) {
            if(abs() >= v.abs()) {
                bigint res = *this;
                for(int i = 0, carry = 0; i < (int)v.z.size() || carry; ++i) {
                    res.z[i] -= carry + (i < (int)v.z.size() ? v.z[i] : 0);
                    carry = res.z[i] < 0;
                    if(carry) {
                        res.z[i] += base;
                    }
                }
                res.trim();
                return res;
            }
            return -(v - *this);
        }
        return *this + (-v);
    }

    void operator*=(int v) {
        if(v < 0) {
            sign = -sign, v = -v;
        }
        for(int i = 0, carry = 0; i < (int)z.size() || carry; ++i) {
            if(i == (int)z.size()) {
                z.push_back(0);
            }
            long long cur = z[i] * (long long)v + carry;
            carry = (int)(cur / base);
            z[i] = (int)(cur % base);
            // asm("divl %%ecx" : "=a"(carry), "=d"(a[i]) : "A"(cur),
            // "c"(base));
        }
        trim();
    }

    bigint operator*(int v) const {
        bigint res = *this;
        res *= v;
        return res;
    }

    friend pair<bigint, bigint> divmod(const bigint& a1, const bigint& b1) {
        int norm = base / (b1.z.back() + 1);
        bigint a = a1.abs() * norm;
        bigint b = b1.abs() * norm;
        bigint q, r;
        q.z.resize(a.z.size());

        for(int i = a.z.size() - 1; i >= 0; i--) {
            r *= base;
            r += a.z[i];
            int s1 = b.z.size() < r.z.size() ? r.z[b.z.size()] : 0;
            int s2 = b.z.size() - 1 < r.z.size() ? r.z[b.z.size() - 1] : 0;
            int d = ((long long)s1 * base + s2) / b.z.back();
            r -= b * d;
            while(r < 0) {
                r += b, --d;
            }
            q.z[i] = d;
        }

        q.sign = a1.sign * b1.sign;
        r.sign = a1.sign;
        q.trim();
        r.trim();
        return make_pair(q, r / norm);
    }

    friend bigint sqrt(const bigint& a1) {
        bigint a = a1;
        while(a.z.empty() || a.z.size() % 2 == 1) {
            a.z.push_back(0);
        }

        int n = a.z.size();

        int firstDigit = (int)sqrt((double)a.z[n - 1] * base + a.z[n - 2]);
        int norm = base / (firstDigit + 1);
        a *= norm;
        a *= norm;
        while(a.z.empty() || a.z.size() % 2 == 1) {
            a.z.push_back(0);
        }

        bigint r = (long long)a.z[n - 1] * base + a.z[n - 2];
        firstDigit = (int)sqrt((double)a.z[n - 1] * base + a.z[n - 2]);
        int q = firstDigit;
        bigint res;

        for(int j = n / 2 - 1; j >= 0; j--) {
            for(;; --q) {
                bigint r1 =
                    (r - (res * 2 * base + q) * q) * base * base +
                    (j > 0 ? (long long)a.z[2 * j - 1] * base + a.z[2 * j - 2]
                           : 0);
                if(r1 >= 0) {
                    r = r1;
                    break;
                }
            }
            res *= base;
            res += q;

            if(j > 0) {
                int d1 =
                    res.z.size() + 2 < r.z.size() ? r.z[res.z.size() + 2] : 0;
                int d2 =
                    res.z.size() + 1 < r.z.size() ? r.z[res.z.size() + 1] : 0;
                int d3 = res.z.size() < r.z.size() ? r.z[res.z.size()] : 0;
                q = ((long long)d1 * base * base + (long long)d2 * base + d3) /
                    (firstDigit * 2);
            }
        }

        res.trim();
        return res / norm;
    }

    bigint operator/(const bigint& v) const { return divmod(*this, v).first; }

    bigint operator%(const bigint& v) const { return divmod(*this, v).second; }

    void operator/=(int v) {
        if(v < 0) {
            sign = -sign, v = -v;
        }
        for(int i = (int)z.size() - 1, rem = 0; i >= 0; --i) {
            long long cur = z[i] + rem * (long long)base;
            z[i] = (int)(cur / v);
            rem = (int)(cur % v);
        }
        trim();
    }

    bigint operator/(int v) const {
        bigint res = *this;
        res /= v;
        return res;
    }

    int operator%(int v) const {
        if(v < 0) {
            v = -v;
        }
        int m = 0;
        for(int i = z.size() - 1; i >= 0; --i) {
            m = (z[i] + m * (long long)base) % v;
        }
        return m * sign;
    }

    void operator+=(const bigint& v) { *this = *this + v; }
    void operator-=(const bigint& v) { *this = *this - v; }
    void operator*=(const bigint& v) { *this = *this * v; }
    void operator/=(const bigint& v) { *this = *this / v; }

    bool operator<(const bigint& v) const {
        if(sign != v.sign) {
            return sign < v.sign;
        }
        if(z.size() != v.z.size()) {
            return z.size() * sign < v.z.size() * v.sign;
        }
        for(int i = z.size() - 1; i >= 0; i--) {
            if(z[i] != v.z[i]) {
                return z[i] * sign < v.z[i] * sign;
            }
        }
        return false;
    }

    bool operator>(const bigint& v) const { return v < *this; }
    bool operator<=(const bigint& v) const { return !(v < *this); }
    bool operator>=(const bigint& v) const { return !(*this < v); }
    bool operator==(const bigint& v) const {
        return !(*this < v) && !(v < *this);
    }
    bool operator!=(const bigint& v) const { return *this < v || v < *this; }

    void trim() {
        while(!z.empty() && z.back() == 0) {
            z.pop_back();
        }
        if(z.empty()) {
            sign = 1;
        }
    }

    bool isZero() const { return z.empty() || (z.size() == 1 && !z[0]); }

    bigint operator-() const {
        bigint res = *this;
        res.sign = -sign;
        return res;
    }

    bigint abs() const {
        bigint res = *this;
        res.sign *= res.sign;
        return res;
    }

    long long longValue() const {
        long long res = 0;
        for(int i = z.size() - 1; i >= 0; i--) {
            res = res * base + z[i];
        }
        return res * sign;
    }

    friend bigint gcd(const bigint& a, const bigint& b) {
        return b.isZero() ? a : gcd(b, a % b);
    }
    friend bigint lcm(const bigint& a, const bigint& b) {
        return a / gcd(a, b) * b;
    }

    void read(const string& s) {
        sign = 1;
        z.clear();
        int pos = 0;
        while(pos < (int)s.size() && (s[pos] == '-' || s[pos] == '+')) {
            if(s[pos] == '-') {
                sign = -sign;
            }
            ++pos;
        }
        for(int i = s.size() - 1; i >= pos; i -= base_digits) {
            int x = 0;
            for(int j = max(pos, i - base_digits + 1); j <= i; j++) {
                x = x * 10 + s[j] - '0';
            }
            z.push_back(x);
        }
        trim();
    }

    friend istream& operator>>(istream& stream, bigint& v) {
        string s;
        stream >> s;
        v.read(s);
        return stream;
    }

    friend ostream& operator<<(ostream& stream, const bigint& v) {
        if(v.sign == -1) {
            stream << '-';
        }
        stream << (v.z.empty() ? 0 : v.z.back());
        for(int i = (int)v.z.size() - 2; i >= 0; --i) {
            stream << setw(base_digits) << setfill('0') << v.z[i];
        }
        return stream;
    }

    static vector<int> convert_base(
        const vector<int>& a, int old_digits, int new_digits
    ) {
        vector<long long> p(max(old_digits, new_digits) + 1);
        p[0] = 1;
        for(int i = 1; i < (int)p.size(); i++) {
            p[i] = p[i - 1] * 10;
        }
        vector<int> res;
        long long cur = 0;
        int cur_digits = 0;
        for(int i = 0; i < (int)a.size(); i++) {
            cur += a[i] * p[cur_digits];
            cur_digits += old_digits;
            while(cur_digits >= new_digits) {
                res.push_back(int(cur % p[new_digits]));
                cur /= p[new_digits];
                cur_digits -= new_digits;
            }
        }
        res.push_back((int)cur);
        while(!res.empty() && res.back() == 0) {
            res.pop_back();
        }
        return res;
    }

    typedef vector<long long> vll;

    static vll karatsubaMultiply(const vll& a, const vll& b) {
        int n = a.size();
        vll res(n + n);
        if(n <= 32) {
            for(int i = 0; i < n; i++) {
                for(int j = 0; j < n; j++) {
                    res[i + j] += a[i] * b[j];
                }
            }
            return res;
        }

        int k = n >> 1;
        vll a1(a.begin(), a.begin() + k);
        vll a2(a.begin() + k, a.end());
        vll b1(b.begin(), b.begin() + k);
        vll b2(b.begin() + k, b.end());

        vll a1b1 = karatsubaMultiply(a1, b1);
        vll a2b2 = karatsubaMultiply(a2, b2);

        for(int i = 0; i < k; i++) {
            a2[i] += a1[i];
        }
        for(int i = 0; i < k; i++) {
            b2[i] += b1[i];
        }

        vll r = karatsubaMultiply(a2, b2);
        for(int i = 0; i < (int)a1b1.size(); i++) {
            r[i] -= a1b1[i];
        }
        for(int i = 0; i < (int)a2b2.size(); i++) {
            r[i] -= a2b2[i];
        }

        for(int i = 0; i < (int)r.size(); i++) {
            res[i + k] += r[i];
        }
        for(int i = 0; i < (int)a1b1.size(); i++) {
            res[i] += a1b1[i];
        }
        for(int i = 0; i < (int)a2b2.size(); i++) {
            res[i + n] += a2b2[i];
        }
        return res;
    }

    bigint operator*(const bigint& v) const {
        vector<int> a6 = convert_base(this->z, base_digits, 6);
        vector<int> b6 = convert_base(v.z, base_digits, 6);
        vll a(a6.begin(), a6.end());
        vll b(b6.begin(), b6.end());
        while(a.size() < b.size()) {
            a.push_back(0);
        }
        while(b.size() < a.size()) {
            b.push_back(0);
        }
        while(a.size() & (a.size() - 1)) {
            a.push_back(0), b.push_back(0);
        }
        vll c = karatsubaMultiply(a, b);
        bigint res;
        res.sign = sign * v.sign;
        for(int i = 0, carry = 0; i < (int)c.size(); i++) {
            long long cur = c[i] + carry;
            res.z.push_back((int)(cur % 1000000));
            carry = (int)(cur / 1000000);
        }
        res.z = convert_base(res.z, 6, base_digits);
        res.trim();
        return res;
    }
};

int k, n;

void read() { cin >> k >> n; }

bigint derangements(int m) {
    if(m == 0) {
        return bigint(1);
    }
    if(m == 1) {
        return bigint(0);
    }

    bigint d_prev2 = 1, d_prev1 = 0;
    for(int i = 2; i <= m; i++) {
        bigint d_next = (d_prev1 + d_prev2) * (i - 1);
        d_prev2 = d_prev1;
        d_prev1 = d_next;
    }

    return d_prev1;
}

bigint factorial(int m) {
    bigint res = 1;
    for(int i = 2; i <= m; i++) {
        res *= i;
    }
    return res;
}

void solve() {
    // Classify each guessed matching by how it permutes the secret pairing.
    // Fixing exactly k of n pairs means choosing k correct positions and
    // deranging the remaining m = n - k (a permutation with no fixed point).
    // The favourable count is C(n, k) * D_m and the total is n!, so the
    // probability simplifies to D_m / (k! * m!).
    //
    // The derangement number is built with the recurrence
    // D_i = (i - 1) * (D_{i-1} + D_{i-2}), D_0 = 1, D_1 = 0. All quantities are
    // computed with the vendored arbitrary-precision bigint since n can reach
    // 100, then the fraction is reduced by its gcd.

    int m = n - k;
    bigint d = derangements(m);

    if(d.isZero()) {
        cout << 0 << '\n';
        return;
    }

    bigint numerator = d;
    bigint denominator = factorial(k) * factorial(m);

    bigint g = gcd(numerator, denominator);
    numerator /= g;
    denominator /= g;

    cout << numerator << '/' << denominator << '\n';
}

int main() {
    ios_base::sync_with_stdio(false);
    cin.tie(nullptr);

    int T = 1;
    // cin >> T;
    for(int test = 1; test <= T; test++) {
        read();
        // cout << "Case #" << test << ": ";
        solve();
    }

    return 0;
}

=================
p356.in1
======================
2 5

=================
p356.in2
======================
9 10

=================
p356.py
======================
from math import gcd, factorial


def derangements(m):
    """
    Count permutations of m elements with no fixed points.
    """
    if m == 0:
        return 1
    if m == 1:
        return 0
    d_prev2, d_prev1 = 1, 0  # D_0, D_1
    for i in range(2, m + 1):
        d_prev2, d_prev1 = d_prev1, (i - 1) * (d_prev1 + d_prev2)
    return d_prev1


def main():
    # Each matching can be classified by a permutation. Let the secret one
    # be denoted as 1, ..., n. We are interested in the number of permutations
    # p1, ..., pn that match with 1, ..., n in exactly k positions. This is a
    # derangements problem: choose k positions to be fixed points, then the
    # remaining n-k must form a derangement (permutation with no fixed points).
    # The count is C(n, k) * D_{n-k}, where D_m is the m-th derangement number.
    # The probability is C(n, k) * D_{n-k} / n! = D_{n-k} / (k! * (n-k)!).
    #
    # Derangements can be counted via inclusion-exclusion. Let A_i be the set of
    # permutations where i is a fixed point. We want n! - |A_1 U ... U A_n|.
    # By inclusion-exclusion:
    # |A_1 U ... U A_n| = C(n,1) * (n-1)! - C(n,2) * (n-2)! + C(n,3) * (n-3)! - ...
    # So D_n = n! - C(n,1) * (n-1)! + C(n,2) * (n-2)! - ...
    #        = sum_{i=0}^{n} (-1)^i * C(n,i) * (n-i)!
    #        = sum_{i=0}^{n} (-1)^i * n! / i!
    #        = n! * sum_{i=0}^{n} (-1)^i / i!
    #
    # Alternatively, the recurrence D_n = (n-1) * (D_{n-1} + D_{n-2}), with
    # D_0 = 1, D_1 = 0. This follows from considering where element 1 maps: if
    # 1 -> j, then either j -> 1 (giving D_{n-2} ways for the rest) or j maps to
    # something else, equivalent to a derangement of n-1 elements as j and 1 can be
    # treated like one element (D_{n-1} ways).

    k, n = map(int, input().split())
    m = n - k
    d = derangements(m)

    if d == 0:
        print(0)
        return

    numerator = d
    denominator = factorial(k) * factorial(m)

    g = gcd(numerator, denominator)
    numerator //= g
    denominator //= g

    print(f"{numerator}/{denominator}")


if __name__ == "__main__":
    main()

=================
statement.txt
======================
356. Extrasensory Perception
Time limit per test: 0.25 second(s)
Memory limit: 65536 kilobytes
input: standard
output: standard



His Royal Highness the King of Berland Berl XV has just decapitated his magician for being politically incorrect. It is a truth universally acknowledged that every king in the possession of a good fortune is in a want of a decent sorcerer. Heralds all over the country have spread the news about the coming Berland-wide wizard competition that is to end up with the announcement of the winner's name. The winner is to join the royal retinue as His Royal Highness' Counselor.

Magician Resool graduated from Bexford University ten years ago and since that time he has been chased by misfortunes and disasters. His wife, an accomplished, still pretty looking lady, believes her husband can do much better than just grumbling about the unfair life. She suggests that he should try his magic powers at the competition and reminds him of his famous trick. The gist of the trick is to see through people's souls.

After two weeks of training Resool makes himself enter the Great Hall of Berland Palace and demonstrate his exceptional abilities. His Royal Highness the King of Berland Berl XV gives him an assignment: N men and N women are presented, and the fact is known that they form N married couples; Resool is to identify a wife for every husband. After the magician casts a spell, burns some magic herbs, and performs ritual dance he still fails to answer correctly. He only manages to guess K married couples out of N. But The King is full of sympathy for the magician and can't just let the magician go away jobless and sorrowful. So, he hires him as a yard-keeper. Now Resool has a lot of free time on his new job and wants to calculate the chances he had at the competition. He wants to know the probability of guessing correctly exactly K married couples when given N men, N women, and a fact that they form N married couples. This will definitely help him during the competition next time.

Input
Input file contains two integers K and N (1 ≤ N ≤ 100; 0 ≤ K ≤ N) separated by a space.

Output
Print the required probability to the output file. If the answer is zero, simply print "0" (without quotes). Otherwise print the answer in a form of irreducible fraction "A/B" (without quotes), where A and B are positive integers without leading zeroes. See examples below for the format of output.

Example(s)
sample input
sample output
2 5
1/6

sample input
sample output
9 10
0

=================
