Inclusion-Exclusion

Why Inclusion-Exclusion Matters

Suppose you need to count the integers from 1 through 30 that are divisible by 2, 3, or 5. There are 15 multiples of 2, 10 multiples of 3, and 6 multiples of 5. Adding gives 31, which already exceeds the size of the entire range. The problem is overlap: 6 belongs to two groups, and 30 belongs to all three.

Inclusion–Exclusion counts the union of overlapping sets. It adds individual set sizes, subtracts pairwise intersections, adds triple intersections, and continues with alternating signs. The result counts every qualifying object exactly once, even when it satisfies several conditions.

This lesson explains the counting principle first, then implements one useful specialization: counting integers in 1…N divisible by at least one number in a short list. You will see every subset of the three-rule example, understand why intersections use the least common multiple, and implement the same algorithm in Python, C++, and Java. The principle also supports probability calculations, restricted arrangements, and counting assignments that avoid forbidden conditions.

Core Concepts and Terminology

Let U be a finite universe and let A₁, …, Aₖ be subsets of U. The union contains objects belonging to at least one set. An intersection contains objects belonging to every selected set. For two sets, the rule is size(A ∪ B) = size(A) + size(B) − size(A ∩ B).

For three sets, add the sizes of A, B, and C; subtract the sizes of A ∩ B, A ∩ C, and B ∩ C; then add size(A ∩ B ∩ C). An object in all three sets is initially counted three times, removed three times, and restored once: 3 − 3 + 1 = 1.

In general, for every nonempty subset S of the set indices, add the size of the corresponding intersection if S has odd size, and subtract it if S has even size. This is a sum over subsets of conditions, not a sum over individual objects.

For our specialization, U = {1, 2, …, N}, with N a nonnegative integer. Each positive divisor dᵢ defines Aᵢ = {x in U: dᵢ divides x}. An intersection for a selected subset S consists of multiples of L = lcm(dᵢ for i in S), so its size is floor(N/L). An empty divisor list describes no qualifying sets and gives a union count of zero.

How Inclusion-Exclusion Works, Step by Step

We enumerate every nonempty subset of the k divisor positions using masks from 1 through 2ᵏ − 1. Bit j selects divisor dⱼ, with bit 0 written at the right of the binary representation. For [2, 3, 5], mask 101 selects 2 and 5. A list may contain equal or redundant divisors; its positions still define separate conditions, and the formula remains valid.

For each mask, start with L = 1 and selected_count = 0. Scan the divisor positions. Whenever a bit is set, increase selected_count and replace L with lcm(L, d). Compute this safely using g = gcd(L, d) and q = L/g. The proposed new LCM is q × d.

Before multiplying, test whether q > floor(N/d). If so, the proposed LCM exceeds N, so the intersection has no elements in 1…N. Skip the entire mask immediately: adding more divisors cannot make an LCM smaller. This comparison also prevents fixed-width multiplication overflow. For N = 0, return zero after validating the inputs.

Otherwise, after scanning the selected positions, calculate term = floor(N/L). Add term when selected_count is odd and subtract it when selected_count is even. The implementations reduce the running total modulo M = 1,000,000,007 after each term. The diagrams and trace show ordinary integer totals, which are all smaller than M in these examples.

flowchart TD
    A["Validate N, positive divisors, and k at most 20"] --> B{"N = 0 or k = 0?"}
    B -->|Yes| R["Return 0"]
    B -->|No| C["total = 0; mask = 1"]
    C --> D{"mask less than 2^k?"}
    D -->|No| Z["Return total modulo M"]
    D -->|Yes| E["L = 1; selected = 0; j = 0"]
    E --> F{"j less than k?"}
    F -->|No| G["term = floor(N / L)"]
    F -->|Yes| H{"Bit j selected?"}
    H -->|No| I["j = j + 1"]
    I --> F
    H -->|Yes| J["selected += 1; q = L / gcd(L, d[j])"]
    J --> K{"q greater than floor(N / d[j])?"}
    K -->|Yes| N["mask = mask + 1"]
    K -->|No| L["L = q times d[j]"]
    L --> I
    G --> O{"selected is odd?"}
    O -->|Yes| P["Add term; normalize modulo M"]
    O -->|No| Q["Subtract term; normalize modulo M"]
    P --> N
    Q --> N
    N --> D

Worked Example with a Full Trace

Take N = 30 and divisors [2, 3, 5]. The seven nonempty masks are listed in numeric order, not grouped by subset size. Therefore the signs appear as +, +, −, +, −, −, +.

1. 001 selects {2}: L changes 1 → 2; term = floor(30/2) = 15. One set means add: total 0 → 15.
2. 010 selects {3}: L changes 1 → 3; term = floor(30/3) = 10. One set means add: total 15 → 25.
3. 011 selects {2, 3}: L changes 1 → 2 → 6; term = floor(30/6) = 5. Two sets mean subtract: total 25 → 20.
4. 100 selects {5}: L changes 1 → 5; term = floor(30/5) = 6. One set means add: total 20 → 26.
5. 101 selects {2, 5}: L changes 1 → 2 → 10; term = floor(30/10) = 3. Two sets mean subtract: total 26 → 23.
6. 110 selects {3, 5}: L changes 1 → 3 → 15; term = floor(30/15) = 2. Two sets mean subtract: total 23 → 21.
7. 111 selects {2, 3, 5}: L changes 1 → 2 → 6 → 30; term = floor(30/30) = 1. Three sets mean add: total 21 → 22.

Grouping the same terms by intersection size gives 15 + 10 + 6 − 5 − 3 − 2 + 1 = 22. The eight nonqualifying numbers are 1, 7, 11, 13, 17, 19, 23, and 29, providing the independent check 30 − 8 = 22.

Implementation in Python, C++, and Java

Each panel is a standalone program. All three enumerate masks in the same order, calculate intersections using a bounded LCM, and return the union count modulo M = 1,000,000,007. The small sample answers below are also the full counts because they are smaller than M. For larger answers, a residue must not be mistaken for the full count.

The input contract is N ≥ 0, strictly positive integer divisors, and at most 20 divisor positions. The limit keeps the mask representation safe and makes the exponential cost explicit; it is a practical guard, not a mathematical limit of Inclusion–Exclusion. C++ and Java accept signed 64-bit values for N and the divisors. Python integers can grow further, with correspondingly larger arithmetic costs.

Python uses arbitrary-precision integers. C++ uses long long and tests the LCM bound before multiplying; its modular accumulator stays within a few multiples of M. Java uses the same safe long LCM calculation and a BigInteger accumulator. All subtraction is normalized into the range 0…M − 1. Divisors equal to 1, duplicates, empty lists, and intersections larger than N are handled directly.

from math import gcd

MOD = 1_000_000_007


def count_union_mod(n: int, divisors: list[int]) -> int:
    if n < 0:
        raise ValueError("N must be nonnegative")
    if len(divisors) > 20:
        raise ValueError("At most 20 divisors are supported")
    if any(d <= 0 for d in divisors):
        raise ValueError("Divisors must be positive")
    k = len(divisors)
    if n == 0 or k == 0:
        return 0

    total = 0
    for mask in range(1, 1 << k):
        lcm_value = 1
        selected = 0
        too_large = False
        for j in range(k):
            if mask & (1 << j):
                selected += 1
                d = divisors[j]
                q = lcm_value // gcd(lcm_value, d)
                if q > n // d:
                    too_large = True
                    break
                lcm_value = q * d
        if too_large:
            continue
        term = (n // lcm_value) % MOD
        if selected % 2 == 1:
            total = (total + term) % MOD
        else:
            total = (total - term) % MOD
    return total


def main() -> None:
    samples = [
        (30, [2, 3, 5]),
        (12, [2, 4]),
        (10, [2, 2]),
        (5, [6, 10]),
        (0, [2, 3]),
        (10, []),
        (100, [1, 7]),
    ]
    print(f"modulus = {MOD}")
    for n, divisors in samples:
        print(f"N={n}, divisors={divisors} -> {count_union_mod(n, divisors)}")


if __name__ == "__main__":
    main()

#include <iostream>
#include <stdexcept>
#include <utility>
#include <vector>

constexpr long long MOD = 1'000'000'007LL;

long long gcd_positive(long long a, long long b) {
    while (b != 0) {
        const long long next = a % b;
        a = b;
        b = next;
    }
    return a;
}

long long count_union_mod(long long n, const std::vector<long long>& divisors) {
    if (n < 0) {
        throw std::invalid_argument("N must be nonnegative");
    }
    if (divisors.size() > 20) {
        throw std::invalid_argument("At most 20 divisors are supported");
    }
    for (long long d : divisors) {
        if (d <= 0) {
            throw std::invalid_argument("Divisors must be positive");
        }
    }
    const int k = static_cast<int>(divisors.size());
    if (n == 0 || k == 0) {
        return 0;
    }

    long long total = 0;
    const int limit = 1 << k;
    for (int mask = 1; mask < limit; ++mask) {
        long long lcm_value = 1;
        int selected = 0;
        bool too_large = false;
        for (int j = 0; j < k; ++j) {
            if ((mask & (1 << j)) != 0) {
                ++selected;
                const long long d = divisors[j];
                const long long q = lcm_value / gcd_positive(lcm_value, d);
                if (q > n / d) {
                    too_large = true;
                    break;
                }
                lcm_value = q * d;
            }
        }
        if (too_large) {
            continue;
        }
        const long long term = (n / lcm_value) % MOD;
        if (selected % 2 == 1) {
            total = (total + term) % MOD;
        } else {
            total = (total - term + MOD) % MOD;
        }
    }
    return total;
}

int main() {
    const std::vector<std::pair<long long, std::vector<long long>>> samples = {
        {30, {2, 3, 5}}, {12, {2, 4}}, {10, {2, 2}},
        {5, {6, 10}}, {0, {2, 3}}, {10, {}}, {100, {1, 7}}
    };
    std::cout << "modulus = " << MOD << '\n';
    for (const auto& sample : samples) {
        std::cout << "N=" << sample.first << ", divisors=[";
        for (std::size_t j = 0; j < sample.second.size(); ++j) {
            if (j != 0) {
                std::cout << ", ";
            }
            std::cout << sample.second[j];
        }
        std::cout << "] -> " << count_union_mod(sample.first, sample.second) << '\n';
    }
}

import java.math.BigInteger;
import java.util.Arrays;

public class Main {
    private static final BigInteger MOD = BigInteger.valueOf(1_000_000_007L);

    private static long gcdPositive(long a, long b) {
        while (b != 0) {
            long next = a % b;
            a = b;
            b = next;
        }
        return a;
    }

    public static BigInteger countUnionMod(long n, long[] divisors) {
        if (n < 0) {
            throw new IllegalArgumentException("N must be nonnegative");
        }
        if (divisors.length > 20) {
            throw new IllegalArgumentException("At most 20 divisors are supported");
        }
        for (long d : divisors) {
            if (d <= 0) {
                throw new IllegalArgumentException("Divisors must be positive");
            }
        }
        int k = divisors.length;
        if (n == 0 || k == 0) {
            return BigInteger.ZERO;
        }

        BigInteger total = BigInteger.ZERO;
        int limit = 1 << k;
        for (int mask = 1; mask < limit; ++mask) {
            long lcmValue = 1;
            int selected = 0;
            boolean tooLarge = false;
            for (int j = 0; j < k; ++j) {
                if ((mask & (1 << j)) != 0) {
                    ++selected;
                    long d = divisors[j];
                    long q = lcmValue / gcdPositive(lcmValue, d);
                    if (q > n / d) {
                        tooLarge = true;
                        break;
                    }
                    lcmValue = q * d;
                }
            }
            if (tooLarge) {
                continue;
            }
            BigInteger term = BigInteger.valueOf(n / lcmValue).mod(MOD);
            if (selected % 2 == 1) {
                total = total.add(term).mod(MOD);
            } else {
                total = total.subtract(term).mod(MOD);
            }
        }
        return total;
    }

    public static void main(String[] args) {
        long[] limits = {30, 12, 10, 5, 0, 10, 100};
        long[][] lists = {
            {2, 3, 5}, {2, 4}, {2, 2}, {6, 10}, {2, 3}, {}, {1, 7}
        };
        System.out.println("modulus = " + MOD);
        for (int i = 0; i < limits.length; ++i) {
            System.out.println("N=" + limits[i] + ", divisors="
                + Arrays.toString(lists[i]) + " -> "
                + countUnionMod(limits[i], lists[i]));
        }
    }
}

Run Python with Python 3.9 or later. Compile C++ with C++14 or later. Save Java as Main.java, compile with javac Main.java, and run with java Main. Each program prints these exact eight lines:

• modulus = 1000000007
• N=30, divisors=[2, 3, 5] -> 22
• N=12, divisors=[2, 4] -> 6
• N=10, divisors=[2, 2] -> 5
• N=5, divisors=[6, 10] -> 0
• N=0, divisors=[2, 3] -> 0
• N=10, divisors=[] -> 0
• N=100, divisors=[1, 7] -> 100

As a separate modular check, N = 2,000,000,020 and divisors [2] give a full count of 1,000,000,010 and a returned residue of 3. To compute full counts, Python and Java can omit modular reduction while retaining arbitrary-precision totals. C++ needs an arbitrary-precision accumulator or a separate bound on all partial sums: even when the final union count fits, intermediate signed sums can be larger.

Complexity, Correctness, and Invariants

There are 2ᵏ − 1 nonempty subsets. The straightforward implementation scans up to k positions per subset and keeps only a mask, an LCM, a selected count, and an accumulator. Its worst-case cost is O(k × 2ᵏ) bounded-width arithmetic/GCD operations and O(1) auxiliary space, excluding the input list. It does not store all intersections or enumerate all N objects.

The distinction between a GCD operation and an elementary arithmetic operation matters. Euclid's algorithm can take O(log D) remainder steps, where D bounds the positive divisors and intermediate LCMs. Counting those steps gives an O(k × 2ᵏ × log D) arithmetic-step upper bound. Arbitrary-precision inputs require a further bit-cost analysis; the fixed-width model describes the common C++/Java domain here.

When no LCM cutoff fires, all masks scan all k positions, giving Θ(k × 2ᵏ) position checks. For positive N, inputs with every divisor larger than N allow each subset to stop at its first selected position and give Θ(2ᵏ) position checks overall. N = 0 returns after Θ(k) input validation, while an empty list takes Θ(1). There is no distribution-independent average-case bound sharper than the worst-case upper bound: average behavior depends on the distribution of divisors and how often intersections become empty.

For correctness, first establish the intersection invariant: after processing selected divisors within a mask, L is their LCM. It is updated by the identity lcm(a, b) = (a/gcd(a, b)) × b. Therefore floor(N/L) counts exactly the objects satisfying every selected divisibility rule. If a partial LCM exceeds N, all extensions of that selected subset also have empty intersections, so skipping the mask loses no contribution.

Next, establish the accumulator invariant: after processing a mask, total equals the signed sum of intersection counts for all masks processed so far, modulo M. At termination, it is the full Inclusion–Exclusion sum. To see why that sum counts the union, consider one object that belongs to exactly t condition sets. For t ≥ 1, its coefficient is C(t, 1) − C(t, 2) + … + (−1)ᵗ⁺¹C(t, t) = 1 − (1 − 1)ᵗ = 1. An object in no sets appears in no nonempty intersection and contributes zero.

The table's arithmetic-operation bounds treat each GCD as one operation; the logarithmic GCD-step factor applies where relevant. Hash-set costs assume expected constant-time hashing, and U in its space entry is the number of distinct union elements.

Variants, Pitfalls, and Optimizations

Use the right intersection formula. For divisibility rules, use the LCM, not the product. The product equals the LCM for pairwise coprime selected divisors, but fails when factors are shared. For general sets, there may be no divisibility formula at all: an intersection-counting procedure must come from the problem's own structure.

Reduce the number of relevant conditions. Remove duplicates. If a divisor a divides another divisor b, every multiple of b is already covered by a, so b can be removed when counting this union. Divisors larger than N describe empty sets and can also be removed. These reductions can cut k dramatically, although the supplied programs retain the original list to make the subset mechanics visible.

Reuse or prune partial intersections. Recursive enumeration can carry the current LCM and parity into child calls, avoiding repeated scans. Once a partial LCM exceeds N, it can prune all supersets containing those already-selected conditions. A subset-DP approach can build each mask's LCM from a smaller mask, trading O(2ᵏ) memory for fewer repeated GCD calculations.

Count the complement when appropriate. The number satisfying none of the conditions is N minus the union count. For modular output, calculate ((N mod M) − union_mod + M) mod M. For probability, the same alternating-intersection identity holds with probabilities in place of cardinalities; independence is not required. Exact rational or appropriate numerical arithmetic replaces the integer modulus.

Use truncation carefully. Bonferroni inequalities give an upper bound from an odd number of complete intersection-size levels and a lower bound from an even number. For example, summing all single sets is an upper bound; subtracting all pair intersections gives a lower bound. An arbitrary prefix of numeric masks is not such a truncation. For larger structural problems, Möbius inversion generalizes the same cancellation idea on partially ordered sets.

Practice Problems and Real-World Use

For each problem, name the universe, identify the unwanted or qualifying sets, and explain how to count every selected intersection. Only then decide whether Inclusion–Exclusion is efficient. A formula with 2ᵏ terms is useful when k is small or when symmetry lets you combine many equal terms; it is not automatically an efficient algorithm for arbitrary large collections of sets.

Begin by checking small inputs with direct enumeration. In the divisibility problem, loop over x from 1 through N and count x once if any divisor divides it. Comparing this simple reference with the subset method catches incorrect signs, wrong intersections, missing triple terms, and duplicate-condition mistakes.

In the bounded-allocation example, unrestricted nonnegative solutions number C(7, 2) = 21. For each variable that is at least 4, subtract C(3, 2) = 3 solutions after shifting that variable down by 4. There are three such violation sets, and no two can occur together because their sum would already exceed 5, giving 21 − 3 × 3 = 12.

For audience counts or reliability estimates, the identity is exact only when the supplied intersection statistics are exact and refer to the same universe. Independence can simplify an intersection probability, but it must come from the model; Inclusion–Exclusion itself does not imply independence. An overlap estimate remains an estimate after applying the formula.

Before finishing, verify empty input, N = 0, a divisor equal to 1, duplicate divisors, divisibility chains such as [2, 4], and every divisor exceeding N. Include a three-condition example with a nonempty triple intersection and a large-value test that exercises the pre-multiplication bound. These checks test the counting argument and its numerical implementation together.