Project Euler - Problem 279

Triangles with integral sides and an integral angle

辺の長さが整数で, 少なくとも 1 つの角が整数(度で計測)な三角形のうち, 周長が 10⁸ を超えないものはいくつあるか.

(日本語訳より)

f:id:my316g:20180207070403p:plain

少なくとも 1 つの辺が整数角 x° である図のような三角形を考えると，余弦定理より，

$\displaystyle \cos \left( \frac{x}{180} \pi \right) = \frac{ a^{2} + b^{2} - c^{2} }{ 2 a b }$

各辺 a, b, c は整数だから，$\displaystyle \cos \left(\frac{x}{180} \pi \right)$ は有理数．

$\displaystyle \cos \left(\frac{x}{180} \pi \right)$ が有理数となるような整数 0 < x < 180 は，実のところ x = 60, 90, 120 の 3 つしかないことがわかる．*1

x = 90 の時は普通にピタゴラス数をの総数を求めればいいので，よくある単位円周上の有理点を求める方法を使う．*2

残りの x = 60, 120 の場合も実は似たような方法が使えて，楕円 $x^{2} \pm xy + y^{2} = 1$ 上の有理点を求める問題になる．

重複しないようにうまく数えるのが大変．

#include "bits/stdc++.h"
using namespace std;
typedef long long ll; 

template<typename T> T GCD(T m, T n) { T k; do { k = m % n; m = n; n = k; } while (k != 0); return m; }

int main()
{
    clock_t start, end;
    start = clock();

    const ll MAX = 100000000;

    ll ans = 0;

    // 90 deg
    for (ll m = 1; 4 * m * m <= MAX; m++) {
        for (ll n = m + 1; 2 * n * (m + n) <= MAX; n += 2) {
            if (GCD(m, n) != 1) continue;
            ll a = 2 * m * n;
            ll b = n * n - m * m;
            ll c = n * n + m * m;
            ll L = a + b + c;
            ans += MAX / L;
        }
    }

    // 60 deg
    for (ll m = 1; 2 * m * m <= 3 * MAX; m++) {
        for (ll n = 2 * m; (m + n) * (2 * n - m) <= 3 * MAX; n++) {
            if (GCD(m, n) != 1) continue;
            ll a = 2 * m * n - m * m;
            ll b = n * n - m * m;
            ll c = n * n - m * n + m * m;
            ll d = GCD(a, GCD(b, c));
            if (d > 1) { a /= d; b /= d; c /= d; }
            ll L = a + b + c;
            if (L <= MAX) {
                ans += MAX / L;
            }
        }
    }

    // 120 deg
    for (ll m = 1; 2 * m * m <= MAX; m++) {
        for (ll n = 2 * m; (m + n) * (2 * n + m) <= 3 * MAX; n++) {
            if (GCD(m, n) != 1) continue;
            ll a = 2 * m * n + m * m;
            ll b = n * n - m * m;
            if (a > b) continue;
            ll c = n * n + m * n + m * m;
            ll d = GCD(a, GCD(b, c));
            if (d > 1) { a /= d; b /= d; c /= d; }
            ll L = a + b + c;
            if (L <= MAX) {
                ans += MAX / L;
            }
        }
    }

    cout << ans << endl;

    end = clock();
    printf("%d msec.\n", end - start);

    return 0;
}