[C/C++] 프로젝트 오일러 #21 10,000 이하의 친화수 찾기.

 

이번 문제는 밀레니엄 문제급은 아니지만, 홀수 완전수가 존재하지 않는다의 증명과 같이 오랫동안 풀리지 않는 문제의 영역입니다.  완전수에 대해서 아시는 분들은 많을거라고 생각합니다.  완전수는 자신을 제외한 약수의 합이 자신이 되는 수를 뜻합니다.  대표적으로는 6, 28 과 같이 짝수만 존재합니다.

 

완전수의 갯수가 무한한가라는 것도 사실 밝혀진 바가 없습니다.  짝수 완전수의 경우에는 반드시 이 형태여야 된다는 것은 이미 증명되었고, 그 증명과정이 어렵지는 않습니다.

 

완전수와 비슷하게 친화수도 자기 자신을 제외한 약수의 합을 통해서 나타냅니다.  그 약수들의 합이 다른 수가 되고, 그 다른 수의 약수들의 합이 자신이 되면, 두 사이 관계는 친화수라고 합니다.  친화수가 수학에서 큰 비중을 가지는 수는 아닙니다.  (완전수 역시 비슷하죠.)  역사적으로는 기원전 피타고라스 학파때로 올라갑니다.  그 유명한 피타고라스 정리가 나왔죠.  그 외에도 피타고라스 학파는 상당히 많은 수학적 성과를 이룩했는데요.  사실, 전 개인적으로 피타고라스 학파를 좋아하는 편은 아닙니다.  (너무 종교적이었고, 꽉 막혀 있는 집단이었다고 생각합니다.  그 시대를 살지 않았으니..)

 

 

 

 

요즘이야 컴퓨터로 계산하면 1초도 안되어서 몇만, 몇십만의 친화수를 계산할 수 있지만요, 당시만 해도 손으로 친화수를 계산한다는 것은 어려웠을겁니다.  친화수 중 몇가지는 공식에 의해서 만들어집니다.  그러나 그 공식으로 도출할 수 있는 친화수도 몇가지 안 됩니다.

 

완전수의 갯수는 겨우 몇십개입니다.  새로운 완전수를 찾는다면, 그 사람은 완전수 발견자이기에 앞서서, 메르센 소수 발견자로 수학사에 남게되겠죠.  홀수 완전수의 증명이 되지 않는 이상은 현재까지는 완전수 발견은 메르센 소수 발견자와 동일합니다.

 

가끔은 비트코인 욕하면서, 비트코인 마이닝 작업에 메르센소수 찾기를 했다면, 메르센소수를 더 찾았을지도 모른다고 생각합니다.  잡설이고요.. ^^

 

친화수를 찾는 로직자체는 그다지 다른 사람과 다르지 않습니다.  단지 소수를 찾는 과정이 추가된 것과 약수의 합을 계산하는 과정이 다릅니다.

 

어떤 수의 약수의 합은 다음과 같이 표현할 수 있습니다.

\[ n = \prod_p p^{a_p} \]

일 때,

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

 

 

이 방식을 따르기 위해서는 모든 소수를 가지고 계산을 해야했습니다.  그래서 소수 리스트가 필요했고요.

 

소수 리스트를 인수로 넘겨줄 수도 있었지만, 전역변수로 사용을 해서 프로그램을 작성했습니다.  그리고 소수 리스트를 이용해서 소인수 분해를 하고, 약수들의 합을 구했습니다.

 

그래서 만들어진 소스는 다음과 같습니다.

 

#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;
    }
}