There is a scene in one of the Tintin books, the Black Island I think, where Tintin scolds Snowy (or Milou if you prefer) for letting his clever ideas run away with them. Wise words.
There is a wonderful book called "Algorithms to live by" by Brian Christian and Tom Griffiths in which, seeking to explain a Hamiltonian circuit, the authors use the observation that Abraham Lincoln was a circuit judge. Being a Gresham Professor, and realising I ought to use examples from my predecessors, I noted that Sir Christopher Wren provided an alternative exemplar. Within the City of London there remain thirteen of Wren's churches untouched by the remodelling efforts of developers or the Luftwaffe. Excellent! We can use them as an example.
Wikipedia lists their coordinates so, I thought, I'll write a script that takes their locations and computes the distance between each pair of churches. There are thirteen churches so that will give me a 13 by 13 distance matrix. And then I'll use a standard TSP algorithm to solve for the optimum route around the churches.
Then I paused; how to find the distance between churches? Well, if I have the coordinates, then I could model the earth as an oblate spheroid and measure the curvilinear distance. Too complicated, and anyway my audience, here I am picturing some svelte intelligentsia who have nothing better to do in a pandemic than yomp around 13 churches, my audience, will want to know the walking distance. No problem - I will use the Google maps walking distance.
Well obviously I am not going to enter, by hand, 56 pairs of churches so I need to access the Google maps API. Ah ha! Google maps has a distance matrix API. Lovely. Hang on, in order to access it I need to have a Google account and I need an API key. How do I get an API key? Well ignoring the slight twinge in my philosophical position (did I in previous lectures point out that Google was not entirely saintly?), thirty-minutes or so later, I have a Google API key. Right then how to interface? Well, the online notes give some examples in Python and C. I'm using Matlab but, not problem, I can rewrite in Matlab and then Matlab can parse JSON objects into a single object. Oops, the API doesn't seem to allow lists of coordinates (It says it does but I get an error code). No problem, I will write a loop. Careful now as my Google account gets charged each time I call these loops. And here we have it:
clear
% Lat and Lon of churches in Decimal format
C = [51.511669, -0.099272,
51.511325, -0.086892,
51.512394, -0.0863,
51.511108, -0.093761,
51.514697, -0.088867,
51.510764, -0.082983,
51.514042, -0.101942,
51.511672, -0.088347,
51.512778, -0.093333,
51.513194, -0.085467,
51.513611, -0.098056,
51.51325, -0.084583,
51.512628, -0.089919];
[N,M] = size(C);
% Adjacency matrix
A = ones(N,N);
A = 0*A +triu(A,1) + tril(A,-1);
% Fully conencted graph
G = graph(A);
% Plot the graph
gp =plot(G);
gp.XData = C(:,2);
gp.YData = C(:,1);
Nodelabels = {'St Benet Paul''s wharf',
'St Clement Eastcheap',
'St Edmund, King and Martyr',
'St James Garlickhythe',
'St Margaret Lothbury',
'St Margaret Pattens',
'St Martin, Ludgate',
'St Mary Abchurch',
'St Mary Aldermary',
'St Michael''s Cornhill',
'St Paul''s Cathedral',
'St Peter upon Cornhill',
'St Stephen Walbrook'};
gp.NodeLabel = Nodelabels;
% Google API key
API_key = 'Your value goes in here';
% Google distance API is called as
%https://maps.googleapis.com/maps/api/distancematrix/outputFormat?parameters
%
% Quick test example - the walking distance between the first two churches in miles
%
% http_str = 'https://maps.googleapis.com/maps/api/distancematrix/json?units=imperial&mode=walking&origins=51.511669,-0.099272&destinations=51.511325%2C-0.086892&key=AIzaSyBr01KnUtOscbNL6vvug-VkJAkRamm6j8k';
%
% data = webread( http_str );
%
% The returned structure will have the distance in
% data.rows.elements.distance
% ans =
% text: '0.7 mi'
% value: 1132
% Now to construct our Distance matrix.
if (~exist('D.mat')),
preamble_str = 'https://maps.googleapis.com/maps/api/distancematrix/json?units=imperial&mode=walking';
term_str = ['key=',API_key];
% Construct a well-formed string of the church-coordinates
coord_str = []
for n = 1:N
origin_str = [num2str(C(n,1)),',',num2str(C(n,2))];
for m = 1:N
dest_str = [num2str(C(m,1)),'%2C',num2str(C(m,2))];
url_str = [preamble_str,'&','origins=',origin_str,'&','destinations=',dest_str,'&',term_str];
data = webread( url_str );
D(m,n) = data.rows.elements.distance.value;
end
end
% Save D - as each of those requests cost a bit!
save D
else
load D
endSimples! Well there is one step left which is to symmeterise D
Dt = round(0.5*(D+D.'));
Why do we need to symmetrise D? Well I'll leave that as an exercise for the reader but it makes the problem a bit more manageable and the asymmetry is mostly zero (two churches have routes between them that differ by 57m depending on whether you walk from A-to-B or B-toA). And we are there.
Now to find TSP. Let's use Google OR tools. Super that seems easy. The webpage explains how to do that in Python, C++, Java and C#. But I'm using Matlab. Option 1 is to work out how to call OR-tools from Matlab. This is all getting a bit out of hand. Maybe I should just save the distance matrix from Matlab and then run a Python program to read in the distance matrix and then return the Hamiltonian cycle? Yes, Option2, that is simpler. Oh hang-on, I just changed computer, so I need to install OR-tools, which means installing homebrew and then OR-tools and then the command line tools and then Xcode and several 100Gb later we are ready.
Well that was great. Now we have the shortest route, time is running a bit short now of course, as all of that takes a a day or two, so let's quickly plot out some of the graphs (using the lovely new Gresham colours but being careful because the Macintosh will colour calibrate the colours if you are not careful).
And that's how we got the example for the lecture.
One teeny weeny problem though. Before my lecture I wandered down from Barnard's Inn to my dentist on the east side of the City of London. Oh this will be fun, I'll pass some of these churches I thought. Hang on! They are all in the wrong place!
Well dear reader, in my haste to plot them out, I interchanged latitude and longitude on all my slides. Of course, the example still works and the shortest path is still the same but, if you want to read my slides as a map, you need to lie on your side. Or, if you really do fancy a Hamiltonian circuit around Wren's thirteen unmarked churches in the City of London here it is, this time the right way up! Enjoy!
And if you want to see the original lecture then it is available at Gresham College. The slides and the notes have Wren's churches on their side. I'm minded to leave them like that and see if anyone notices.