Problem 179: Consecutive Positive Divisors

View on Project Euler

Project Euler Problem 179 Solution

EulerSolve provides an optimized solution for Project Euler Problem 179, Consecutive Positive Divisors, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary We must count the integers \(n\) with \(1<n<10^7\) such that \(n\) and \(n+1\) have the same number of positive divisors. If \(d(m)\) denotes the divisor-count function, the required quantity is $$A=\sum_{n=2}^{10^7-1}\mathbf{1}_{d(n)=d(n+1)}.$$ Consecutive integers are always coprime, so this is not a question about two numbers sharing prime factors. Instead, it asks when two different factorizations happen to produce the same divisor count. The implementations therefore avoid factoring one pair at a time and compute the whole divisor-count table in a single structured pass. Mathematical Approach The problem is governed by the arithmetic function \(d(m)\), also written \(\tau(m)\). If $$m=\prod_{i=1}^{r} p_i^{a_i},$$ then the number of positive divisors is $$d(m)=\prod_{i=1}^{r}(a_i+1).$$ This explains why neighboring integers can still tie. For example, \(14=2^1\cdot 7^1\) and \(15=3^1\cdot 5^1\), so both have \(4\) divisors even though \(\gcd(14,15)=1\). The identity behind the sieve Factoring each integer separately would repeat the same work many times....

Detailed mathematical approach

Problem Summary

For a positive integer \(m\) let \(d(m)\) be the number of its positive divisors. We must count the integers \(n\) with \(1\lt n\lt 10^7\), that is \(2\le n\le 10^7-1\), for which \(n\) and \(n+1\) have the same number of positive divisors, \(d(n)=d(n+1)\). The statement gives one example: \(14\) has the divisors \(1,2,7,14\) and \(15\) has \(1,3,5,15\), so \(d(14)=d(15)=4\) and \(n=14\) is counted. [1]

The whole solution in one sentence: we let every number \(k\) add \(1\) to the table entries of all its multiples, which leaves every divisor count in the table after about \(1.6\cdot10^8\) additions, and then count the equal neighbouring entries.

Mathematical Approach

1. Small cases

The divisor counts of \(1,2,\dots,16\) are

$$\begin{aligned}&1,\ 2,\ 2,\ 3,\ 2,\ 4,\ 2,\ 4,\\&3,\ 4,\ 2,\ 6,\ 2,\ 4,\ 4,\ 5.\end{aligned}$$

Comparing neighbours, only \(d(2)=d(3)=2\) and \(d(14)=d(15)=4\) agree, so for \(1\lt n\lt 15\) the count is \(2\). Going on to \(99\), exactly \(15\) numbers \(n\) below \(100\) have a partner:

$$\begin{aligned}&2,\ 14,\ 21,\ 26,\ 33,\\&34,\ 38,\ 44,\ 57,\ 75,\\&85,\ 86,\ 93,\ 94,\ 98.\end{aligned}$$

The common value is \(2\) for \(n=2\), \(6\) for \(n=44\), \(75\) and \(98\), and \(4\) for the other eleven. Larger ranges give \(118\) matches below \(10^3\), \(1119\) below \(10^4\), \(10585\) below \(10^5\) and \(102093\) below \(10^6\), about one number in ten.

2. Divisor counts from the prime factorisation

By the fundamental theorem of arithmetic every \(m\ge2\) is, up to the order of the factors, a unique product \(m=p_1^{a_1}\cdots p_r^{a_r}\) of distinct primes. [3] If \(m=kj\), multiplying the factorisations of \(k\) and \(j\) gives the factorisation of \(m\), so a divisor \(k\) contains no prime other than \(p_1,\dots,p_r\), and \(p_i\) at most \(a_i\) times. Conversely, every product \(p_1^{b_1}\cdots p_r^{b_r}\) with \(0\le b_i\le a_i\) divides \(m\), and different choices of the exponents give different numbers, again by uniqueness. Each \(b_i\) can be chosen independently in \(a_i+1\) ways, so [2]

$$d(m)=(a_1+1)\cdots(a_r+1),$$

and \(d(1)=1\). For example \(14=2\cdot7\) and \(15=3\cdot5\) both have \((1+1)(1+1)=4\) divisors. So do \(26=2\cdot13\) and \(27=3^3\), two different shapes with the same product, while \(44=2^2\cdot11\) and \(45=3^2\cdot5\) share \((2+1)(1+1)=6\).

3. What a match looks like

Two neighbours \(n\) and \(n+1\) never share a prime factor, because a common divisor of both divides their difference \(1\). A match therefore means that the products of Section 2 agree although the two numbers have no prime factor in common, and two more facts explain the list of Section 1.

  • The common value is even. The divisors of \(m\) come in pairs \(k\), \(m/k\), which are different unless \(k^2=m\), so \(d(m)\) is odd exactly when \(m\) is a perfect square. [2] The squares \(k^2\) and \((k+1)^2\) with \(k\ge1\) differ by \(2k+1\ge3\), so \(n\) and \(n+1\) are never both squares, as an odd common value would require.
  • Only \(n=2\) has the common value \(2\). A number with two divisors is a prime, and for \(n\ge3\) one of \(n\), \(n+1\) is an even number greater than \(2\).

Matches can overlap: \(33=3\cdot11\), \(34=2\cdot17\) and \(35=5\cdot7\) all have \(4\) divisors and give the two matches \(n=33\) and \(n=34\), since each \(n\) is compared only with \(n+1\). Heath-Brown proved in 1984 that there are infinitely many matches, [2] and bounds of Hildebrand and of Erdős, Pomerance and Sárközy show that their proportion tends to \(0\), but extremely slowly. [4]

4. Turning the question around

Computing each \(d(m)\) on its own, for example by testing every \(k\le\sqrt m\) and counting the pairs \(k\), \(m/k\), costs about \(\sqrt m\) divisions per number, in total \(\sum_{m\le10^7}\lfloor\sqrt m\rfloor\)\({}\approx2.1\cdot10^{10}\). It is much cheaper to ask the opposite question: not which \(k\) divide \(m\), but which \(m\) a given \(k\) divides. The answer is the multiples \(k,2k,3k,\dots\), and up to a bound \(M\) there are \(\lfloor M/k\rfloor\) of them. Counting the pairs \((k,m)\) with \(k\mid m\) and \(m\le M\) once by \(m\) and once by \(k\) gives [5]

$$\sum_{m=1}^{M}d(m)=\sum_{k=1}^{M}\left\lfloor\frac{M}{k}\right\rfloor.$$

For \(M=16\) both sides are \(50\), the sum of the sixteen counts in Section 1.

5. The divisor sieve

The right-hand side is an algorithm. Start with a table \(c[1],\dots,c[M]\) of zeros and, for \(k=1,2,\dots,M\), add \(1\) to \(c[j]\) for every multiple \(j\) of \(k\) up to \(M\). The sieve of Eratosthenes walks through multiples in the same way, but only for primes, and crosses numbers off instead of counting. [6]

Claim. After the passes \(k=1,\dots,K\), the entry \(c[m]\) equals the number of divisors of \(m\) that are at most \(K\). In particular, after all \(M\) passes, \(c[m]=d(m)\) for every \(1\le m\le M\).

Proof. Pass \(k\) changes only the entries of multiples of \(k\), and it raises each of them by exactly \(1\), because the walk \(k,2k,3k,\dots\) meets each multiple once. So after the passes \(1,\dots,K\), the entry \(c[m]\) counts the \(k\le K\) with \(k\mid m\). Every divisor of \(m\) is at most \(m\le M\), so after the last pass all of them have been counted. \(\square\)

Pass \(k\) makes \(\lfloor M/k\rfloor\) additions, so the whole sieve makes exactly the right-hand side of Section 4.

6. Worked example: the limit 15

The programs build the table up to \(M=L+1\) for the limit \(L\), so for \(L=15\) they use \(M=16\). Pass \(1\) adds \(1\) to all \(16\) entries, pass \(2\) to the \(8\) even ones, pass \(3\) to \(3,6,9,12,15\), pass \(4\) to \(4,8,12,16\), pass \(5\) to \(5,10,15\), the passes \(6\), \(7\) and \(8\) to two entries each, and the passes \(9\) to \(16\) only to their own entry: \(50\) additions in all. Entry \(12\) is raised by the passes \(1,2,3,4,6,12\) and ends at \(6=d(12)\). Entry \(14\) is raised by \(1,2,7,14\) and entry \(15\) by \(1,3,5,15\), so both end at \(4\). The final table is the list of Section 1. The scan over \(n=2,\dots,14\) compares \(c[n]\) with \(c[n+1]\) and finds the matches \(n=2\) and \(n=14\), so the count for \(L=15\) is \(2\). The entry \(c[16]=5\) is computed but never read, since \(n+1\le L\).

7. The cost is a harmonic sum

Let \(U(M)=\sum_{k\le M}\lfloor M/k\rfloor\) be the number of additions. Since \(M/k-1\lt\lfloor M/k\rfloor\le M/k\), it lies between \(MH_M-M\) and \(MH_M\), where \(H_M=1+\frac12+\dots+\frac1M\) is the \(M\)-th harmonic number. Comparing \(1/k\) with the area under \(1/x\) from \(k\) to \(k+1\), or from \(k-1\) to \(k\), gives \(\ln(M+1)\le H_M\le1+\ln M\), [7] so

$$\begin{aligned}U(M)&\gt M\ln M-M,\\U(M)&\le M\ln M+M.\end{aligned}$$

The sieve thus needs about \(M\ln M\) additions, only a factor \(\ln M\approx16\) more than the \(M\) entries it fills. Dirichlet proved the sharper estimate [8]

$$\begin{aligned}U(M)={}&M\ln M+(2\gamma-1)M\\&+O\bigl(\sqrt M\bigr),\end{aligned}$$

where \(\gamma\approx0.5772\) is the Euler–Mascheroni constant and \(O(\sqrt M)\) stands for an error of at most a fixed constant times \(\sqrt M\). For the problem, \(M=10^7+1\) and the sieve makes exactly \(162\,725\,368\) additions, about \(16.3\) per entry, while Dirichlet's two main terms give about \(162\,725\,287\). That is about \(130\) times fewer steps than the \(2.1\cdot10^{10}\) trial divisions of Section 4.

8. How large the entries become

The table must hold the largest divisor count up to \(M\). Two numbers are decisive:

$$\begin{aligned}8648640&=2^6\cdot3^3\cdot5\cdot7\cdot11\cdot13,\\1081080&=2^3\cdot3^3\cdot5\cdot7\cdot11\cdot13.\end{aligned}$$

A scan of the full table shows that the largest count is \((6+1)(3+1)\cdot2^4=448\), reached only by \(8648640\). The number \(1081080\) has \(4\cdot4\cdot2^4=256\) divisors and is highly composite, which means that it has more divisors than every smaller number, [9] so it is the first number with more than \(255\) divisors. An \(8\)-bit counter, which holds at most \(255\), would therefore overflow there, and \(1374\) entries up to \(M\) need more than \(8\) bits. A \(16\)-bit counter, with room up to \(65535\), is ample.

9. Faster variants

Pairing each divisor \(k\le\sqrt m\) with \(m/k\) roughly halves the work, to \(81\,364\,265\) additions for the problem: only the passes \(k\le\sqrt M\) run, each starting at \(k^2\), where it adds \(1\), and adding \(2\) at every larger multiple. A linear sieve, a form of Euler's sieve that reaches every composite number exactly once, finds the smallest prime factor \(p\) of every \(m\le M\) in \(O(M)\) steps, [6] and Section 2 gives \(d(m)=(e+1)\,d(m/p^e)\) for the highest power \(p^e\) dividing \(m\). The programs keep the simple version; the linear sieve serves as an independent check.

How the Code Works

In the C++ program, [10] parse_arguments accepts --skip-checkpoints and --limit= with decimal digits, which parse_int_after_prefix reads into an int without an overflow check, so an oversized limit is undefined behaviour; [11] the default limit \(10^7\) is set in Options. Any other argument, or a limit below \(3\), gives exit code \(1\). Unless switched off, run_checkpoints requires solve(100) to return \(15\), the count of Section 1, with exit code \(2\) otherwise, and main then prints solve(options.limit).

solve creates div_count, a std::vector<std::uint16_t> of limit + 2 zeros for the indices \(0\) to limit + 1, where index \(0\) stays unused. The loop over d, the \(k\) of Section 5, from \(1\) to limit + 1, with the inner loop over m = d, 2d, … up to limit + 1 that executes ++div_count[m], is the sieve of Section 5 with \(M=L+1\). The final loop over n from \(2\) to limit - 1 counts the \(n\) with div_count[n] == div_count[n + 1]. As in Section 6 the last entry is never read; it costs only a handful of additions. By Section 8 the \(16\)-bit entries cannot overflow, but from a limit of about \(1.07\cdot10^9\) on, the int step m += d could.

The Python and Java programs fix the limit at \(10^7\), with no options or checkpoint, and run the same loops: Python in solve with a list div_count and range(d, limit + 2, d), returning the count as a string, and Java in main with an int[] array dc. The element type decides the size of the table: about \(20\) MB in C++ with \(2\) bytes per entry, \(40\) MB in Java with \(4\) bytes, and in 64-bit CPython about \(80\) MB for the \(8\)-byte references of the list alone.

Complexity and Verification

By Section 7 the sieve makes \(U(L+1)=O(L\log L)\) additions, \(162\,725\,368\) for \(L=10^7\), the scan makes \(L-2\) comparisons, and the memory is one table of \(L+2\) entries. Timed from start to finish on a laptop with an Apple M4 Max chip, a full run takes about \(0.1\) s in C++ (g++ -O2), \(0.2\) s in Java and \(11\) s in CPython, whose interpreter executes each addition as bytecode.

For this article the C++ program was also run with --limit=L for every \(L\) from \(3\) to \(3000\), and each result agreed with a direct count by trial division. An independent C program computed \(d(m)\) for all \(m\le10^7+1\) with the linear sieve of Section 9, agreed with the additive sieve entry by entry, confirmed the numbers of Sections 1, 7 and 8, and gave the same result for the problem as the three published programs. These checks are not part of the published solutions.

Footnotes and References

The references below supply the background facts used above; the sieve, its cost and the worked example are derived in this article.

  1. Project Euler 179 — Consecutive Positive Divisors. The official statement asks for the number of integers 1 < n < 10⁷ for which n and n + 1 have the same number of positive divisors, with the example of 14 (divisors 1, 2, 7, 14) and 15 (divisors 1, 3, 5, 15).
  2. Wikipedia — Divisor function. The number of divisors of p₁^a₁ ⋯ p_r^a_r is (a₁ + 1) ⋯ (a_r + 1). Divisors pair off, so only a square, whose root has no partner, has an odd number of them. In 1984 Heath-Brown proved that σ₀(n) = σ₀(n + 1) for infinitely many n.
  3. Wikipedia — Fundamental theorem of arithmetic. Every integer greater than 1 is a prime or a product of primes in exactly one way, up to the order of the factors, and multiplying two numbers adds the exponents of each prime.
  4. OEIS A005237 — Numbers k such that k and k+1 have the same number of divisors. The sequence begins 2, 14, 21, 26, 33, 34, 38; the entry notes that Heath-Brown proved it infinite and that Hildebrand and Erdős, Pomerance and Sárközy showed n·√(log log n) ≪ a(n) ≪ n·(log log n)³ for its n-th term.
  5. Wikipedia — Double counting (proof technique). Two expressions are equal when they count the same set in two ways; Section 4 counts the pairs (k, m) with k dividing m.
  6. Wikipedia — Sieve of Eratosthenes. The sieve marks the multiples of each prime and needs O(n log log n) operations; in Euler's sieve, described there too, each composite number is eliminated exactly once, which takes linear time.
  7. Wikipedia — Harmonic number. The article compares Hₙ with the area ln(n + 1) under 1/x from 1 to n + 1 and gives Hₙ = ln n + γ + 1/(2n) − ⋯ with γ ≈ 0.5772156649.
  8. Wikipedia — Divisor summatory function. D(x) = d(1) + ⋯ + d(⌊x⌋) is the sum of ⌊x/k⌋ over k ≤ x, and D(x) = x log x + x(2γ − 1) + O(√x), first established by Dirichlet in 1849.
  9. Wikipedia — Highly composite number. A highly composite number has more divisors than every smaller positive integer; the table lists 720720, with 240 divisors, followed by 1081080, with 256.
  10. C++ — Euler179.cpp. The linked immutable revision contains parse_arguments, solve and run_checkpoints. It is the implementation record for this article.
  11. cppreference — Arithmetic operators. When a signed integer operation overflows, its result not fitting in the type, the behaviour is undefined.

Problem 179 source code

C++

#include <cstdint>
#include <iostream>
#include <string>
#include <vector>

namespace {

struct Options {
    int limit = 10000000;
    bool run_checkpoints = true;
};

bool parse_int_after_prefix(const std::string& arg, const std::string& prefix, int& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }

    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    int parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<int>(c - '0');
    }

    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);
        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_int_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }

    return options.limit >= 3;
}

int solve(const int limit) {
    std::vector<std::uint16_t> div_count(static_cast<std::size_t>(limit + 2), 0);

    for (int d = 1; d <= limit + 1; ++d) {
        for (int m = d; m <= limit + 1; m += d) {
            ++div_count[static_cast<std::size_t>(m)];
        }
    }

    int answer = 0;
    for (int n = 2; n < limit; ++n) {
        if (div_count[static_cast<std::size_t>(n)] == div_count[static_cast<std::size_t>(n + 1)]) {
            ++answer;
        }
    }

    return answer;
}

bool run_checkpoints() {
    if (solve(100) != 15) {
        std::cerr << "Checkpoint failed for limit=100" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }

    std::cout << solve(options.limit) << '\n';
    return 0;
}

Python

def solve():
    limit = 10000000
    div_count = [0] * (limit + 2)
    for d in range(1, limit + 2):
        for m in range(d, limit + 2, d):
            div_count[m] += 1
    answer = 0
    for n in range(2, limit):
        if div_count[n] == div_count[n + 1]:
            answer += 1
    return str(answer)

if __name__ == '__main__':
    print(solve())

Java

public class Euler179 {
    public static void main(String[] args) {
        int limit = 10000000;
        int[] dc = new int[limit + 2];
        for (int d = 1; d <= limit + 1; d++)
            for (int m = d; m <= limit + 1; m += d)
                dc[m]++;
        int ans = 0;
        for (int n = 2; n < limit; n++)
            if (dc[n] == dc[n + 1])
                ans++;
        System.out.println(ans);
    }
}