A Holiday Brain Teaser: [Diophantine quadruples]
From The Sunday Times, 4th August 1957 [link]
The series below is peculiar in that if any two of the numbers are multiplied together and added to unity, the result is a perfect square:
1, 3, 8
Let A[4], be a fourth number in this series. There is a further series with the same characteristic:
8, B[2], 190, B[4]
in which not only is B[1], the same as A[3], but B[2] is also the same as A[4].
What is B[4]?
(Note: zero is excluded)
This is one of the occasional Holiday Brain Teasers published in The Sunday Times prior to the start of numbered Teasers in 1961. A prize of £5 was offered for a solution to the puzzle, and a further prize of £25 for finding A[5] or B[5].
This puzzle was originally published with no title.
[teaser-1957-08-04] [teaser-unnumbered]
Frits, 
GeoffR, and 1 other are discussing. 
S2T2 
Jim Randell 11:28 am on 17 October 2021 Permalink |
I am assuming we are looking for a sequence of increasing integers.
With computers we can easily calculate the next terms.
Even a simple program runs in less than a second:
from enigma import (irange, inf, is_square, printf) # find the smallest number which when multiplied by any of ns, and 1 # is added to the result, give a square number def extend(ns): for x in irange(1, inf): if all(is_square(x * n + 1) for n in ns): return x # extend the series 1, 3, 8, ... A4 = extend([1, 3, 8]) printf("A4 = {A4}") # extend the series 8, A4, 190, ... B4 = extend([8, A4, 190]) printf("B4 = {B4}")But we can be a little bit cleverer, and this gives us a program that runs in just a few milliseconds:
Run: [ @replit ]
from enigma import (irange, inf, div, is_square, printf) # find the smallest number which when multiplied by any of ns, and 1 # is added to the result, give a square number def extend(ns): n0 = ns[0] # look for squares, z^2, such that: z^2 = n0 * x + 1 for z in irange(1, inf): x = div(z * z - 1, n0) if x and all(is_square(n * x + 1) for n in ns[1:]): return x # extend the series 1, 3, 8, ... A4 = extend([1, 3, 8]) printf("A4 = {A4}") # extend the series 8, A4, 190, ... B4 = extend([8, A4, 190]) printf("B4 = {B4}")Solution: A[4] = 120. B[4] = 730236.
See OEIS A030063 [@oeis], which says that the sequence cannot be extended further (although there is a further term in rationals).
The following was published in the 18th August 1957 edition of The Sunday Times [link]:
LikeLike
Jim Randell 12:10 pm on 17 October 2021 Permalink |
Euler found an infinite family of such 4-tuples:
(Which means that any pair of values from one quadruple can be used to seed another).
The A sequence of the puzzle is given by: a = 1; b = 3; r = 2.
And the B sequence is given by: a = 8; b = 120; r = 31.
from enigma import (div, printf) def extend(ns): (a, b, c) = ns r = div(c - a - b, 2) if r is not None and r * r == a * b + 1: return 4 * r * (r + a) * (r + b) # extend the series 1, 3, 8, ... A4 = extend([1, 3, 8]) printf("A4 = {A4}") # extend the series 8, A4, 190, ... B4 = extend([8, A4, 190]) printf("B4 = {B4}")In 1969 it was shown that the {1, 3, 8, 120} quadruple cannot be extended with another value, and in 2016 (a mere 59 years after the puzzle was set) it was proved that no Diophantine quintuple exists (see: [ @wikipedia ]). So the £25 prize was safe.
LikeLike
Jim Randell 5:48 pm on 20 October 2021 Permalink |
And here is a solution using Pell’s equations [revised to use the pells.py module from the py-enigma-plus repository] (see: Teaser 2994), which means fewer candidate numbers are considered when looking for solutions.
The following Python program runs in 63ms. (Internal runtime is 334µs).
from enigma import (div, is_square, peek, subsets, diff, printf) import pells # in general to extend a pair (a, b) we need x, such that: # # a.x + 1 = A^2 # b.x + 1 = B^2 # # so, eliminating x: # # B^2 = (b/a)(A^2 - 1) + 1 # B^2 = (b/a)A^2 - (b/a) + 1 # B^2 - (b/a)A^2 = 1 - b/a # # which, if b/a is an integer, is a form of Pell's equation # extend diophantine tuple where b/a is an integer def extend(a, b, *rest): D = div(b, a) for (B, A) in pells.diop_quad(1, -D, 1 - D): x = div(A * A - 1, a) if x and all(is_square(x * n + 1) for n in rest): yield x # extend 1, 3, 8 ... A4 = peek(extend(1, 3, 8)) printf("A4 = {A4}") # choose a viable (a, b) pair ns = (8, A4, 190) (a, b) = peek((a, b) for (a, b) in subsets(sorted(ns), size=2) if b % a == 0) B4 = peek(extend(a, b, *diff(ns, [a, b]))) printf("B4 = {B4}")LikeLike
GeoffR 3:09 pm on 17 October 2021 Permalink |
The best run-time I could get was 12.81 sec on an I9 processor. The Geocode solver took much longer. I had to use the answer to fix variable upper bounds.
% A Solution in MiniZinc include "globals.mzn"; predicate is_sq(var int: y) = exists(z in 1..ceil(sqrt(int2float(ub(y))))) ( z*z = y); % 1st series A1, A2, A3, A4 int: A1 == 1; int: A2 == 3; int: A3 == 8; var 1..200: A4; constraint is_sq(A1 * A2 + 1) /\ is_sq(A1 * A3 + 1) /\ is_sq(A1 * A4 + 1) /\ is_sq(A2 * A3 + 1) /\ is_sq(A2 * A4 + 1) /\ is_sq(A3 * A4 + 1); % 2nd series B1, B2, B3, B4 int: B1 == 8; var 9..200: B2; int: B3 == 190; var 191..800000: B4; constraint B2 == A4; constraint is_sq(B1 * B2 + 1) /\ is_sq(B1 * B3 + 1) /\ is_sq(B1 * B4 + 1) /\ is_sq(B2 * B3 + 1) /\ is_sq(B2 * B4 + 1) /\ is_sq(B3 * B4 + 1); solve satisfy; output [ "A1, A2, A3, A4 = " ++ show([A1, A2, A3, A4]) ++ "\n" ++ "B1, B2, B3, B4 = " ++ show([B1, B2, B3, B4])]; % A1, A2, A3, A4 = [1, 3, 8, 120] % B1, B2, B3, B4 = [8, 120, 190, 730236]LikeLike
A Holiday Brain Teaser | PuzzlingInPython 7:48 pm on 17 October 2021 Permalink |
[…] From The Sunday Times, 4th August 1957 [link] […]
LikeLike
GeoffR 10:38 pm on 17 October 2021 Permalink |
A Python solution ran in 258 msec.
from enigma import is_square as is_sq # 1st sequence is 1, 3, 8, A4 for A4 in range (9, 200): if not is_sq(A4 * 1 + 1): continue if not is_sq(A4 * 3 + 1): continue if not is_sq(A4 * 8 + 1): continue print('A4 = ', A4) # 2nd sequence is 8, A4, 190, B4 for B4 in range(191, 800000): if not is_sq(B4 * 8 + 1): continue if not is_sq(B4 * A4 + 1):continue if not is_sq(B4 * 190 + 1): continue print('B4 = ', B4) # A4 = 120 # B4 = 730236LikeLike
Frits 12:33 am on 18 October 2021 Permalink |
from enigma import div, first, inf, irange, is_square as is_sq # sqrt(8 * 190 + 1) = 39 A4s = [A4 for i in range(3, 40) if is_sq((A4 := i * i - 1) * 8 + 1) and is_sq(A4 * 3 + 1)] # 1st sequence is 1, 3, 8, A4 for A4 in A4s: # 2nd sequence is 8, A4, 190, B4 # sqrt(190 * 190 + 1) = 190.xx B4 = first(B4 for i in irange(191, inf) if (B4 := div(i * i - 1, 190)) and is_sq(B4 * A4 + 1) and is_sq(B4 * 8 + 1) ) print(f"A4 = {A4}, B4 = {B4[0]}")LikeLike
GeoffR 9:17 am on 18 October 2021 Permalink |
@Frits:
Neat solution.
The combined use of
first()in an iterator andinffrom the enigma.py library solves the problem I had in fixing the upper bound in a search for value B4.LikeLike
Frits 9:57 am on 18 October 2021 Permalink |
@GeoffR: the other way of achieving it is to break out of a while True loop for the first correct B4.
Also possible is using a wheel [2, 40, 152, 190] so that (i * i – 1) always is a multiple of 190.
B4 = first(B4 for k in irange(0, inf) for i in [189 + (190 * k) + x for x in [2, 40, 152, 190]] if is_sq((B4 := (i * i - 1) // 190) * A4 + 1) and is_sq(B4 * 8 + 1))LikeLike