123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081 |
- #!/usr/bin/ruby
- # Find all the primitive solutions to the inverse Pythagorean equation:
- # 1/x^2 + 1/y^2 = 1/z^2
- # such that x <= y and 1 <= x,y,z <= N.
- # It can be shown that all the primitive solutions are generated from the following parametric form:
- #
- # x = 2*a*b*(a^2 + b^2)
- # y = a^4 - b^4
- # z = 2*a*b*(a^2 - b^2)
- #
- # where gcd(a, b) = 1 and a+b is odd.
- # See also:
- # https://oeis.org/A341990
- # https://math.stackexchange.com/questions/2688808/diophantine-equation-of-three-variables
- func S(N) {
- var solutions = []
- var limit = N.iroot(4)
- for a in (1..limit), b in (1 ..^ a) {
- is_odd(a+b) || next
- is_coprime(a,b) || next
- var x = (2*a*b * (a**2 + b**2))
- var y = (a**4 - b**4)
- var z = (2*a*b * (a**2 - b**2))
- x <= N || next
- y <= N || next
- z <= N || next
- solutions << [x,y,z]
- assert_eq(1/x**2 + 1/y**2, 1/z**2)
- }
- solutions.sort
- }
- var N = 1e4
- say <<"EOT"
- :: Primitve solutions (x,y,z) with 1 <= x,y,z <= #{N} and x <= y, to equation:
- 1/x^2 + 1/y^2 = 1/z^2
- EOT
- for x,y,z in (S(N)) {
- say "(#{x}, #{y}, #{z})"
- }
- __END__
- :: Primitve solutions (x,y,z) with 1 <= x,y,z <= 10000 and x <= y, to equation:
- 1/x^2 + 1/y^2 = 1/z^2
- (20, 15, 12)
- (136, 255, 120)
- (156, 65, 60)
- (444, 1295, 420)
- (580, 609, 420)
- (600, 175, 168)
- (1040, 4095, 1008)
- (1484, 2385, 1260)
- (1640, 369, 360)
- (2020, 9999, 1980)
- (3060, 6545, 2772)
- (3504, 4015, 2640)
- (3640, 2145, 1848)
- (3660, 671, 660)
- (6540, 9919, 5460)
- (6984, 6305, 4680)
- (7120, 3471, 3120)
- (7140, 1105, 1092)
|