As part of my Numberphile series of posts:
here is another one about a recent video called The Square-Sum Problem
The question is: if you have the integers 1 through n, can you arrange that list in such a way that every two adjacent ones sum to a square number. As seen in the video and the extra footage.
We can easily check this in the Wolfram Language:
Let's see which number can pair up to a pair:
SquareEdges[n_Integer?Positive]:=Reap[Do[If[IntegerQ[Sqrt[i+j]],Sow[{i,j}]],{i,n-1},{j,i+1,n}]][[2,1]]
Now let's try for 15, as in the main video:
n = 15;
poss = SquareEdges[n];
gr = Graph[TwoWayRule @@@ poss, VertexLabels -> Automatic];
path = FindHamiltonianPath[gr, PerformanceGoal :> "Speed"]
HighlightGraph[gr, BlockMap[Rule @@ # &, path, 2, 1]]
giving:
{9, 7, 2, 14, 11, 5, 4, 12, 13, 3, 6, 10, 15, 1, 8}
In the extra footage, it is revealed that they found the solution for up to n=299. Can we do better? Yes we can! Changing n to 300 in the above code and rerunning gives us the solution in 0.28 sec on my laptop.
{289,35,65,259,30,294,67,257,32,292,69,100,44,125,71,154,135,189,211,113,248,8,281,119,205,195,166,158,283,6,250,191,133,156,285,4,252,277,12,244,117,207,193,168,273,16,240,160,164,236,20,269,131,94,230,59,197,92,232,57,199,90,234,22,267,217,224,137,152,73,123,46,150,75,121,48,148,77,179,110,214,270,19,237,163,161,239,17,272,128,41,103,297,27,262,62,227,97,99,190,210,114,175,50,146,79,177,112,212,188,253,3,286,155,134,266,23,233,91,198,58,231,93,196,60,229,95,130,159,165,276,13,243,118,206,194,167,274,15,241,288,1,255,186,138,223,218,143,181,108,88,201,55,170,86,203,53,172,84,37,107,182,142,299,25,264,220,221,140,184,216,225,64,36,85,171,54,202,87,169,56,200,89,235,21,268,132,157,284,5,251,278,11,245,116,208,192,249,7,282,247,9,280,204,120,136,153,72,124,45,151,74,122,47,149,76,180,109,215,185,139,222,219,265,24,300,141,183,106,38,83,173,52,144,81,40,104,296,28,261,63,226,98,127,42,102,298,26,263,61,228,96,129,271,18,238,162,279,10,246,115,209,275,14,242,287,2,254,187,213,111,178,78,147,49,176,80,145,51,174,82,39,105,295,29,260,101,43,126,70,291,33,256,68,293,31,258,66,34,290}
and a completely mess of a graph:
Can we go beyond? Let's optimize a code a bit, and let's find the solution for larger n:
Let's store our intermediate results in the association db:
SetDirectory[NotebookDirectory[]];
$HistoryLength=1;
db=If[FileExistsQ["squaresumdb.mx"],
Import["squaresumdb.mx"]
,
<||>
];
And now our main code:
ClearAll[SquareEdges,SquareEdges2,CheckSol,TryToFind]
SquareEdges[n_Integer?Positive]:=Reap[Do[If[IntegerQ[Sqrt[i+j]],Sow[{i,j}]],{i,n-1},{j,i+1,n}]][[2,1]]
SquareEdges2[n_Integer]:=Module[{tmp},
tmp=Table[
{i,#}&/@(Range[Ceiling[Sqrt[2 i]],Floor[Sqrt[i+n]]]^2-i)
,
{i,1,n-1}
];
tmp=Join@@tmp;
Select[tmp,Less@@#&]
]
CheckSol[l_List]:=Sort[l]===Range[Length[l]]\[And](And@@BlockMap[IntegerQ@*Sqrt@*Total,l,2,1])
TryToFind[n_Integer?Positive]:=Module[{edges,out},
If[!KeyExistsQ[db,n],
edges=SquareEdges2[n];
If[Union[Flatten[edges]]===Range[n],
edges=TwoWayRule@@@edges;
edges=RandomSample[edges];
Do[
out=TimeConstrained[FindHamiltonianPath[Graph[edges],PerformanceGoal:>"Speed"],5+i,$Failed];
If[out=!=$Failed,
If[Length[out]==0,
Print[Style["No solution for ",Red],n];
,
status=Row[{"Found solution for ",n,":",i}];
];
AssociateTo[db,n->out];
Break[]
];
Print["Failed ",n,":",i];
edges=RandomSample[edges];
,
{i,5}
]
,
Print["Edges are not connected for ",n];
AssociateTo[db,n->{}]
]
]
]
Let's scan the first 1000:
Dynamic[status]
status = "";
Do[TryToFind[k], {k, 3, 1000}]
Export["squaresumdb.mx", db];
Note that if finding the Hamiltonian path takes too long I mix the edges and try again, sometimes, seemingly random, it then finds the solution quickly.
I can tell you now that all of them have a solution. In fact I went up to larger numbers and found that all up to 2667 have a solution, and possibly beyond. I attached the notebook and the solutions in form of a mx file.
Attachments:
