p441.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;
};

const int mod = 2007;

int64_t n;
int k;

using matrix = vector<vector<int64_t>>;

matrix multiply(const matrix& a, const matrix& b) {
    int sz = a.size();
    matrix c(sz, vector<int64_t>(sz, 0));
    for(int i = 0; i < sz; i++) {
        for(int j = 0; j < sz; j++) {
            for(int l = 0; l < sz; l++) {
                c[i][j] = (c[i][j] + a[i][l] * b[l][j]) % mod;
            }
        }
    }
    return c;
}

matrix mat_pow(matrix m, int64_t p) {
    int sz = m.size();
    matrix result(sz, vector<int64_t>(sz, 0));
    for(int i = 0; i < sz; i++) {
        result[i][i] = 1;
    }
    while(p > 0) {
        if(p & 1) {
            result = multiply(result, m);
        }
        m = multiply(m, m);
        p >>= 1;
    }
    return result;
}

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

void solve() {
    // We can see that N is very large while K is small, so we should think of
    // either a closed form. Let's start with a simplification by making the
    // friends indistinguishable - we can think of this as sorting them by their
    // lowest album ID. After making this reduction, we can notice that each
    // sorted configuration is mapped to exactly K! of the original. The main
    // reason to make this change is that it let's us forget the ordering of the
    // friends and makes the transitions easier.
    //
    // We will start giving albums from 1 until N. Let's say we have P people
    // that have already received an album. We now have only two options -
    // give the current album to any of the people that already have albums in P
    // ways, or give it to someone that still hasn't received an album (bool)(K
    // - P) ways. This means we have a direct dp[N][P] solution and furthermore
    // each state depends only on the previous one. This is already enough to
    // create a O(K^3 log N) solution, because we can directly express dp[N][.]
    // as a linear combination of dp[N - 1][.] as the coefficients don't depend
    // on N. Then we can do binary exponentiation.
    //
    // There is also an alternative approach based on inclusion exclusion. Let's
    // drop the symmetry division by K!, and instead solve the original problem.
    // If we didn't have the "at least 1 album" constraint, the solution would
    // have been K^n. However, we can do inclusion-exclusion here. Particularly,
    // there are C(K, j) ways to have at least j people receive 0 albums, and
    // then (K - j)^n to distribute the rest. This leads us to the closed form
    // of sum_{j=0}^{K} C(K, j) * (-1)^j * (K - j)^n. With a precompute of the
    // factorials up to K, we can solve this in just O(K log N).
    //
    // People familiar with number theory might notice the problem here asks for
    // k! * S(n, k), or the Stirling numbers of the second kind. The two
    // approaches above are precisely the ways to compute it: due to a recursive
    // relation, and with inclusion-exclusion. Here we implement the first
    // solution.

    if(n < k) {
        cout << 0 << "\n";
        return;
    }

    matrix trans(k + 1, vector<int64_t>(k + 1, 0));
    for(int p = 0; p <= k; p++) {
        if(p > 0) {
            trans[p][p] = p % mod;
        }
        if(p < k) {
            trans[p + 1][p] = 1;
        }
    }

    matrix result = mat_pow(trans, n);
    cout << result[k][0] << "\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;
}

=================
p441.in1
======================
3 2

=================
statement.txt
======================
441. Set Division
Time limit per test: 0.25 second(s)
Memory limit: 262144 kilobytes
input: standard
output: standard



Ruslan has K friends. And all of them have birthday tomorrow. He has bought N different photo albums yesterday, and wants to present these photo albums to his friends. Of course, he can't give less than 1 photo album to anybody. Your task is to calculate how many possible ways is there to do it.

All photo albums are different. Two distributions are considered the same if they differ only by order of albums in the gifts or by the persons receiving gifts.

Input
In the only line of the input there are two numbers separated by a space — , .

Output
Output should contain one number — the number of possible ways to distribute the albums modulo 2007.

Example(s)
sample input
sample output
3 2 
3

=================
p441.ans1
======================
3

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