[PYTHON] I ran all over Japan at the shortest distance at 10km / h.

How many days can it be achieved?

I was studying optimization problems and came up with a problem. __ How long does it take to travel all over Japan in the shortest distance? __ In order to answer the question, I calculated a little with python.

Since the purpose is to hands-on the optimization problem with python, please forgive the sweetness of the problem setting. Also, as I wrote in this article, I couldn't solve the problem of creating a circular route other than the favorite route. It is an approximate solution, not an optimal solution. If you are familiar with the solutions around here, please comment! !!

Problem setting

This time I will use Takashikun, which everyone loves (I don't know if it's true) that often appears in math problems.

Takashi-kun specs

・ __ Takashi-kun is strong __ so you can cross the sea. ・ __ Takashi-kun is strong __ so you can run at a constant speed regardless of whether there are mountains or valleys. ・ __ Takashi-kun is strong __ so you can run for 24 hours.

Configuration

・ __ Takashi-kun has a lot of intellectual curiosity, so I thought about going around Japan at a speed of 10km / h. ・ "All over Japan" is defined as "the location of all prefectures" this time. ・ Since __ Takashi-kun is impatient, all movements will be done in a straight line.

I tried it

Acquisition of prefectural office location data

The longitude and latitude of the prefectural capital are quoted from this article. After downloading Excel, it will be molded only to the location name, longitude, and latitude data.

import pandas as pd
import numpy as np
df = pd.read_excel("city_list.xlsx")
display(df)

The data looks like this. lati is latitude and long is longitude. image.png

Distance calculation between each location

Calculate the distance between 47 locations. Originally, I think it would be good to narrow down the combinations, thinking, for example, "Okinawa-Hokkaido will not move." This time, ~~ the selection work was troublesome, and ~~ __ I had a faint expectation that something miracle might happen __, so I calculated all the combinations.

This article was used for distance calculation using longitude and latitude. It is a function that takes the longitude and latitude of two points as arguments (note that each is a list type!) And returns the distance between them in km.

def cal_rho(lon_a,lat_a,lon_b,lat_b):
    ra=6378.140  # equatorial radius (km)
    rb=6356.755  # polar radius (km)
    F=(ra-rb)/ra # flattening of the earth
    rad_lat_a=np.radians(lat_a)
    rad_lon_a=np.radians(lon_a)
    rad_lat_b=np.radians(lat_b)
    rad_lon_b=np.radians(lon_b)
    pa=np.arctan(rb/ra*np.tan(rad_lat_a))
    pb=np.arctan(rb/ra*np.tan(rad_lat_b))
    xx=np.arccos(np.sin(pa)*np.sin(pb)+np.cos(pa)*np.cos(pb)*np.cos(rad_lon_a-rad_lon_b))
    c1=(np.sin(xx)-xx)*(np.sin(pa)+np.sin(pb))**2/np.cos(xx/2)**2
    c2=(np.sin(xx)+xx)*(np.sin(pa)-np.sin(pb))**2/np.sin(xx/2)**2
    dr=F/8*(c1-c2)
    rho=ra*(xx+dr)
    return rho

After that, I will put all the combinations into this function.

#Prepare a container for the list as an argument
loc_A = []
lon_A = []
lat_A = []

loc_B = []
lon_B = []
lat_B = []

#Intently
for i in range(47):
    for j in range(47):
        loc_A.append(i)
        lon_A.append(df["long"][i])
        lat_A.append(df["lati"][i])
        
        loc_B.append(j)
        lon_B.append(df["long"][j])
        lat_B.append(df["lati"][j])

#Calculate!
rho=cal_rho(lon_A,lat_A,lon_B,lat_B)

#Put the result in a data frame
combi_df = pd.DataFrame([loc_A,loc_B,lon_A,lat_A,lon_B,lat_B,rho]).T
combi_df.columns = ["loc_A","loc_B","long_A","lati_A","long_B","lati_B","Dist/km"]
combi_df = combi_df.fillna(0)
combi_df.head(5)

As a result, we were able to calculate the distance between each point. image.png

Calculation of shortest path

This is almost the same implementation as Previous article and Original article, but only the size of drawing by matplotlib is changed (Reference article). ) Click here for the drawing function.

import random, numpy as np, pandas as pd, networkx as nx, matplotlib.pyplot as plt
from itertools import chain, combinations
from pulp import *

def draw(g):
    plt.figure(dpi=150)
    rn = g.nodes() #Node list
    pos = nx.spring_layout(g) #Node position
    """drawing"""
    nx.draw_networkx_labels(g, pos=pos)
    nx.draw_networkx_nodes(g, node_color='w', pos=pos)
    nx.draw_networkx_edges(g, pos=pos)
    plt.show()

Click here for the optimization calculation code. For Networkx, this article is easy to read! This time, I solved it as a discrete optimization problem with possible values ​​of 0 and 1 ([Reference article]([Reference article]) https://docs.pyq.jp/python/math_opt/pulp.html))。

#Define 47 points
n = 47 #Number of nodes
g = nx.random_graphs.fast_gnp_random_graph(n, 1, 8)

#Input the distance between 47 points
for a,b in g.edges():
    col_idx = a * 47 + b
    g.edges[a, b]["dist"] = combi_df.loc[col_idx,"Dist/km"]

#Calculation
source, sink = 0, 46 #start point,end point
r = list(enumerate(g.edges()))
m = LpProblem() #Mathematical model
x = [LpVariable('x%d'%k, lowBound=0, upBound=1,cat="Integer" ) for k, (i, j) in r] #variable(Whether to enter the road)
m += lpSum(x[k] * g.edges[i,j]["dist"] for k, (i, j) in r) #Objective function
for nd in g.nodes():
    m += lpSum(x[k] for k, (i, j) in r if i == nd) \
        + lpSum(x[k] for k, (i, j) in r if j == nd) == {source:1, sink:1}.get(nd, 2) #Constraint
m.solve()

print([(i, j) for k, (i, j) in r if value(x[k]) == 1])
print(LpStatus[m.status])
print(value(m.objective))

The first print is the display of the route, the second is the display of the optimization calculation, and the third is the display of the optimum solution. Click here for the output result. image.png Since it is optimal, it seems that the optimal solution could be found! !! And the shortest route is __4304km! !! __ So Takashi-kun seems to be able to travel all over Japan in __430 hours, that is, about 17 days and 22 hours! __

__ I did it! I was able to calculate! !! !! __ It was a moment when I thought it was __. __

Circular path problem

What kind of route did you take? I thought, I tried drawing. Click here for the code to check the optimization result.

G = nx.Graph()
for k,(i,j) in r:
    if value(x[k]) ==1:
        nx.add_path(G, [i, j])
draw(G)

Here is the drawing result. The __circle generation problem __, which was also mentioned in the previous article (https://qiita.com/canonno/items/711201febd5b2bf746c4), occurred. This is a problem that arises from the lack of constraint that "when focusing on a certain point, two roads extend at the waypoint and one road extends at the start and end points."

I borrowed a blank map from this site and plotted it. There are places in the Tohoku, Kanto, and Chugoku regions that go round and round.

From here, I will play with the route a little and make it a single route. Since arbitrary operations are added, it is difficult to know whether it will be the optimum solution, but let's think that it will be close.

Akita-Niigata

Try to make the length of the path you want to cut extremely long. It's like this.

pair_list = [(4,14),(14,4)]
for a,b in pair_list:
    col_idx = a * 47 + b
    combi_df.loc[col_idx,"Dist/km"] = 1000000

Now let's do the optimization calculation. The result is only a map. Successful integration of the Tohoku region! !! I did it. Even in this case, it is about 4309km, so changing only one route does not seem to have that much effect.

Cut between Niigata and Nagano

I was worried about what to do with the integration of the Kanto region. So, I thought about a route that goes from Niigata to Gunma without going to Nagano, joins the Kanto circle, and then goes to Nagano. So I cut the route between Niigata and Nagano.

The result is here.

I went from Niigata to Tochigi and was able to integrate well. The distance has increased a little, but it is now 4326km.

Cut between Shimane and Okayama

Finally, what should we do with the integration of the Chugoku region? You can see the movement between Shimane and Okayama, but if you want to move between Okayama and Shimane, it seems better to go via Tottori or Hiroshima. So I will cut between Shimane and Okayama.

__ I got involved in Kagawa ... __ __ How much do you want to be small there? __

Cut between Shimane and Kagawa

So I will cut the Shimane-Kagawa section. The result is here! Although it is warping from Shimane to Kochi, it is finally one! !! The distance is 4361km.

I was a little worried about the warp between Shimane and Kochi, so I cut it off here as well. This time it warps from Kochi to Oita, and Yamaguchi seems to go to Kyushu and then pick it up. The distance is about the same as 4361km. It seems difficult to write a single stroke around Chugoku and Shikoku.

Conclusion

__ Distance to run: Approximately 4361km__ __ Running time: Approximately 18 days 4 hours __

at the end

Optimization problem interesting! I will continue to experiment. Thank you for reading to the end, and it will be encouraging to have LGTM!

Recommended Posts

I ran all over Japan at the shortest distance at 10km / h.
I want to find the shortest route to travel through all points