# Search for solutions to quartic diophantine equation

I wrote this code to find the variables that will satisfy the diophantine equation (i.e. a^4 + b^4 + c^4 = d^4), for all variables belong to: (0 , 430000), and it takes forever to find the answer (I left the code running for two hours and got nothing). I know that this is caused by the algorithm or the approach I took, but I really couldn't think of any better solution...

here's the code:

#include <iostream>
using namespace std;
void main()
{
//Here I use double because max value for int variables is too small for what I need.
double n = 0, LHS, RHS, valA, valB, valC, valD;
double *arrA, *arrB, *arrC, *arrD;
arrA = new double[430000]; //array to store all the powers of 4 for var A
arrB = new double[430000]; //array to store all the powers of 4 for var B
arrC = new double[430000]; //array to store all the powers of 4 for var C
arrD = new double[430000]; //array to store all the powers of 4 for var D
for (int count = 0; count < 430000; count++)
{
arrA[count] = n*n*n*n;
arrB[count] = n*n*n*n;
arrC[count] = n*n*n*n;
arrD[count] = n*n*n*n;
n++;
}
//this loop is a counter, it increments "b" by 1 when "a" reaches max value, and so on for each variable...
/*
N.B:
The variables that satisfy the equation should be the following:
a = 95800
b = 217519
c = 414560
d = 422481
Plug it in as initializers for the variables to check that the program is working fine.
*/
int a = 0, b = 0, c = 0, d = 0;
while (d < 430000)
{
if (a == 430000)
{
a = 0;
b++;
}
if (b == 430000)
{
a = 0;
b = 0;
c++;
}
if (c == 430000)
{
a = 0;
b = 0;
c = 0;
d++;
}
if (d == 430000)
{
break;
}
LHS = arrA[a] + arrB[b] + arrC[c];
RHS = arrD[d];
if (LHS == RHS && (arrA[a] != arrB[b] != arrC[c] != arrD[d]))
{
valA = a; valB = b; valC = c; valD = d;
break;
}
a++;
}

if (d == 430000){
cout << "Couldn't find an answer! \n" }
else
{
cout << "value of A is:  " << valA << endl;
cout << "value of B is:  " << valB << endl;
cout << "value of C is:  " << valC << endl;
cout << "value of D is:  " << valD << endl;
}
}


If anyone could tell me how to speed up the code or could think of a better algorithm. There's nothing wrong (IMO) in the code's syntax or logic, but it is very very slow...

N.B: here's a link to an article regarding the equation for further reading: https://www.jstor.org/stable/2008781?seq=1#page_scan_tab_contents

I'm using a 2013 macbook pro "running windows 10", (intel i5 2.4GHz, 8GB RAM), MS Visual Studio 2017...

• – JS1
Jan 2, 2018 at 23:17
• @JS1 Thank you so much for referencing this! :) Jan 2, 2018 at 23:37
• Welcome to Code Review! I changed the title so that it describes what the code does per site goals: "State what your code does in your title, not your main concerns about it.". Feel free to give it a different title if there is something more appropriate. Jan 3, 2018 at 0:08
• You should take a look at the code profiling tools provided by your IDE. They can run your application and tell you which parts are taking the longest to process. That can help you figure out where you should start when optimizing. Jan 3, 2018 at 4:41
• Too small for an answer, but you don't need to check every combination (i.e. a=1,b=0,c=0 is equivalent to a=0,b=1,c=0 and to a=0,b=0,c=1) of values. Change your looping mechanism or cache the ones you've already done. Jan 4, 2018 at 23:20

This looks like a bug:

if (LHS == RHS && (arrA[a] != arrB[b] != arrC[c] != arrD[d]))


The operator != returns a boolean result. Its Associativity is also left to write. so the above is equivalent to:

if (LHS == RHS && (((arrA[a] != arrB[b]) != arrC[c]) != arrD[d]))


So:

arrA[a]         != arrB[b]     =>  true or false
(true or false) != arrC[c]     =>  true or false  (so usually true)
true            != arrD[d]     =>  true or false  (so usually true)

LHS == RHS && true             =>  Looking for LHS == RHS


This is a bad way of writing a nested loop:

while (d < 430000)
{
if (a == 430000)
{
a = 0;
b++;
}
if (b == 430000)
{
a = 0;
b = 0;
c++;
}
if (c == 430000)
{
a = 0;
b = 0;
c = 0;
d++;
}
if (d == 430000)
{
break;
}
LHS = arrA[a] + arrB[b] + arrC[c];
RHS = arrD[d];
// stuff
}


Just write:

for(int a = 0; a < 430000; ++a) {
for(int b = 0; b < 430000; ++b) {
for(int c = 0; c < 430000; ++c) {
for(int d = 0; d < 430000; ++d) {
LHS = arrA[a] + arrB[b] + arrC[c];
RHS = arrD[d];
// stuff
}
}
}
}


There is no point looping over all the values.

if a^4 + b^4 + c^4 > 430000^4 then there is no point in running through the d loop as you will never find the value.

There is no point in checking all values of d. If you must use a search arrD then do a binary search of the array. This will vastly reduce the number of items checked in the innermost loop. Alternatively just find the cube root of the result using a maths function.

 double root = std::pow(LHS, 1.0/4.0);
double base = std::floor(root);
if (floor * floor * floor * floor == LHS) {
break;
}


I am not convinced that calculating and saving all the powered values is going to save time. Also you don't need to save the same result four times. Rather than have arrA and arrB and arrC and arrD have a single array arr. It is the index into the array that is important.

Also don't manually allocate this space with new use a vector<double>.

The question is a memory lookup faster than or slower than doing 4 multiplications (or callling std::pow(x, 4)). That is something you will need to time and find out.

• Why should I user vector instead? I tried using vectors and found out that the loops take much longer time...maybe it's because of the way I wrote the loop? anyways, I'll take in account what you've mentioned and re-write the code. many thanks! Jan 3, 2018 at 20:12
• Vector is no slower than bare C-Arrays. If you find any timing difference there is something wrong with your timing (a lot of work has gone into making sure that there is no speed difference). stackoverflow.com/a/3664349/14065 Jan 3, 2018 at 21:14
• The reason you should use vector is that it is RAII. There will never be any memory leaks (so your program executes correctly). Jan 3, 2018 at 21:15

First, a bit of quick analysis. The solutions with $a=b=0$ are trivial: $c=d$. There are no solutions with $a=0$ other than those trivial solutions, by Wiles' theorem. Therefore we may as well ignore the cases where the variables are zero.

You don't need to loop over four variables ($N^4$ tuples) to solve $a^4 + b^4 + c^4 = d^4$. As Loki Astari has already pointed out, given $a,b,c$ you just need to test $a^4 + b^4 + c^4$ to see whether it is a fourth power: so $N^3$ tuples times a lookup (binary chop is a reasonable suggestion). Combined with Dannnno's suggestion in a comment on the question that without loss of generality $a \le b \le c$ you can reduce that further to approximately $\frac{1}{6}N^3$ tuples.

But actually you can take the loop reduction a step further: rewrite the equation as $a^4 + b^4 = d^4 - c^4$ (with $c < d$). Then (at an abstract level) you can write two parallel loops, one over $(a, b)$ and the other over $(c, d)$, processing a total of about $N^2$ tuples.

The implementation detail is that $N$ is too large to store the tuples from one of the loops. A standard approach to get around that is to generate the two sequences in order using priority queues. To keep it simple, I would initialise the ab queue with $(a, a)$ for each $1 \le a \le N$, and the cd queue with $(d-1, d)$ for each $2 \le d \le N$. Then the main loop is (Pythonesque pseudocode):

while not ab.isEmpty() and not cd.isEmpty():
(a, b), (c, d) = ab.peek(), cd.peek()
if a^4 + b^4 + c^4 == d^4:
print "Solution", a, b, c, d

if a^4 + b^4 <= d^4 - c^4:
ab.pop()
if b < N:
ab.push((a, b+1), a^4 + (b+1)^4)
else:
cd.pop()
if c > 1:
cd.push((c-1, d), d^4 - (c-1)^4)