[C/C++] Project Euler #21 Amicable Numbers

This problem, while not at the level of a Millennium Problem, belongs to the realm of long-standing unsolved mathematical challenges, similar to the proof that no odd perfect number exists. Many people might already be familiar with the concept of a perfect number. A perfect number is defined as a number that equals the sum of its proper divisors (excluding itself). Examples include 6 and 28, which are all even.

It is still unknown whether the set of perfect numbers is infinite. In the case of even perfect numbers, it has already been proven that they must follow a specific form, and the process of proving this is not particularly difficult.

Similar to perfect numbers, there are also amicable numbers. These numbers are defined through the sum of proper divisors. If the sum of the proper divisors of one number equals another number, and vice versa (the sum of the second number’s proper divisors equals the first number), the two numbers are said to be amicable. Amicable numbers do not hold a significant position in mathematics (neither do perfect numbers, for that matter). Historically, they date back to the Pythagorean school in ancient times. This is the same school known for the famous Pythagorean theorem. In addition to that, the Pythagoreans achieved numerous mathematical milestones. Personally, I don’t have much admiration for the Pythagorean school (as I find it to have been overly religious and rigid, but I wasn’t there to witness it firsthand, so who knows…).

These days, with computers, it takes less than a second to calculate tens or even hundreds of thousands of amicable numbers. However, back in the day, calculating amicable numbers by hand must have been extremely challenging. Some amicable numbers can be generated using formulas, but even those formulas only yield a limited number of examples.

The number of known perfect numbers is only a few dozen. If someone were to discover a new perfect number, they would go down in mathematical history as the discoverer of a Mersenne prime first and foremost, as finding a new perfect number is currently equivalent to finding a new Mersenne prime unless the proof for the existence of odd perfect numbers is established.

Sometimes I joke about how, instead of mining Bitcoin, if all that computational effort had been directed toward finding Mersenne primes, we might have discovered more of them. But that’s just a random thought. 😊

The logic for finding amicable numbers isn’t much different from what others use. It just involves an additional step of finding prime numbers and calculating the sum of divisors.

The sum of divisors of a number can be expressed as follows:

If
\[ n = \prod_p p^{a_p} \]
then
\[ \sigma(n) = \prod_p \frac{p^{a_p + 1 } - 1}{p-1} \]

To use this method, one must calculate with all primes, which necessitated a list of prime numbers.

Although it was possible to pass the prime list as an argument, I opted to use it as a global variable when writing the program. Using the prime list, I performed prime factorization and calculated the sum of divisors.

 

amicable numbers

The resulting code is as follows.

//------------------------------------------------
//    Project Euler #21 - Amicable Numbers
//        - by Aubrey Choi
//        - created at 2015-01-17
//------------------------------------------------
#include <stdio.h>

int sigma(int n);
void getprimes(int limit);

int primes[10000];
int pcount = 0;

int main()
{
    int sum = 0;

    getprimes(10000);
    for( int i = 1 ; i <= 10000 ; i++ )
    {
        int s = sigma(i)-i;
        if( s <= i || s > 10000 ) continue;
        if( i != sigma(s)-s ) continue;
        sum += s + i;
        printf("(%d, %d)\n", i, s);
    }
    printf("Ans = %d\n", sum);
}

int sigma(int n)
{
    int sum = 1;
    int t = n;
    for( int i = 0 ; t > 1 && i < pcount ; i++ )
    {
        int r = 1;
        while( t%primes[i] == 0 ) { r = r*primes[i] + 1; t /= primes[i]; }
        sum *= r;
    }
    return sum;
}

void getprimes(int limit)
{
    pcount = 0;
    primes[pcount++] = 2;
    for( int i = 3 ; i < limit ; i += 2 )
    {
        int j;
        for( j = 1 ; j < pcount ; j++ )
            if( i%primes[j] == 0 ) break;
        if( j == pcount ) 
            primes[pcount++] = i;
    }
}

 

Main Components:
1. sigma Function
• Calculates the sum of all divisors of a given number n.
• It uses prime factorization with the list of primes (primes[]) and applies the divisor sum formula:

\[ \sigma(n) = \prod_p \frac{p^{a_p + 1} - 1}{p - 1} \]

• The function iteratively factors n using primes from the primes[] array.
2. getprimes Function
• Generates all prime numbers up to 10,000 and stores them in the primes array.
• It uses a method similar to the Sieve of Eratosthenes to efficiently find primes.
3. Main Function
• Iterates through numbers from 1 to 10,000, computing the sum of proper divisors of each number using the sigma function.
• For each number i, it calculates s = sigma(i) - i, the sum of proper divisors excluding i.
• If s is within the valid range and sigma(s) - s equals i, then i and s are amicable numbers.
• Adds i and s to the total sum if they form an amicable pair.

Notes:
• To avoid double-counting amicable pairs (e.g., (i, s) and (s, i)), only pairs where s > i are considered.
• The final output is the sum of all amicable numbers below 10,000.