#include #include #include #include #include using namespace std; // Enumerates (n, A004114(n)) for 0 < n < count. // Relies on the fact that 2(A004114(n))^2 = A081324(n + 1) for n >= 1 void a004114(const uint64_t count) { uint64_t n = 1; // Index in Sequence uint64_t a_n; // A081324(n+1) uint64_t e; // sqrt(a_n/2) uint64_t iSqr, jSqr; // Tests for iSqr + jSqr = a_n. uint64_t i; // sqrt(iSqr) uint64_t j; // sqrt(jSqr) bool isSquareSum; // Is 2e^2 a sum of two perfect squares i^2 and j^2 such that i != j? // Skips A081324(1) = 0 by initializing e to 1. for (e = 1; n <= count; e++) { isSquareSum = false; a_n = 2 * e * e; // For counting the number of integer pairs (i, j) such that i^2 + j^2 == 2e^2 // and such that i > j (to avoid repeats), i > e. Furthermore, i^2 != 2e^2, as // 2e^2 cannot be a perfect square, and i^2 <= 2e^2, as i^2 > 2e^2 would imply // that j^2 < 0, which is obviously impossible, so i^2 < 2e^2. for (i = e + 1, iSqr = i * i; iSqr < a_n; i++, iSqr = i * i) { jSqr = a_n - iSqr; // Solves for j^2. j = (uint64_t)(sqrtl((long double)jSqr) + 0.5); // Calculates j rounded to the nearest integer. if (j * j == jSqr) { // If jSqr is a perfect square isSquareSum = true; break; } } if (!isSquareSum) { cout << n << ' ' << e << endl; n++; } } }