http://projecteuler.net/index.php?section=problems&id=9
いつものようにピタゴラス数の式を使うと、
a + b + c = 2lm(m + n)
これが1000になる。約数を出すオブジェクトを作ればよいが、これがなかなか大変だった。
#include <cmath> #include "itertools.h" using namespace std; using namespace itertools; typedef tuple<int,int> Fact; class factorize : public cGenerator<Fact> { int n; Fact f; factorize *fs; bool first; public: factorize(int n) : n(n) { first = true; for(int p = 2; p * p <= n; p++) { int e = 0; while(n % p == 0) { n /= p; e++; } if(e > 0) { f = Fact(p, e); if(n > 1) fs = new factorize(n); else fs = NULL; return; } } f = Fact(n, 1); fs = NULL; } Fact next() { try { if(first) { first = false; return f; } else { if(fs == NULL) throw(1); return fs->next(); } } catch(int e) { delete fs; throw(1); } } factorize *copy() const { return new factorize(n); } }; typedef vector<Fact> Facts; ostream& operator<<(ostream& os, const Facts& f); class divisors : public cGenerator<tuple<Facts,Facts>> { Facts::iterator begin; Facts::iterator end; divisors *g; int p, e; int e_max; tuple<Facts,Facts> f; public: divisors(Facts::iterator b, Facts::iterator e) : begin(b), end(e) { if(begin == end) { g = NULL; } else { g = new divisors(begin + 1, end); p = fst(*begin); e_max = snd(*begin); } this->e = 0; } tuple<Facts,Facts> next() { if(g == NULL) { if(e == 0) { e++; return tuple<Facts,Facts>(Facts(),Facts()); } else { throw(1); } } else { if(e == 0) f = g->next(); else if(e == 1) get<0>(f).insert(fst(f).begin(), Fact(p, 1)); else get<1>(fst(f).front()) = e; if(e == 0) get<1>(f).insert(snd(f).begin(), Fact(p, e_max)); else if(e == e_max) get<1>(f).erase(snd(f).begin()); else get<1>(snd(f).front()) = e_max - e; if(e == e_max) e = 0; else e++; return f; } } divisors *copy() const { return new divisors(begin, end); } }; int value_f(const Facts& f) { int v = 1; for(int i = 0; i < (int)f.size(); i++) { int p = get<0>(f[i]); int e_max = get<1>(f[i]); for(int e = 0; e < e_max; e++) v *= p; } return v; } typedef tuple<int,int,int> triplet; class divs : public cGenerator<triplet> { divisors *g1; divisors *g2; int n1; public: divs(int n) { Facts f = list(new factorize(n)); g1 = new divisors(f.begin(), f.end()); next_f2(); } divs(divisors *g1, divisors *g2, int n1) : g1(g1), g2(g2), n1(n1) { } tuple<int,int,int> next() { try { tuple<Facts,Facts> ff = g2->next(); return triplet(n1, value_f(fst(ff)), value_f(snd(ff))); } catch(int e) { delete g2; next_f2(); tuple<Facts,Facts> ff = g2->next(); return triplet(n1, value_f(fst(ff)), value_f(snd(ff))); } } divs *copy() const { return new divs(g1, g2, n1); } private: void next_f2() { tuple<Facts,Facts> ff = g1->next(); n1 = value_f(fst(ff)); Facts f2 = snd(ff); g2 = new divisors(f2.begin(), f2.end()); } }; template<class T> T gcd(T a, T b) { if(a % b == 0) return b; else return gcd(b, a % b); } bool is_valid(triplet d) { int l = get<0>(d); int m = get<1>(d); int n = get<2>(d) - m; return 0 < n && n < m && (m + n) % 2 == 1 && gcd(m, n) == 1; } int abc(triplet d) { int l = get<0>(d); int m = get<1>(d); int n = get<2>(d) - m; int a = l * (m * m - n * n); int b = 2 * l * m * n; int c = l * (m * m + n * n); return a * b * c; } const int N = 1000; int main() { cout << head(new map<int,triplet>(abc, new filter<triplet>(is_valid, new divs(N / 2)))) << endl; }