Project Euler 問題9（特別なピタゴラス数）を解いてみる

Project Euler で紹介されている問題9を解いてみました。素直なプログラム例とピタゴラス数が満たす性質を使った解き方を紹介します。

問題9（特別なピタゴラス数）

解答例

方針としては、$a<b<c$ かつ $a+b+c$ を満たす $a$、$b$、$c$ を生成して、それが三平方の式を満たすか調べることにします。

検算をするために解答に加えて、求めたピタゴラス数も出力しています。C 言語での実装例です。

#include <stdio.h>

#define SUM 1000

int main(void)
{
	int a, b, c;

	for (a = 1; a <= (SUM/3); ++a) {
		for (b = a + 1; b <= (SUM/2); ++b) {
			c = SUM - (a + b);
			if (b < c) {
				if (a * a + b * b == c * c) {
					printf("Problem009: %d\n", a * b * c);
					printf("Problem009: a:%d b:%d c:%d\n", a, b, c);
				}
			}
		}
	}

	return 0;
}

プログラムは、条件をみたすピタゴラス数をすべて出力します。プログラムを動作させると、問題に記載されている通り、1通りしか解が出力されないことも確認できます。

考察

この問題は、工夫すれば、手計算でも解を求めることができます。

原始ピタゴラス数

$a$、$b$、$c$ の最大公約数が 1 である、つまり $a$、$b$、$c$ が互いに素であるピタゴラス数を原始ピタゴラス数（primitive Pythagorean）と呼びます。

例えば、$3$、$4$、$5$ は原始ピタゴラス数です。明らかに原始ピタゴラス数を $k$ 倍した数もピタゴラス数になります。$3$、$4$、$5$ を2倍した、$6$、$8$、$10$ は、$6^2+9^2={10}^2$ を満たします。

次の事実（定理）が古代ギリシャですでに知られていました。証明は省略します。

問題に適用する

つまり、すべてのピタゴラス数は、原始ピタゴラス数を $k$ 倍（$k＝1$ の場合も含む）した数となるため、すべてのピタゴラス数は、以下の形になります。

$$a=k(m^2–n^2),b=2kmn,c=k(m^2+n^2)$$

これから $a+b+c$ を計算すると、

$a+b+c=k(m^2–n^2)+2kmn+k(m^2+n^2)$
$=2km(m+n)$

一方、$a+b+c=1000=2\times 2\times 2\times 5\times 5\times 5$ となるため

$$500=2\times 2\times 5\times 5\times 5=km(m+n)$$

ここで、$m+n$ は、明らかに $m$ より大きいため、素因数から、考えられる $m$ と $m+n$ をすべて抜き出して、$m>n$ を満たすか確認します。

解が2つでたように表から読み取れますが、$m=4=2\times 2$、$m+n=5$ の場合、$k=25=5\times 5$となります。もう一方の$m=20=2\times 2\times 5$、$m+n=25=5\times 5$ となるため、$k=1$ となります。この解は、以下を計算して一致していることが分かります（この問題では、$a<b$ となることを求めているため、気になるかたは、以下の $a$ と $b$ を入れ替えて読んでください）。

$$a=25\times (4^2－1^2),b=25\times 2\times 4\times 1,c=25\times (4^2+1^2)$$

$$a={20}^2－5^2,b=2\times 20\times 5,c={20}^2+5^2$$

ピタゴラス数がもつ性質によって、解答を求めることができました。

解答案（改善版）

考察したロジックをプログラムに落とし込みます。手計算では素因数から候補を求めましたが、プログラムの場合は、2重ループで、$k$ の値が自然数になるか確認しています。プログラムは、$m$、$n$、$k$ も出力しています。

#include <stdio.h>

#define SUM 1000

int main(void)
{
	int m, n, k;
	int a, b, c;

	for (n = 1; n * n <= SUM; ++n) {
		for (m = n + 1; m * m <= SUM; ++m) {
			if ((SUM / 2)%(m *(m + n)) == 0) {
				k = (SUM / 2)/(m *(m + n));
				if ((m * m - n * n)<= 2 * m * n) {
					a = k *(m * m - n * n);
					b = k * 2 * m * n;
				} else {
					a = k * 2 * m * n;
					b = k *(m * m - n * n);
				}
				c = k *(m * m + n * n);
				printf("Problem009: %d\n", a * b * c);
				printf("Problem009: a:%d b:%d c:%d\n", a, b, c);
				printf("Problem009: m:%d n:%d k:%d\n", m, n, k);
			}
		}
	}

	return 0;
}

素直なプログラムは、1000に対して、1000/2、1000/3 まで調べています（工夫すれば、もう少し小さくできます）。一方、考察に従ったプログラムは、$\sqrt{1000}$ まで調べることになるため、和の数（問題9の場合、1000）が大きくなっていくと、考察に従ったプログラムの方が速く動作するようになります。

最後に

素直なプログラムを提示したあと、ピタゴラス数を生成する式から解答を求めてみました。

今回の問題を解くためには使いませんでしたが、高校数学の数学A「整数の性質」の範囲で、余り（合同式）を使う証明として、以下の例がよくでてきます（a、b、c はピタゴラス数とします）。

• a、b のうち、少なくともひとつは、3の倍数である。
• a、b のうち、少なくともひとつは、4の倍数である。
• a、b、c のうち、少なくともひとつは、5の倍数である。

ピタゴラス数は、昔から知られていますが、興味深い性質が多くあるようです。

引き続き、Project Euler の問題を紹介していきます。