Inclusion-exclusion in code

Why this matters when you build with AI

Many counting questions ask for the number of arrangements that have none of several bad features. A model that has seen the textbook formula for two bad features often tries to extend it to more by guessing. The code compiles, the sample passes, and the answer is wrong when two bad features overlap or the list grows.

Our runs show the kind of failure to expect. A version of the subset code that treated a repeated forbidden pair as two separate bad features returned 0 for n = 2 with the pair (0, 0) listed twice. The correct answer is 1, because only the arrangement [1, 0] avoids the pair. Nothing reported an error. A second version used the C++ remainder on a negative number and printed -7 for -7 % 1000000007. Both bugs pass a sample that has no repeats and no negative intermediate values.

This lesson covers the principle, the bitmask enumeration that implements it, the sign rules when the answer is taken modulo a prime, and three classic shapes: derangements, surjections and permutations with forbidden pairs. The code is C++17, compiled with g++ -std=c++17 -O2 -Wall -Wextra using GCC 15.2 in our scratch folder. The outputs below are from those runs. Timings depend on the machine, and the counts do not.

The principle

Inclusion-exclusion counts the objects that have none of m bad properties. Call the set of objects that have bad property i the event A_i. We want the number of objects that lie in no event.

Start with two events. Take 10 students who sit in 10 seats. Event A_1 is "student 1 sits in seat 1", and event A_2 is "student 2 sits in seat 2". Without restrictions there are 10! = 3628800 seatings. Each single event fixes one student, so 9! = 362880 seatings are in A_1, and the same number are in A_2. Both events together fix two students, so 8! = 40320 seatings are in both.

If we subtract the two single events from the total, the seatings in both events were removed twice. We add them back once:

total - |A_1| - |A_2| + |A_1 and A_2| = 3628800 - 362880 - 362880 + 40320 = 2943360

That is the whole idea. For m events, the rule is a sum over every subset S of the events. The subset S with no events is the total. Each subset adds or subtracts the number of objects that lie in all events of S. The sign is minus when S has an odd number of events and plus when it has an even number:

answer = sum over all subsets S of (-1)^|S| * (number of objects in every event of S)

You can check the rule on one object that lies in exactly two events, A_1 and A_2, and no others. The subsets of that contain only events the object lies in are: the empty set (plus 1), (minus 1), (minus 1), and (plus 1). The total is 1 - 1 - 1 + 1 = 0, so the object is not counted, which is what we want. An object in no event is counted only by the empty set, so it is counted once.

Subsets as bitmasks

With m events there are 2^m subsets. We store each subset as an integer. Bit i of the integer is 1 when event i is in the subset. The loop is then a plain loop over the integers from 0 to 2^m - 1:

for (unsigned long long mask = 0; mask < (1ULL << m); mask++) {
    // bit i of mask is 1 when event i is in this subset
    bool inSubset = (mask >> i) & 1;
    int size = __builtin_popcountll(mask);  // number of events in the subset
    int sign = (size % 2 == 0) ? +1 : -1;   // (-1)^|S|
}

The cost is 2^m subsets times the work for each subset. For m = 20 that is about one million subsets, which is fast. For m = 30 it is about a billion, which is too slow in most contests. Mask enumeration is for a few properties, up to about 20 or 25. When m is larger, the event structure has to be used, for example with a DP over the events that stays small.

Use 1ULL << i and not 1 << i once i can reach 31. In a signed 32-bit int, 1 << 31 is undefined and 1 << 40 is worse.

Sign handling modulo a prime

Counting problems often ask for the answer modulo p = 1000000007. Here the subtraction is the risky part. Three rules keep it correct.

Keep the running answer in the range from 0 to p - 1. When a term is subtracted, write ans = (ans - term + MOD) % MOD. The + MOD keeps the value nonnegative before the %.
In C++, the % operator keeps the sign of the left operand. So -7 % 1000000007 is -7, not 1000000000. Our run printed -7. The safe form is ((x % MOD) + MOD) % MOD, which printed 1000000000 for that value.
A product of two values below p is below about 10^18, which fits in a signed 64-bit integer. Take % after each multiplication, not at the end of a long chain.

Inclusion-exclusion needs no division, so no modular inverse is needed for the sign. Closed forms of some classic answers do contain divisions, such as n! / j!. The recurrences in the next section avoid those divisions.

Three classic shapes

Derangements

A derangement of n items is a permutation with no fixed point. That means p[i] differs from i for every i. Let the event A_i be "p[i] = i". A set S of j events fixes j items, so the intersection has (n - j)! permutations. There are C(n, j) sets of size j. So

D(n) = sum over j from 0 to n of (-1)^j * C(n, j) * (n - j)!

The same sum is sum over j of (-1)^j * n! / j!. A simpler way to compute it is the recurrence D(n) = (n - 1) * (D(n - 1) + D(n - 2)) with D(0) = 1 and D(1) = 0. Our demo printed:

derangements D(0..10): 1 0 1 2 9 44 265 1854 14833 133496 1334961

Surjections

A surjection from an n-element set onto a k-element set is a function where every target value is used at least once. Let the event B_t be "target t is not used". A set of j missed targets leaves k - j available values, so the functions that miss all of them number (k - j)^n. Hence

S(n, k) = sum over j from 0 to k of (-1)^j * C(k, j) * (k - j)^n

The formula needs 0^0 = 1 for n = 0. Here is the demo table for rows n = 0 to 5 and columns k = 0 to 4. Check one entry by hand: S(4, 2) should be 2^4 - 2 = 14, and the table has 14.

surjections S(n,k), rows n=0..5, columns k=0..4:
   1   0   0   0   0
   0   1   0   0   0
   0   1   2   0   0
   0   1   6   6   0
   0   1  14  36  24
   0   1  30 150 240

Permutations with forbidden pairs

This is the shape of the worked example below. We are given a list of forbidden pairs (pos, val). A permutation p of n items is bad when p[pos] = val for some pair. Let the event A_i be "p[pos_i] = val_i".

A set S of pairs can all happen together only when no two of its pairs share a position and no two share a value. A position holds one value, and a value sits in one position. If the pairs of S conflict, the intersection is empty and the term is 0. If S has c pairs with no conflicts, the c positions and c values are fixed, and the other n - c positions take the other n - c values in (n - c)! ways. So

answer = sum over conflict-free subsets S of (-1)^|S| * (n - |S|)!

In code, two bitmasks record which positions and values are already used by the subset. A conflict is a bit that is already set.

Worked example: the workshop assignment

Seven jobs must go to seven machines, numbered 0 to 6. Each job goes to exactly one machine, and each machine takes exactly one job. A list of m forbidden job-machine pairs says which pairs cannot be used. Count the assignments that use no forbidden pair, modulo 1000000007.

Input: a first line with n m (n at most 60, m at most 20). Then m lines, each with a job and a machine. Output: one number.

The full solution is below. Save it as solution.cpp and compile it with g++ -std=c++17 -O2 -Wall -Wextra solution.cpp -o solution.

#include <bits/stdc++.h>
using namespace std;
using ll = long long;
using ull = unsigned long long;
 
const ll MOD = 1000000007LL;
 
ll power(ll b, ll e) {
    b %= MOD;
    if (b < 0) b += MOD;
    ll r = 1;
    while (e > 0) {
        if (e & 1) r = r * b % MOD;
        b = b * b % MOD;
        e >>= 1;
    }
    return r;
}
 
// Factorials and inverse factorials up to N. MOD is prime, so inverses exist.
struct Comb {
    vector<ll> f, inv;
    explicit Comb(int N) : f(N + 1), inv(N + 1) {
        f[0] = 1;
        for (int i = 1; i <= N; i++) f[i] = f[i - 1] * i % MOD;
        inv[N] = power(f[N], MOD - 2);
        for (int i = N; i > 0; i--) inv[i - 1] = inv[i] * i % MOD;
    }
    ll C(int n, int k) const {
        if (n < 0 || k < 0 || k > n) return 0;
        return f[n] * inv[k] % MOD * inv[n - k] % MOD;
    }
};
 
// D(n): permutations of 0..n-1 with no fixed point.
// Recurrence: D(n) = (n - 1) * (D(n - 1) + D(n - 2)), with D(0) = 1, D(1) = 0.
ll derangements(int n) {
    ll a = 1, b = 0;  // D(0), D(1)
    if (n == 0) return a;
    for (int i = 2; i <= n; i++) {
        ll c = (ll)(i - 1) * ((a + b) % MOD) % MOD;
        a = b;
        b = c;
    }
    return b;
}
 
// Surjections from an n-set onto a k-set: sum_j (-1)^j * C(k, j) * (k - j)^n.
ll surjections(int n, int k, const Comb& cb) {
    ll ans = 0;
    for (int j = 0; j <= k; j++) {
        ll term = cb.C(k, j) * power(k - j, n) % MOD;
        ans = (j % 2 == 0) ? (ans + term) % MOD : (ans - term + MOD) % MOD;
    }
    return ans;
}
 
// Permutations p of 0..n-1 with p[pos] != val for every forbidden pair (pos, val).
// Inclusion-exclusion over subsets of the forbidden pairs. Requires n <= 63.
ll count_avoiding(int n, vector<pair<int, int>> bad, const Comb& cb) {
    // Keep only pairs that a permutation can use, and drop repeated pairs.
    // A repeated pair is one event, not two events.
    vector<pair<int, int>> ev;
    for (auto [p, v] : bad)
        if (0 <= p && p < n && 0 <= v && v < n) ev.push_back({p, v});
    sort(ev.begin(), ev.end());
    ev.erase(unique(ev.begin(), ev.end()), ev.end());
 
    int m = ev.size();
    ll ans = 0;
    for (ull mask = 0; mask < (1ULL << m); mask++) {
        // Pairs in this subset must not share a position or a value.
        // Otherwise no permutation makes all of them true, and the term is 0.
        ull usedPos = 0, usedVal = 0;
        int cnt = 0;
        bool ok = true;
        for (int i = 0; i < m && ok; i++) {
            if (!((mask >> i) & 1)) continue;
            auto [p, v] = ev[i];
            if (((usedPos >> p) & 1) || ((usedVal >> v) & 1)) {
                ok = false;
            } else {
                usedPos |= 1ULL << p;
                usedVal |= 1ULL << v;
                cnt++;
            }
        }
        if (!ok) continue;
        // Fixing cnt pairs leaves n - cnt free positions: (n - cnt)! permutations.
        ll term = cb.f[n - cnt];
        if (__builtin_popcountll(mask) % 2 == 1)
            ans = (ans - term + MOD) % MOD;
        else
            ans = (ans + term) % MOD;
    }
    return ans;
}
 
// Input: n m, then m lines "pos val". Output: the count modulo 1000000007.
int main() {
    int n, m;
    cin >> n >> m;
    vector<pair<int, int>> bad(m);
    for (auto& [p, v] : bad) cin >> p >> v;
    Comb cb(max(n, 1));
    cout << count_avoiding(n, bad, cb) << "\n";
    return 0;
}

On the sample input below, it printed:

A brute force that tries all 5040 permutations of seven items and keeps the ones with no forbidden pair also gives 2119. The brute force is in Python, and it is a check, not part of the answer.

Testing against brute force

A sample is one input. We also ran a stress test, which generates random inputs, runs the solution and a brute force, and compares them. The brute force tries every permutation for n up to 8, so it is only for small inputs. The stress harness includes solution.cpp with a macro that renames its main, so the same functions are tested:

// Brute force: try every permutation and count the ones with no forbidden pair.
ll brute_avoiding(int n, const vector<pair<int, int>>& bad) {
    vector<int> p(n);
    iota(p.begin(), p.end(), 0);
    ll cnt = 0;
    do {
        bool ok = true;
        for (auto [pos, v] : bad)
            if (0 <= pos && pos < n && 0 <= v && v < n && p[pos] == v) { ok = false; break; }
        if (ok) cnt++;
    } while (next_permutation(p.begin(), p.end()));
    return cnt % MOD;
}

The harness checks four things. The derangement recurrence, the mask count on the diagonal pairs (i, i), and the brute force agree for n from 0 to 9. The surjection formula agrees with a brute force over all functions for 35 pairs of n and k. The mask count agrees with the brute force on 3000 random inputs with n up to 8, up to 12 pairs, out-of-range pairs and repeated pairs. The output was:

derangements n=0..9: recurrence, mask count and brute force agree
surjections n=0..6, k=0..4: 35 cases agree
random forbidden-pair sets: 3000 cases agree with brute force

The largest allowed case, with m = 20 and n = 60, ran in 33 ms on our machine. The answer there was 60077476, but it has no brute force check, because 60! permutations cannot be enumerated.

The pitfalls, in one demo

The demo program prints the three tables above and two pitfalls. It uses the same functions. Its output was:

worked example (n=7, 6 pairs): 2119
pair (0,0) listed twice, n=2: correct 1, without dedupe 0
-7 % 1000000007 in C++ = -7, fixed = 1000000000

The second line shows what happens when repeated pairs are counted as two events. The two identical events conflict with themselves in the subset check, so the subsets with both copies give 0, and the sum comes out to 0 instead of 1.

The Python snippet below checks the derangement formula against brute force. It uses exact integers, so there is no modulus. It runs without input.

Python · editable

Its output was:

0 1 1
1 0 0
2 1 1
3 2 2
4 9 9
5 44 44
6 265 265
7 1854 1854
8 14833 14833

Python is too slow for the mask version at m = 20, so the subset enumeration is written in C++ only. The Python snippet above is fast enough because it stops at n = 8.

Complexity

The mask version takes O(2^m * m) time for the subsets plus O(n) for the factorial table. Memory is O(n + m). With m = 20 and n = 60 it ran in 33 ms on our machine. The time will be different on yours.

Common mistakes

Repeated events counted twice. The symptom is a result that is too small when a pair is listed twice. In our run, n = 2 with the pair (0, 0) listed twice gave 0 instead of 1. The fix is to sort the pairs and remove duplicates before the mask loop.
A wrong sign. The empty subset is the total, and it is added. Single events are subtracted. Use popcount % 2 to choose the sign, and check the result on a known value. For the diagonal pairs with n = 5, the answer must be D(5) = 44.
A negative remainder. The symptom is a negative number on the screen where a count should be. -7 % MOD is -7 in C++. Add MOD before the %, as in the sign rules above.
Shifts on a 32-bit int. 1 << i with i at 31 or more is undefined. Use 1ULL << i. Positions and values must also be at most 63 for a 64-bit mask, which is why the solution asks for n at most 63.
Out-of-range pairs. A pair with a position of -1 or a value of n indexes past an array. In our first stress harness, the brute force read p[-1] and stopped with the assertion __n < this->size() from the vector's indexing check. That bug was in the harness, not in the solution, which drops out-of-range pairs at the start. Test the edges before you trust the code.
Too many events. A list of 40 forbidden pairs gives 2^40 subsets, which will not finish. The fix is a different method, not a faster loop. Check m against the limit in the statement before you pick the method.

Checklist

Name the bad events in words first. Then decide which sets of events can happen together.
Remove repeated events before enumerating.
Give each subset the sign (-1) to the number of events in it.
Keep every running value in the range from 0 to MOD - 1, and add MOD before you subtract.
Use 1ULL for masks, and check that n and m are within the limits of the statement.
Check a known small value, such as D(5) = 44, before you run anything bigger.
Compare with a brute force on random small inputs, including repeats and out-of-range values.

Exercise

Write a program that counts permutations of n items avoiding a list of forbidden (position, value) pairs. Use mask enumeration with n at most 9 and m at most 12. Write a brute force that tries all permutations. Run a stress loop of 500 random inputs and confirm that the two agree on every case.

Your deliverables are the solution, the brute force, and the loop. Use a fixed seed so a failing input can be reproduced.

Here is a prompt for an AI tool that can help with the harness. Read its output before you run it, because the tool may generate inputs outside the limits.

Write a C++17 program that counts permutations of 0..n-1 that avoid a list
of forbidden pairs (pos, val), using inclusion-exclusion over subsets of the
pairs. Use bitmasks for the used positions and values. Remove repeated pairs
and ignore pairs outside the range. Keep all values modulo 1000000007 and
add MOD before subtracting. Then write a separate brute force that tries
every permutation for n up to 8. Write a random generator with a fixed seed
that makes n from 0 to 8, up to 12 pairs, some out-of-range values and some
repeated pairs. The generator must print the input that fails. Do not guess
the outputs. Show the commands and the real output.

A good answer has the repeat removal before the mask loop, a sign taken from the popcount, and the generator printing the failing input. Check that the brute force really tries every permutation and that the generator uses the same limits as the statement.

Show the answer

The answer is the program described in the worked example, plus the brute force and the random loop. Run the loop with for seed in $(seq 1 500), and pass the seed to the generator. A correct run prints no mismatch lines. If any line says mismatch, print the input and recheck the repeat removal, the sign, and the range checks.

Practice problems

Compute D(6) by the recurrence, without a computer.

Show the answer

Use D(0) = 1, D(1) = 0. D(2) = 1, D(3) = 2, D(4) = 9, D(5) = 44, and D(6) = 5 * (44 + 9) = 265. The demo printed 265 for D(6).
Compute S(5, 3), the number of surjections from a 5-element set onto a 3-element set, with the formula.

Show the answer

The formula gives 3^5 - C(3, 1) * 2^5 + C(3, 2) * 1^5 = 243 - 96 + 3 = 150. The demo table row for n = 5 has 150 in column k = 3.
Count permutations of 0 to 3 that avoid the forbidden pairs (0, 0) and (1, 1).

Show the answer

Hint: each single event has 3! = 6 permutations, and both together have 2! = 2. The answer is 4! - 3! - 3! + 2! = 24 - 6 - 6 + 2 = 14.
Count permutations of 0 to 5 that avoid the forbidden pairs (0, 0) and (0, 1). Both pairs use position 0.

Show the answer

Hint: the two events cannot happen together, because position 0 holds one value. So the intersection of the two events is empty. The answer is 6! - 5! - 5! + 0 = 720 - 120 - 120 = 480.
In C++, a program keeps its answer in unsigned long long and writes ans = ans - term;. Say what happens when term is larger than ans, and how to fix it.

Show the answer

Hint: unsigned arithmetic wraps around modulo 2^64, not modulo p. The result is a huge number, not a negative one. Keep the answer in a signed type and reduce with (ans - term + MOD) % MOD, or use the same fix in an unsigned type by adding MOD before the subtraction.

Summary

Inclusion-exclusion turns "none of these bad properties" into a signed sum over subsets of properties. Bitmasks enumerate the subsets, and the sign is the parity of the subset size. Modulo a prime, the sign needs a nonnegative remainder, and repeated events must be removed first.

Reference for the general method: Richard P. Stanley, Enumerative Combinatorics, Volume 1, second edition, Cambridge University Press, 2011. The edition details come from the publisher listing. The derangement and surjection formulas above are checked against brute force in our runs, and the reference is for the general theory.