ϐλΰϥεΛܭࢉ͠Α͏ɿࠜੑฤ ࠜੑͰܭࢉͩʂ from itertools import product for a, b, c in product(range(1, 20), repeat=3): if a ** 2 + b ** 2 == c ** 2: print(a, b, c) Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 10 / 39
طϐλΰϥεΛܭࢉ͠Α͏ ͪΐͬͱͨ͠ from itertools import combinations from math import gcd for a, b, c in combinations(range(1, 50), 3): if a ** 2 + b ** 2 == c ** 2 and gcd(a, b) == 1: print(a, b, c) ͪΐͬͱͨ͠ › a; b; c ͦΕͧΕҧ͏Ͱ͋Δ › a; b ޓ͍ʹૉͰ͋Δ Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 13 / 39
طϐλΰϥεΛͬͱ্ख͘ܭࢉ͠Α͏ ิ (͋ͱͰͬ͘͡Γಡ͏) طϐλΰϥε (a; b; c) ͷ a; b ͷҰํۮͰ͏ҰํحͰ͋Δɻ ·ͨɺc ৗʹحͰ͋Δɻ ূ໌ʢ͋ͱͰͬ͘͡Γಡ͏ʣ. ࣗવͷ (a; b; c) Λطϐλΰϥεͱ͢Δɻ ·ͣɺa; b ͕ڞʹۮͳΒɺϐλΰϥεͷఆ͔ٛΒ c ۮͱͳΓɺ (a; b; c) ڞ௨Ҽ 2 Λ࣋ͭ͜ͱʹͳΓԾఆʹ͢Δɻ࣍ʹɺa; b ͕ڞ ʹحͳΒɺϐλΰϥεͷఆ͔ٛΒ c ۮͱͳΔɻ͜ͷͱ͖ɺ a = 2x + 1; b = 2y + 1; c = 2z ͱ͢Δͱ a2 + b2 = c2 ) 2(x2 + x + y2 + y) + 1 = 2z2 ͱͳΓໃ६ɻ Αͬͯɺa; b ͷɺҰํۮͰ͏ҰํحͰ͋Δɻ·ͨɺϐλΰϥ εͷఆ͔ٛΒ c حͱͳΔɻ Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 16 / 39
طϐλΰϥεΛͬͱ্ख͘ܭࢉ͠Α͏ ্खͨ͘͠ from itertools import combinations from math import gcd for t, s in filter( lambda x: gcd(*x) == 1, combinations(range(1, 10, 2), r=2), ): print( s * t, (s ** 2 - t ** 2) // 2, (s ** 2 + t ** 2) // 2, ) Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 20 / 39
A2 + B2 = Mp ͳΔ (A; B; M) Λ୳͢ ߹ಉํఔࣜ x2 ” `1 (mod p) ͷղΛࢉग़͢ΔʢఱԼΓʣ from random import randint from sympy.ntheory import isprime, legendre_symbol def solve_quadratic(p): if not (isprime(p) and p % 4 == 1): raise ValueError while True: a = randint(1, p - 1) b = pow(a, (p - 1) // 4, p) if legendre_symbol(a, p) == -1 and 1 <= b < p: return b Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 29 / 39
A2 + B2 = Mp ͳΔ (A; B; M) Λ୳͢ (A; B; M) Λ୳͢ def sums_of_two_squares(p): """4 Λ๏ͱͯ͠ 1 ͱ߹ಉͳૉΛ 2 ͭͷฏํͷͰද͢""" if not (isprime(p) and p % 4 == 1): raise ValueError A, B = solve_quadratic(p), 1 M = divmod(A ** 2 + B ** 2, p)[0] Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 30 / 39
a2 + b2 = Mr Λ୳͢ a2 + b2 = Mr Λ୳͢ ҎԼͷ݅Λຬͨ͢Α͏ʹ a; b ΛબͿɻ a ” A (mod M) b ” B (mod M) ` 1 2 M » a; b » 1 2 M r ͱ M ͷؔͷൿີ r = a2 + b2 M » “ M 2 ”2 + “M 2 ”2 M = M 2 < M: Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 31 / 39
a2 + b2 = Mr Λ୳͢ a2 + b2 = Mr Λ୳͢ from sympy.ntheory.modular import solve_congruence def sums_of_two_squares(p): """4 Λ๏ͱͯ͠ 1 ͱ߹ಉͳૉΛ 2 ͭͷฏํͷͰද͢""" if not (isprime(p) and p % 4 == 1): raise ValueError A, B = solve_quadratic(p), 1 M = divmod(A ** 2 + B ** 2, p)[0] a = solve_congruence((A, M), symmetric=True)[0] b = solve_congruence((B, M), symmetric=True)[0] r = divmod(a ** 2 + b ** 2, M)[0] Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 32 / 39
2 ͭͷฏํͷͰදͤΔૉ 2 ͭͷฏํͷͰදͤΔૉ def sums_of_two_squares(p): A, B = solve_quadratic(p), 1 M = divmod(A ** 2 + B ** 2, p)[0] while True: a = solve_congruence((A, M), symmetric=True)[0] b = solve_congruence((B, M), symmetric=True)[0] r = divmod(u ** 2 + v ** 2, M)[0] s = abs((a * A + b * B) // M) t = abs((b * A - a * B) // M) if r == 1: print(f"${s}^2 + {t}^2={p}$") return else: A, B, M = s, t, r Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 35 / 39
2 ͭͷฏํͷͰදͤΔૉ 2 ͭͷฏํͷͰදͤΔૉ from sympy import primerange primes_1_mod_4 = ( p for p in primerange(1, 100) if p % 4 == 1 ) for p in primes_1_mod_4: sums_of_two_squares(p) Hayao (Hiroshima 2019) Number Theory with Python October 12, 2019 36 / 39