-
Notifications
You must be signed in to change notification settings - Fork 0
Home
When I was 10-year-old young boy I read the book "Tic-Tac-Toe", Russian translation of "Wheels, life and other mathematical amusements" book of a brilliant American science popularizer Martin Gardner. I was particularly impressed by the problem of a perfect cuboid fascinatingly described by him. I was surprised that it has not been solved yet, while there are many known almost suitable solutions. 30 years have passed since then but the idea of finding the perfect cuboid still does not leave me.
So let's start from well-known for dedicated definitions, properties, equations and notations just to be on the same wavelength.
A perfect cuboid (also called a perfect Euler brick, a perfect box) is a rectangular cuboid whose 3 edges, 3 face diagonals and the body diagonal all have integer lengths. The existance of a perfect cuboid is one of unsolved problems in mathematics. The definition of a perfect cuboid in geometric terms is equivalent to a solution to the following system of Diophantine equations:
a^2 + b^2 = d^2
a^2 + c^2 = e^2
b^2 + c^2 = f^2
a^2 + b^2 + c^2 = g^2
A primitive perfect cuboid is a perfect cuboid whose edge lengths are relatively prime.
Some facts are known about properties that must be satisfied by a primitive perfect cuboid, if one exists, based on modular arithmetic:
- One edge, two face diagonals and the body diagonal must be odd, one edge and the remaining face diagonal must be divisible by 4, and the remaining edge must be divisible by 16.
- Two edges must have length divisible by 3 and at least one of those edges must have length divisible by 9.
- One edge must have length divisible by 5.
- One edge must have length divisible by 7.
- One edge must have length divisible by 11.
- One edge must have length divisible by 19.
- One edge or space diagonal must be divisible by 13.
- One edge, face diagonal or space diagonal must be divisible by 17.
- One edge, face diagonal or space diagonal must be divisible by 29.
- One edge, face diagonal or space diagonal must be divisible by 37.
Most computer attempts to find a perfect cuboid were directed to attack the problem from the side of known edges (smallest or odd). Exhaustive computer searches show that, if a perfect cuboid exists, as of May 2017,
- the odd edge must be greater than 2.5 · 10^13
- the smallest edge must be greater than 5 · 10^11
In 1992 Ivan Korec proposed to change the direction of attack from edges to the body diagonal. He used PASCAL, a software not designed for high-precision integers, so he reported that the body diagonal of a perfect cuboid must exceed 8 billion. In 2006 Terry Raines wrote his own UBASIC program and raised Korec's lower bound from 8 billion to 120 billion. In 2017 we wrote C (pure C, not C++) program based on Korec's ideas, and probably got rid of previous programs bottlenecks and therefore tens of thousands times more efficient, which allowed to raise the lower bound to 2^53, more than 9 · 10^15.
We used several important theorems, lemmas and formulas from Number Theory, which help us on our way to algorytm approach. We leave them without proof because they date back to Euler times and maybe even to ancient Greeks and need no derivation.
Girard's Criteria: A natural number can be represented as the sum of two squares if and only if its prime divisors of the form p(i) = 4k+3 are presented with even degrees.
Girard's Theorem: Primes p ≡ 1 mod 4 could be expressed as the sum of two squares in essentially only one way.
Dirichlet's formula: If n can be represented in the form of the sum of two squares, the amount of representations is equal to [(∏(a(j)+1)+1)/2].
If the amount of factors is equal to 0, then the product is assumed to be 1. Representations do not differ by the order of factors.
To simplify the following exposition let us call the representations of a natural number as the sum of two rational squares as Girard representations.
Main corollaries from the formula:
- all degrees of primes of the form 4k+3 must be even for n has Girard representations;
- the amount of Girard representations of n entirely depends on degrees of primes of the form 4k+1;
- the presence or absence 2 number in the factorization of n do not affect to the amount of Girard representations.
Korec's idea is based on crucial for computer search property of a primitive perfect cuboid:
Lemma 1: If p is a prime divisor of the body diagonal g then p ≡ 1 mod 4.
A trivial but important corollary is
Lemma 2: The body diagonal of a primitive perfect cuboid g ≡ 1 mod 4.
Lemma 3: if n is a product of two natural numbers: n = p(1)·p(2) and each of these numbers has its own Girard representation:
p(1) = a^2 + b^2
p(2) = c^2 + d^2
, the general formula for expressing a product of the sum of two squares in Girard representation form is the following:
(a^2 + b^2)(c^2 + d^2) = (|ac ± bd|)^2 + (|bc ∓ ad|)^2
Note, in general case pairs (a, b) and (c, d) can produce 2 different pairs of Girard representations: (ac + bd, |bc - ad|) and (|ac - bd|, bc + ad), if ac + bd ≠ bc + ad. Otherwise only 1 unique pair is produced.
Lemma 4: The body diagonal of a primitive perfect cuboid cannot be a prime number.
Let's return to Diophantine equations:
a^2 + b^2 = d^2
a^2 + c^2 = e^2
b^2 + c^2 = f^2
a^2 + b^2 + c^2 = g^2
Edges a, b, c must be different, so the square of the body diagonal g^2 must have at least 3 different Girard representations:
g^2 = a^2 + f^2
g^2 = b^2 + e^2
g^2 = c^2 + d^2
If g is a prime number of the form 4k+1, due to Dirichlet's formula g^2 has only 2 different representations, so prime g cannot be a candidate for the body diagonal of a perfect cuboid.
The program searches for a perfect cuboid in a given by lower and higher bounds range for the body diagonal g among the numbers of the form 4k+1. Giving a certain number g first we should check whether it is a product ∏ p(i) of prime divisors p(i) ≡ 1 mod 4. The existance of just one divisor p(i) ≡ 3 mod 4 leads to exclusion of g from the subset of candidates for the body diagonal of a perfect cuboid.
After the sieving (and, actually, the factorization) we got only right candidates. In general they take ~4% of numbers from the range, or ~16% if we take into account the fact that we move through the range with step 4. After that it does batch decomposition of factors of numbers-candidates into the sum of two squares and per each candidate the following routine is proceeded:
- combine all possible Girard representations of g^2 from known decompositions of all g factors;
- search for a perfect cuboid with the body diagonal g based on all known representations of g^2
We experimented with different methods of integer factorization (SQUFOF, Brent methods), but stopped at old-school Eratosthenes sieve. It does sieve and factorization simultaneously. The program received a powerful acceleration after the transition to the batch sieving. It means a given range is divided into fixed size smaller sections and every section entirely undergoes a procedure of sieving. An optimal size of the section is a subject of a separate study. Theoretically section size increasing should reduce program execution time since it needs less sieving loops, but in fact batch size growth causes also linear growth of memory usage and from a particular level of batch size the benefit of fewer sieving loops is leveling down by increasing requirements to consumed by the program memory. As far as we planned to publish single-core application, charge all available cores and do not confine the project by only those who has computers with tens of Gigabytes memory onboard, the task was do not have too high demands to the host's installed memory. So we chose 100000 as an optimal batch size.
We implemeted simplified Cornacchia's algorithm to compute x and y such that p = x^2 + y^2 when prime p ≡ 1 mod 4.
There are two stages of this algorithm:
- find z with z^2 ≡ -1 mod p
- use z to find x and y, where p = x^2 + y^2
If a is a quadratic nonresidue modulo p then a^(p-1)/2 ≡ -1 mod p so we can take z ≡ a^(p-1)/4 mod p (and that is easily computed by modular exponentiation). How does one find quadratic nonresidues? In fact, exactly half of the numbers 1, ..., p-1 have the desired property. We do not know a deterministic process for finding such a number, but they appear to be distributed randomly so we can just try numbers at random until we find one and each trial has a 50% chance of success. Even if we are trying 1000 digit numbers it is unlikely that we will not find a solution before the 5000th trial, so for computational purposes, it is probably worth keeping a list of the first 6542 primes. It is precisely so many primes under 2^16=65536.
Compute the greatest common divisor of p and z+i using the Euclidean algorithm for the Gaussian integers. The answer will be x+yi where x^2+y^2=p.
- set x <─ p, y <─ p mod z, s <─ [√p].
- while y>s, set r <─ x mod y, x <─ y, y <─ r.
Given g = p(1)^a(1)·p(2)^a(2)·...·p(n)^a(n), where p(i) are distinct primes and p(i) ≡ 1 mod 4. How many distinct primes g can have as its factors? It's not difficult to prove that 11 is the largest amount of different prime factors of the form p=4k+1 for numbers g ≡ 1 mod 4 under 2^63. For instance 8418894903232764925 = 5^2 · 13 · 17^2 · 29 · 37 · 41 · 53 · 61 · 73 · 89 · 97 has 11 different factors. Numbers 6437978455413290825^2 and 8418894903232764925^2 have the largest amount of Girard representations among g^2 for numbers g ≡ 1 mod 4 under 2^63. According to Dirihlet's formula: [((4+1)(4+1)(2+1)(2+1)(2+1)(2+1)(2+1)(2+1)(2+1)(2+1)(2+1)+1)/2] = 246038.
Probably we have to turn aside and explain why 2^63 was chosen as a higher limit for investigated numbers. The reason is very simple: as we know, most of modern computer CPUs are built on x86_64 architecture, the 64-bit version of the x86 instruction set, so 64-bit integer arithmetic and logical operations are their primary architecture features and they can operate directly on 64-bit integers. C compiler gcc supports also unsigned 128 bit integers. As we will see below, we shall perform square root computing of the sum of two squares, so 2^63 was chosen to prevent owerflow of 128 bit integers. In fact 2^63 should be enough for our purposes. With the current program speed they will need more than 1 million years on single CPU core to reach 2^63. Ok, now let's return to Girard representations.
Since we know the amount of factors of g and their unique and only one Girard representation, we will apply the following routine to construct all Girard representations of g^2:
- double degrees of g factors, since we are searching for Girard representations of g^2
- start the construction from the pair (0, 1)
- loop through factors
-
- loop through factor degrees
-
-
- construct representations of ∏ p(1)^a(1)·p(2)^a(2)·...·p(i)^a(j) based on representations from the previous step by adding Girard representation of p(i) with the help of formula (2).
-
-
-
- clean array from repetitions
-
- sort all representations by the smallest square (x<y).
- compute and store not only pairs (x, y) from g^2=x^2+y^2, but also values of x^2 and y^2: (x, x^2, y, y^2).
Forasmuch we had known the results of previous searches and did not expect to find a perfect cuboid in ranges at least up to 5 · 10^13, we decided to search not only a perfect cuboid, but also so-called almost-perfect cuboids.
An almost-perfect cuboid has one of the 7 lengths irrational, the other 6 lengths are rational. We had known about 3 types of almost-perfect cuboids, called Body, Edge and Face cuboids.
Body cuboids have irrational body diagonal and therefore were not achievable for our application because of chosen direction of attack.
For the Edge cuboid, one of the edges a, b, c is irrational. The Face cuboid has just one of the face diagonals d, e, f irrational. We had found that Edge and Face are completely achievable for the program. And, it is important, the search completely discovers them. It means any primitive Face and Edge cuboid, like a perfect cuboid, must have the body diagonal as a product of primes of the form 4k+1. We found who is engaged in such kind of research, he was Randall L. Rathbun. As of May 2017 he published 155,151 found by him cuboids with the smallest integer edge under 157,000,000,000. We compared his results with our own and found they are absolutely identical in ranges where they are comparable.
After that Randall inspired us to add also additional, more specific type of almost-perfect cuboids: a cuboid with irrational Gaussian edge. We found that we also can add such type of cuboids without affecting the performance of the program, but unfortunately their search will not exhaustive. For instance [(a, b, c), (d, e, f), g] = [(√-426400, 108, 725), (644i, 315, 733), 333] is not achievable for our program. In response we noticed that Randall's cuboid table is incomplete too, since we found a way to generate new ones cuboids with some lengthes in complex numbers from known one. For instance, Face cuboid [(a, b, c), (d, e, f), g] = [(104,153,672), (185,680,√474993), 697] can produce at least 11 new different cuboids:
- [(104i,185,680), (153,672,√496625), 697]
- [(697,153i,672i), (680,185,√-474993), 104]
- [(185,680,697i), (√496625,672i,153i), 104]
- [(697,104i,672i), (√474993,185,680i), 153]
- [(672,680i,185), (104i,697,√-428175), 153]
- [(697,104i,153i), (√474993,680,185i), 672]
- [(153,185i,680), (104i,697,√428175), 672]
- [(104,680i,697), (672i,√496625,153), 185]
- [(672i,680,153), (104,√-428175,697), 185]
- [(185i,104,697), (153i,672,√496625), 680]
- [(153i,185,672), (104,√428175,697), 680]
If we take into account the fact that we can multiply all lengthes of 12 above cuboids on the imaginary unit and get 12 new different cuboids:
- [(104i,153i,672i), (185i,680i,√474993), 697i]
- [(104,185i,680i), (153i,672i,√-496625), 697i]
- [(697i,153,672), (680i,185i,√474993), 104i]
- [(185i,680i,697), (√-496625,672,153), 104i]
- [(697i,104,672), (√-474993,185i,680), 153i]
- [(672i,680,185i), (104,697i,√428175), 153i]
- [(697i,104,153), (√-474993,680i,185), 672i]
- [(153i,185,680i), (104,697i,√-428175), 672i]
- [(104i,680,697i), (672,√-496625,153i), 185i]
- [(672,680i,153i), (104i,√428175,697i), 185i]
- [(185,104i,697i), (153,672i,√-496625), 680i]
- [(153,185i,672i), (104i,√-428175,697i), 680i]
All 24 cuboids are different, it's obvious and also it is not difficult to prove that they can not be produced by some other almost-perfect cuboid in real numbers, because for all of them max integer part from all rational lengthes is equal to 697. All cases are obtained combinatorially and do not produce additional integer parts of lengthes, so even if, let's imagine, some other almost-perfect cuboid has the body diagonal equal to 697, it has different from this cuboid lengthes.
Note, the combinatorial method for obtaining cuboids nevertheless generates 2 new under square root expressions: √428175 and √496625. Until the opposite is proven we assumed that these radicands can turn out to be a full squares, so potentially almost-perfect cuboid can produce so-called Perfect Complex cuboids.
We found that for cuboids in complex numbers their ability to produce new derivative cuboids does not depend on what is irrational side, an edge or a face diagonal, but whether cuboid contains face diagonals in complex numbers. That is why we had decided to change Randall's notation of almost-perfect cuboids in complex numbers and called them:
- Imaginary cuboid --- cuboid whose only edge(s) is(are) complex number(s), such as:
[(a, b, c), (d, e, f), g] = [(√-3344, 60, 63), (16, 25, 87), 65].
- Twilight cuboid --- cuboid whose only edge(s) and face diagonal(s) are complex numbers, such as:
[(a, b, c), (d, e, f), g] = [(60i, √3344, 65), (16i, 25, 87), 63].
- Midnight cuboid --- cuboid whose space diagonal is complex number, such as:
[(a, b, c), (d, e, f), g] = [(60i, 63i, √3344), (16i, 25i, 87i), 65i].
After series of researches we found a way to generate additional almost-perfect cuboids in complex numbers based on the next rules:
- The existence of any Face cuboid entails the existence of two different Imaginary cuboids or maybe Perfect Complex cuboid(s) if under root expression is a full square. For example:
Face cuboid [(104, 153, 672), (185, 680, √474993), 697] entails:
-
Imaginary cuboid [(104i, 185, 680), (153, 672, √496625), 697]
-
Imaginary cuboid [(153i, 185, 672), (104, √428175, 697), 680]
-
The existence of any Body, Edge, Face or Imaginary cuboid entails the existence of tree different Twilight cuboids. For example:
Face cuboid [(104, 153, 672), (185, 680, √474993), 697] entails:
-
Twilight cuboid [(153i, 104i, 697), (185i, 680, √474993), 672]
-
Twilight cuboid [(672i, 104i, 697), (680i, 185, √474993), 153]
-
Twilight cuboid [(672i, 153i, 697), (√-474993, 153, 672), 104]
-
The existence of any Body, Edge, Face, Imaginary or Twilight cuboid entails the existence of Midnight cuboid by multiply of all its sides on the imaginary unit.
-
The existence of Perfect cuboid in natural numbers entails the existence of 7 different Perfect cuboids in complex numbers which are constructed in this way: Perfect cuboid [(A, B, C), (D, E, F), G] entails:
-
Perfect Complex cuboid [(Bi, Ci, G), (Fi, E, D), A]
-
Perfect Complex cuboid [(Ai, Ci, G), (Ei, F, D), B]
-
Perfect Complex cuboid [(Bi, Ai, G), (Di, E, F), C]
-
Perfect Midnight cuboid [(Ai, Bi, Ci), (Di, Ei, Fi), Gi]
-
Perfect Midnight cuboid [(B, C, Gi), (F, Ei, Di), Ai]
-
Perfect Midnight cuboid [(A, C, Gi), (E, Fi, Di), Bi]
-
Perfect Midnight cuboid [(B, A, Gi), (D, Ei, Fi), Ci]
We decided do not search for Midnight cuboids because of they always have their "sunny brother", so we focused on cuboids with the body diagonal in real numbers.
The application is easily scalable, the width of each task can be arbitrary set. So the application was prepared for voluntary Distributed Computing. September 3, 2017 we have launched the sub-project Perfect Cuboid of the project http://www.rechenkraft.net/yoyo/. We chose 2 as a minimum quorum for workunit validation because some higher ranges were expected to be empty, without almost-perfect cuboids at all. So at least 2 clients must return not only identical found cuboids lists, but also their task check sums must be equal.
December 16, 2017, after 3.5 months we had successfully completed the 1st Batch and reached 2^50 = 1,125,899,906,842,624 threshold. Totally we had spent 84.38 CPU Years and had contributed 63,860,755,968 GFlops (63.86 quintillion floating-point operations).
As expected a perfect cuboid had not been found, but we had discovered:
- 194,652 Edge cuboids
- 350,778 Face cuboids
- 843,594 Imaginary cuboids
- 4,167,072 Twilight cuboids
All almost-perfect cuboids were verified and confirmed.
After a short discussion we had found a way to accelerate the program at the expense of refusal almost-perfect cuboids search. Our new version was more than 4 times faster then the previous one, so we launched it at the range 2^50-2^51 and completed in 24 days with the result:
- 6,665 Face cuboids
- 13,330 Imaginary cuboids
- 59,985 Twilight cuboids
As you can see, all Imaginary and Twilight cuboids was derivated from Face cuboids in 2:1 and 9:1 ratios respectively. And there we should declaim that 6,665 Face cuboids is not a full amount of Face cuboids with body diagonal from the range 2^50-2^51. As a result of program acceleration we refused most of almost-perfect cuboids checks but found possible to save findings of some specific views of Face cuboids without significant impact on application performance... and also for fun. It's always better to know the program finds at least something than receive empty result from each workunit.
We prolonged the search further and covered the range 2^51 (2,251,799,813,685,248) - 2^53 (9,007,199,254,740,992) in almost 9 months and discovered a total of 208'896 different cuboids:
- 17,408 Face cuboids
- 34,816 Imaginary cuboids
- 156,672 Twilight cuboids
and the capacious abstract became:
If a perfect cuboid exists, its body diagonal exceeds 2^53
October 16, 2018, the 3rd Batch was completed. Totally we had spent 227.66 CPU Years and had contributed 172,304,137,008 GFlops (172.3 quintillion floating-point operations).
After that we decided to suspend the project for application code review and since you have already here, we kindly offer you to take part in algorithm and code analysis to be sure our efforts were not in vain.
Honestly, very unlikely a perfect cuboid exists at all, so we rather expect the result of further project progress would become the only short abstract: "If a perfect cuboid exists, its body diagonal exceeds N". The only question is how high this N will.