Today I'd like to introduce a guest blogger, Mike Garrity, who is a graphics software developer at The MathWorks. Although his degree was in computational fluid dynamics, he has spent most of the last 25 years designing graphics hardware, rendering pipelines, and visualization software.
Contents
Introduction
Some numbers like 323 are palindromes. Other numbers like 124 are not. But look what happens when I add that number to a reversed copy of itself.
124 +421 ---- 545
Let’s try another.
150 +051 ---- 201
No, that didn’t work, but what if I keep going?
201 +102 ---- 303
There, it became a palindrome again. Let’s try another one.
291 +192 ---- 483
483 +384 ---- 867
867 +768 ---- 1635
1635 +5361 ----- 6996
That took several turns, but it did become a palindrome. I wonder what the pattern is. Let’s write a program to find out. I'll have it return the number of iterations and the final number.
The Code
Here's the code for producing palindrome numbers.
type pal196
function [n p] = pal196(in)
% Compute the number of iterations of
% reverse and add which are required
% to convert a number into a palindrome.
%
% Here are some fun ones.
% - 147996
% - 1000689
% - 9008299
%
n = 0;
maxn = 10000;
a = num_to_array(fix(in));
while ~all(a == fliplr(a))
a = add_rev(a);
n = n+1;
if (n > maxn)
error('Palindrome:TooManyDigits', ...
['Giving up on ', num2str(in), ' at ' num2str(n), ...
' iterations, value = ' array_to_str(a)]);
end
end
if (nargout == 2)
p = array_to_str(a);
end
end
function a = num_to_array(in)
% Convert a number into an array of digits
%
a = [];
while (in >= 1)
a = [mod(in, 10), a];
in = floor(in/10);
end
end
function out = add_rev(in)
% Add an array of digits to a reversed copy
% of itself.
%
out = in + fliplr(in);
for idx=size(out,2):-1:1
carry = 0;
if (out(idx) >= 10)
carry = 1;
out(idx) = mod(out(idx),10);
end
if (carry > 0)
if (idx > 1)
out(idx-1) = out(idx-1) + carry;
else
out = [carry out];
end
end
end
end
function str = array_to_str(a)
% Format an array of digits as a text string.
%
len = length(a);
lead = mod(length(a)-1, 3)+1;
str = sprintf('%d', a(1:lead));
if len > 3
str = [str, sprintf(',%d%d%d', a(lead+1:end))];
end
end
Explore the First Few Numbers
If I run that for all of the numbers from 1 to 100 and plot the number of iterations, I get this plot.
counts = zeros(1,100); for idx=1:100 counts(idx) = pal196(idx); end bar(counts);
That’s interesting. It looks like there’s a pattern, but I don’t quite get it. Let's keep going.
counts=zeros(1,200); firstfailed = []; try for idx=1:200 counts(idx) = pal196(idx); end catch ME firstfailed = idx; end if ~isempty(firstfailed) disp(['Failure occurred for n = ', num2str(idx)]) end if exist('ME', 'var') ME ME.message(1:45) [ME.message(46:72),'..........', ME.message(end-27:end)] end
Failure occurred for n = 196
ME =
MException
Properties:
identifier: 'Palindrome:TooManyDigits'
message: [1x5591 char]
cause: {}
stack: [6x1 struct]
ans =
Giving up on 196 at 10001 iterations, value =
ans =
6,643,864,828,137,164,591,..........,691,854,727,307,295,584,356
For the number 196, it quit at 10,000 iterations. At that point, the number had 4,159 digits. Wow, that's a big number!
Does It Ever Become a Palindrome?
So 196 doesn't work. Or does it? What if we ran it longer before giving up? Do you think that it would converge on a palindrome if we went further?
How far can you run it out?
Before we try to go too far, perhaps we should optimize the code a bit. What could you change to make it run faster and use less memory?
References
Here's a few interesting references to follow up on this topic.
Other Values?
What other numbers can you find that behave like 196? Also, some of the numbers that don't seem to converge seem to fall into a small number of common patterns. Can you identify any of these (by numeric experimentation preferably, not web look up)? Please add your observations here.
Get
the MATLAB code
Published with MATLAB® 7.8
Well, with all due respectr, improving on that code wasn’t too difficult. The way you convert numbers to arrays and then implement your own addition isn’t very efficient. I used a trick Loren mentioned a while back to use the string manipulation functions.
I also removed all error handling to clarify the actual algorithm. This code (with error checking) runs about twice as fast as the code above. As almost all time is spent in the reverseNum(), I believe the speed can be halved once again by calling it only once per iteration and cacheing the result.
function [n p] = pal196(in) n = 0; a = fix(in); while ~all(a == reverseNum(a)) a = a + reverseNum(a); n = n+1; end end function [ a ] = reverseNum( b ) string = num2str( b ); a = str2double(string(end:-1:1)); endIt seems that just below every multiple of 100 there is at least one non-converging. Numbers just above an even multiple of 100 are either palindrome themselves or become so fast.
Yes, your version is faster, but you need to worry about the fact that the numbers can get very large when you’re doing this. For example, try this starting value: 1000689. It should take 78 turns to converge on 796,589,884,324,966,945,646,549,669,423,488,985,697. If you store that in a double precision number, you’ll only get the first 16 digits. You can’t check whether it’s a palindrome without having all 39 digits.
It is possible to get your approach to work, but you’ll need to detect when you’re about to start losing precision and break your number into pieces. I found it easier to just start out with the number in pieces. Of course my version is using a double precision value to represent each digit. That’s obviously overkill. You really only need ~3.3 bits per digit rather than 64.
I think John D’Errico’s variable precision integer toolbox could come in handy for large numbers!
@Mike
Yeah, my actual code has a test in it to detect precision problems. For starting numbers less than 10000, precision was not a problem.
BTW, shouldn’t one be able to speed things up by alot if one used not one digit per piece but as many as can be fit in a 32-bit number. Then you just check if you have a carry on the most significant digit… just a thought.
Another way to work with these numbers is via the Symbolic Math Toolbox.
–Loren
Some numbers are clearly of this class. From 196, then 691, and 196+691=887. Then numbers, which can reach 887 are of this class, i.e. 196, 295, 394, 493, 592, 691, 788, 790 and 887. Further more, 689 and 986 this is because 689 + 986 = 1675 = 887 + 788. However, I do find another pair numerically, 879 and 978, which does not converge to (887, 788) pairs quickly.
Yes, 196, 295, 394, et al are called “kin numbers” and they converge on the same “thread”. The lowest number that leads to a thread is called the “seed”, so that would be 196 in this case. I think that the next seed after 879 is 1997.
I’ve been trying to come up with a good visual representation of the threads. It’s tough because they’re a very sparse set in a very big, linear space, but there are clearly some interesting patterns in how the threads are interwoven.
Hi! I know this is for MATLAB but I was very interested in doing this a different way. I’m using PHP to do the same thing. I have a web page up to play with. I used bcadd() and strrev() to find the palindrome and do the sums with infinite precision.
http://www.mrsankey.com/pal.php
I don’t feel like escaping the code right now. I’m tired. It wasn’t that tricky. Add and check, add and check, blah, blah…
We cannot assume that any of these numbers will not converge eventually, no matter how improbable. The number 10715 converges after 547 iterations. I haven’t found a number that doesn’t converge until you get above 1000 iterations. There seemed to be a ceiling around 760. Maybe we should offer a prize?