%% Hardy's taxi number 1729 to see Ramanujan in hospital 
%  10.8.2016 Heikki Apiola, HardyRamanujan.m
%
%% Hardy and Ramanujan
%
% Here's a list of references:
%
% * <http://www.imdb.com/title/tt0787524/> The man who new infinity, Movie
% * <http://matematiikkalehtisolmu.fi/2008/1/ramanujan.pdf> Article about
% Ramanujan in "Solmu", in Finnish
% * G.H. Hardy: A mathematician's apology, Cambridge U.P. 1940, Finnish
% translation, Terra Cognita 1997.
% 
% The problem treated here is the classical one mentioned in all the above
% references
% " ...the famous taxicab incident in which Hardy says the taxi license plate 
% has a very dull number, but Ramanujan says "it is a very interesting 
% number"- the lowest number expressed as the sum of two cubes in two 
% different ways ..."
%
% Our starting point is to compute the smallest number that can be so
% expressed without knowing it in advance. So it is kind of converse to the
% incredible "sight" of Ramanujan. On the other hand we will be able to
% solve as many as the computing power will allow, to get more and more
% such numbers. Here we will just show the next such number: 4104.
%
% Also, in the Solmu-reference above was a continuation: Hardy would have
% asked if Ramanujan knew the answer if cubes were replaced by power 4.
% Ramanujan would have thought a while and then said: "It must be a very
% large number". I used near to one hour CPU-time (lunch + coffee) on my 
% laptop to show that it must be more than one million.
%
% The first attempt works with p=3 but is too inefficient for p=4.
% Handling of the case p=4 forced us to design a more efficient way
% of handling the problem. Vectorization (elination of for-loops) 
% and re-thinking the algorithm led us to write the function:
%%
%    *H_RTaxi_Nr*, which is called at the and of this script.
% 
%% Form an M x M grid of pairs of numbers (in the meshgrid-way)
% 
%% 
%        .                                  .
%
% *If you need help in |meshgrid|, here are two links:*
%
%%
% <meshscript.html>
%
% <https://math.aalto.fi/~apiola/matlab/opas/lyhyt/grafiikka.html#surf>
close all
clear
M=16;
i=1:M;
j=1:M;
[I,J]=meshgrid(i,j);
%
% Compute a matrix of sums of cubes for each mesh-point:
%
qubesums=I.^3 + J.^3;
%% Build a matrix N which counts the nr. of times each element of qubesums appears in the matrix. 
%  
for k=1:M
    for l=1:M
        N(k,l)=sum(sum(qubesums(k,l)==qubesums));
    end
end    
%% If there are 2 representations, by symmetry, there are 4 duplicates. 
%
[k,l]=find(N==4);
values=k.^3 + l.^3
[k l values]
Taxinr=min(values)
%%
% This shows that 9^3  +  10^3   =     1729
%  and
%                  12^3 + 1^3      =   1729
%  Also:
%           9^3+15^3  = 2^3 + 16 ^3   = 4104
%
%% Some graphics for illustration:
surfc(i,j,qubesums), colorbar
%%
figure
hold on
contour(i,j,qubesums,20)
[c,h] = contour(i,j,qubesums,[1729 4104]);
clabel(c,h);
%% 
plot(9,10,'or','MarkerSize',10)
plot(1,12,'or','MarkerSize',10)
plot(9,15,'ob','MarkerSize',7)
plot(2,16,'ob','MarkerSize',7)
grid on
grid minor
xlim([0 M])
title('Ramanujan-Hardy, sum of two cubes')

%% Forth powers, rethinking the algorithm.
% 
% Better write it as a function: *HRTaxiNr*
%
%%
help HRTaxiNr
%% Here's the whole code for HRTaxiNr
%
type HRTaxiNr
%%
% Looking at the code, in finding the number of occurrencies, I wrote a 
% worksheet of various approaches:
%%
% <countoccurrencies.html>
%%
[k,l,mindu]=HRTaxiNr(2)   % Sum of squares
%%
[k,l,mindu]=HRTaxiNr(3)   % Sum of qubes
%%
% Here's the taxi number
%%
[k,l,mindu]=HRTaxiNr(4)   % Sum of fourth powers
%
% So we have the answer (635318657) to Hardy's question, where Ramanujan could say:
% 'It must be a very large number'




