|S&SMC#14: My Original Solutions and Comments [LONG]|
Message #27 Posted by Valentin Albillo on 13 Mar 2006, 6:13 a.m.,
in response to message #1 by Valentin Albillo
Thanks to all of you who were interested in this S&SMC#14, and most specially
to the outstanding ones (Marcus von Cube, Eamonn, Gerson W. Barbosa, Mr. Hohmann) who actually cared to take the time to develop and post their very clever solutions, I really hope you enjoyed it and consider the afforded time well spent. These are my original solutions, with comments.
Conjecture 1: Well-done
This is a particular case of a plausible conjecture issued by Euler in 1769,
which resisted every effort to prove or disprove it for almost two centuries until it was finally disproved in 1966 when
the first counterexample was found by Lander and Parkin, using a then quite
powerful CDC 6600 computer.
Now, we're asked to duplicate that feat, namely finding a non-trivial solution
=E5 in non-zero integers, using our favorite handhelds.
A straightforward, dumb brute-force search would require excessive running
times (flattening several sets of batteries in the process), but a few
refinements here and there will cut down the task by several orders of
magnitude. However, these refinements need a significant amount of RAM to
implement, which means only advanced models like the HP-41CX, HP42S,
HP-71B, and HP48/49 can be profitably used.
Here's my original 7-line, 251-byte solution for the HP-71B
1 DESTROY ALL @ OPTION BASE 0 @ K=150 @ DIM P(K)
2 FOR I=0 TO K @ P(I)=I*I*I*I*I @ NEXT I @ R=.2 @ L=2^R @ M=3^R @ N=4^R
3 FOR E=K TO 1 STEP -1 @ T=P(E) @ DISP E
4 FOR D=E-1 TO E/N STEP -1 @ F=T-P(D) @ U=F^R
5 FOR C=INT(U) TO U/M STEP -1 @ G=F-P(C) @ V=G^R @ FOR B=INT(V) TO V/L STEP -1
6 IF NOT FP((G-P(B))^R) THEN A=(G-P(B))^R @ DISP A;B;C;D;E @ END
7 NEXT B @ NEXT C @ NEXT D @ NEXT E
upon running, it eventually outputs:
which is the sought-for counterexample, as indeed
27 84 110 133 144
275 + 845 + 1105 + 1335 = 61917364224 = 1445
The running time is a short 90 min. in an actual, physical HP-71B, and only
13 seconds under Emu71 in my IBM laptop. To achieve such short times
the following techniques are used to speed the search, they key of which is to
iteratively reduce the problem to a similar yet easier one:
- first of all, a table of 5th-powers for integer arguments from 0 to 150 is
precomputed, so that during the search the extremely frequent
need of raising values to the 5th power is reduced to indexing
an array element, which is significantly faster.
- some other needed constants are precomputed as well, which saves
additional time when they're used inside the search's nested
loops, namely 1/5, 21/5, 31/5,
- the outer loop traverses all possible values for E in reverse
order, from the largest possible value to the smallest.
- for earch value of E, we must try to decompose E5 as the sum of
four 5th-powers. Assuming D is the
we must search for
a decomposition of E5-D5 in three
5th-powers where, as D is the largest element by definition, we only need to
search for values of D in the limited range:
(E5/4)1/5 <= D
< (E5)1/5 i.e. E/41/5 <= D < E
which we traverse in reverse order as well. For each value of D,
we must now try to decompose E5-D5
as the sum of three 5th-powers,
so we've reduced our problem to a similar but simpler one,
and the same argument is repeated till a single summand is
left, which we simply check to ascertain whether it is a
5th-power or not (i.e.: the fractional part of its fifth root equals 0).
When this condition is met, we've found a counterexample so
the program outputs it and stops.
Conjecture 2: Medium
We're asked for the smallest integer number from 6 to 2000 that can't be
decomposed as the sum of a prime and a non-zero power.
A brute-force search would again require excessive running time, but
we can trade RAM for speed once more. This is my original 9-line,
259-byte solution for the HP-71:
1 DESTROY ALL @ OPTION BASE 1 @ M=2000 @ D=LOG(M) @ N=INT(SQR(M))
2 INTEGER P(M),Q(310) @ FOR B=2 TO N
3 FOR E=2 TO INT(D/LOG(B)) @ P(B^E)=1 @ NEXT E @ NEXT B @ P(1)=1 @ Q(1)=2 @ I=2
4 FOR N=3 TO M-1 STEP 2 @ FOR D=3 TO INT(SQR(N)) STEP 2 @ IF MOD(N,D)=0 THEN 6
5 NEXT D @ Q(I)=N @ I=I+1
6 NEXT N @ T=I-1 @ FOR N=6 TO M
7 FOR D=1 TO T @ E=N-Q(D) @ IF E<=0 THEN DISP N @ END ELSE IF P(E) THEN 9
8 NEXT D @ DISP N @ END
9 NEXT N
upon running, it outputs:
which is the smallest counterexample. The running
time is a fast 30 min.
in a physical HP-71B, and only 6 seconds under Emu71 in my IBM laptop,
thanks to these simple yet efficient techniques:
- first of all, this challenge's actually far easier than it seems at first.
Actually, for the search to proceed fast and smoothly, we need but two things: (1) to have
all relevant prime numbers instantly at hand, and (2) to be able to test
extremely quickly whether a given integer is
a perfect power or not.
This is achieved as follows:
- all primes up to 2000 are precomputed and stored in an array so
that they can be retrieved very fast during the search
- in order to quickly check whether a given integer is a perfect power or
not, we stablish a 2000-element 'flag' array where each element
has the value 1 or 0 depending on whether the index is a perfect
power or not, i.e., P(7)=0, P(8)=1 ( as 8=23), etc. This could be achieved using much less RAM by packing eight such 'flags' per byte in a string, for instance, but here speed is our top priority so we simply define an integer array and let each individual element act as a 'flag'.
- with primes and very fast power-detection available, the main
search just traverses all possible values from 6 to 2000, in
ascending order. For each, every prime is subtracted in turn
and the result is checked to see if it is a perfect power,
skipping to the next value if it is. When no prime subtracted
will result in a perfect power remaining, the corresponding
value is output as a valid counterexample.
- by the way, the prime array is dimensioned to have 310 elements
because there are some 300+ prime numbers in the range from 2
to 2000 (actually, 303). If you've got INTEGRATE capabilities
in your HP model, you can very quickly compute an approximate
value for the required number of elements as follows:
approx. #primes up to M = INT(INTEGRAL(2,M,1,1/LN(IVAR))
so that the line:
2 INTEGER P(M),Q(310) @ ...
2 INTEGER P(M),Q(INTEGRAL(2,M,1,1/LN(IVAR))) @ ...
and this allows you to expand the search to numbers greater than
2000 by simply changing the value assigned to M at line 1,
without having to consult tables for the proper size of
array Q. As the integration's uncertainty is specified as 1, (i.e: FIX 0
or SCI 0 for models that use the display setting instead), the integral
is computed very quickly.
Now, you'll agree with me that a line which dimensions an
array with a size defined by a non-elementary integral
isn't that frequent a sight.
Conjecture 3: Rare
This is a well-known conjecture as yet unproved or disproved,
though everybody and his uncle feels that 196 is a true counterexample, as
it has failed to produce a palindrome after hundreds and hundreds of millions
of cycles. For instance, after 670,000,000 cycles you're dealing with numbers
280,000,000-digit long(!!), yet no palindrome in sight ...
In order to implement this challenge, multiprecision addition
is mandatory, which can be done with arrays or else
with strings, which is the technique I've used in my original
8-line, 364-byte solution for the HP-71B:
1 DIM M$,N$ @ INPUT "#Cmax=";M @ FOR I=0 TO 200 @ N$=STR$(I) @ K=0
2 CALL MADD(N$,M$,K) @ IF M$=REV$(M$) THEN 4 ELSE IF K<M THEN N$=M$ @ GOTO 2
3 DISP I;"fails (";M$;" after";M;"cycles)"
4 NEXT I @ DISP "OK" @ END
The algorithm is completely straightforward and no special techniques are
needed, though both for speed and to avoid obscuring the inner works with
trivial utility routines, I've made use of several string-handling
keywords (REV$, RPT$, LTRIM$) which are available in a number
of very common LEX files and ROMs (STRNGLEX, REVLEX, RPTLEX, JPC ROM, etc).
The program works as follows:
5 SUB MADD(A$,C$,K) @ DIM B$ @ L=LEN(A$) @ E=10^11 @ C$="" @ C=0
6 B$=REV$(A$) @ FOR I=L TO 1 STEP -11 @ C=VAL(A$[I-10,I])+VAL(B$[I-10,I])+C
7 D$=STR$(MOD(C,E)) @ C$=RPT$("0",11-LEN(D$))&D$&C$ @ C=C DIV E @ NEXT I
8 C$=LTRIM$(STR$(C)&C$,"0") @ K=K+1 @ END SUB
- the multiprecision results will be stored in strings (which in the
case of the HP-71B would allow us to handle up to 65000-digit numbers!),
initially dimensioned to hold up to 500-digit numbers, more than enough
for up to 1000 cycles.
- to perform the multiprecision addition of a number and its mirror
image, a call to the MADD subprogram is made. This subprogram takes
the number as one of its string arguments pased by value, and returns
the result of the addition as another string argument passed by
reference. For instance:
>CALL MADD("78",M$,0) @ M$because 78+87 = 165 and 8263485213753244473212+2123744423573125843628
= 10387229637326370316840. Actually, the subprogram's code could
be inserted directly in the main program proper, so we could get
rid of the subprogram altogether and get a slightly faster, shorter
program (7 lines instead of 8) but 'outsourcing' particular tasks
to subprograms encourages modular programming and makes for clearer,
cleaner code and provides additional functionality as well.
>>CALL MADD("8263485213753244473212",M$,0) @ M$
After calling the subprogram, the main program just checks if the new
result is palindromic,
or if we've exhausted the number of cycles, and iterates as needed.
Upon running, it produces the following, for assorted maximum number of cycles
(10, 20, 40, 100, 200, 1000):
So it's clear that 196 is a very firm candidate
for a counterexample. Of course, it's not the only one, other candidates
include (up to N=2000): 295, 394, 493, 592, 689, 691, 788, 790, 879, 887, 978, 986, 1495,
1497, 1585, 1587, 1675, 1677, 1765, 1767, 1855, 1857, 1945, 1947, and 1997.
89 fails (8872688 after 10 cycles)
98 fails (8872688 after 10 cycles)
167 fails (17050517 after 10 cycles)
177 fails (17794887 after 10 cycles)
187 fails (17735476 after 10 cycles)
196 fails (18211171 after 10 cycles)
89 fails (93445163438 after 20 cycles)
98 fails (93445163438 after 20 cycles)
187 fails (176881317877 after 20 cycles)
196 fails (70446464506 after 20 cycles)
196 fails (13305261530450734933 after 40 cycles)
196 fails (44757771534490515617290699271561508443627774644 after 100 cycles)
196 fails (910449546741765655298269802255629632301207255281
2103235826563197972803556567037646054008 after 200 cycles)
196 fails (35346644392413689785837714402912114362859098083
41012030342772858788740429334664452 after 1000 cycles)
That's all. Thanks again and
Best regards from V.