[PYTHON] ○○ Solving problems in the Department of Mathematics by optimization

what is this

["Takeshi no Koma University Mathematics Department"](https://ja.wikipedia.org/wiki/%E3%81%9F%E3%81%91%E3%81%97%E3%81%AE%E3% 82% B3% E3% 83% 9E% E5% A4% A7% E6% 95% B0% E5% AD% A6% E7% A7% 91) I tried the one that seems to be solved by optimization. It was.

Preparing to run in Python

`python`


import scipy.optimize, numpy as np, networkx as nx
from itertools import product
from pulp import *
from ortoolpy import addvar, addvars, addbinvar, addbinvars

Import it first.

Assignment problem

Assign baseball players A to E to positions to maximize the strength of the team (the lower the number, the better).

player	1st base	2nd base	3rd base	SS	outfield
A	6	5	6	5	6
B	8	7	6	8	7
C	4	5	4	4	5
D	6	7	6	4	7
D	10	8	10	7	10

Minimum weight perfect matching problem That's right.

`python`


w = np.array([[6,5,6,5,6],
              [8,7,6,8,7],
              [4,5,4,4,5],
              [6,7,6,4,7],
              [10,8,10,7,10]])
N = len(w)
g = nx.Graph()
for i,j in product(range(N),range(N)):
    g.add_edge(i, j+N, weight=w.max()-w[i][j])
r = nx.max_weight_matching(g)
for i in range(N):
    print('ABCDE'[i], ['1st base','2nd base','3rd base','SS','outfield',][r[i]-N])
>>>
A outfield
B 3rd base
C 1st base
D SS
E 2nd base

I got the same answer.

Meeting place problem

Where is the total travel distance minimized when 1 to 5 people gather in one place?

[Nonlinear optimization problem](bfbf4c185ed7004b5721 #% E9% 9D% 9E% E7% B7% 9A% E5% BD% A2% E6% 9C% 80% E9% 81% A9% E5% 8C% 96% E5% 95% It will be 8F% E9% A1% 8C).

Considering the distance in L2 norm, [Weber problem](http://www.orsj.or.jp/~wiki/wiki/index.php/%E3%82%A6%E3%82%A7%E3%83 % BC% E3% 83% 90% E3% 83% BC% E5% 95% 8F% E9% A1% 8C) That's right.

`python`


pos = np.array([[1,2,1],[3,7,5],[5,3,2],[6,5,2],[7,1,3]]) # x,y,Number of people
def func(x):
    return sum(p[2]*np.linalg.norm(p[:2]-x) for p in pos)
print(scipy.optimize.fmin(func, [0,0]))
>>>
[ 4.80539237  4.38561398]

Since it is a grid pattern, I will also calculate the L1 norm.

`python`


m = LpProblem()
vx = addvar(cat=LpInteger)
vy = addvar(cat=LpInteger)
vp = addvars(pos.shape[0])
vq = addvars(pos.shape[0])
m += lpDot(pos[:,2], vp) + lpDot(pos[:,2], vq)
for i in range(pos.shape[0]):
    m += vp[i] >=   pos[i,0]-vx
    m += vp[i] >= -(pos[i,0]-vx)
    m += vq[i] >=   pos[i,1]-vy
    m += vq[i] >= -(pos[i,1]-vy)
m.solve()
print(value(vx),value(vy))
>>>
5.0 5.0

Chinese postal delivery problem

Depart from the post office, deliver mail to all homes, and return to the post office again, finding a route that minimizes the distance traveled during this time.

Problems around all sides [Chinese postal delivery problem](https://ja.wikipedia.org/wiki/%E4%B8%AD%E5%9B%BD%E4%BA%BA%E9%83 % B5% E4% BE% BF% E9% 85% 8D% E9% 81% 94% E5% 95% 8F% E9% A1% 8C). One-stroke writing (Euler cycle) If possible, the smallest is self-evident. Since you can write with one stroke by eliminating the points of odd order, you can solve the Minimum weight perfect matching problem by weighting the distance only with the points of odd order. The Eulerian cycle is found in nx.eulerian_circuit. The program is omitted.

Shoelace problem

Find the shortest length (the length to the knot) when you put the laces through the shoes with 16 holes, 8 in each row, so that you can wear them properly! Each hole must be connected to the hole in the opposite row.

Make a graph.

`python`


N = 8
g = nx.Graph()
for i in range(N):
    g.add_edge(i,i+N)
    if i<N-1:
        g.add_edge(i,i+N+1)
    if i:
        g.add_edge(i,i+N-1)
        g.add_edge(i-1,i)
        g.add_edge(i+N-1,i+N)
nx.draw_networkx(g)

I will try to solve it.

`python`


def len_(i,j):
    k =abs(i-j)
    return 1 if k==1 else 2 if k==N else 2.236
m = LpProblem()
#Variable creation
for i,j in g.edges():
    g.edge[i][j]['Var'] = addbinvar()
m += lpSum([len_(i,j)*g.edge[i][j]['Var'] for i,j in g.edges()]) #Objective function
for i in range(2*N):
    #There are in and out points
    m += lpSum(dc['Var'] for dc in g.edge[i].values()) == 2
    m += lpSum(dc['Var'] for j,dc in g.edge[i].items()
               if abs(i-j) >= N-1) >= 1 #Connect with the other side
for k in range(1,N-2):
    #Connect
    m += lpSum(g.edge[i][j]['Var'] for i,j in
               [[k,k+1],[k,k+N+1],[k+N,k+1],[k+N,k+N+1]]) == 2
m.solve()
print(value(m.objective))
for i,j in g.edges():
    if value(g.edge[i][j]['Var']) > 0:
        print((i,j), end=', ')
>>>
25.416000000000004
(0, 8), (0, 1), (8, 9), (9, 2), (1, 10), (10, 11),
(2, 3), (11, 4), (3, 12), (12, 13), (4, 5), (13, 6),
(5, 14), (14, 15), (6, 7), (15, 7),

I got the same answer.

Kakuro

There are 5 cards with numbers from 0 to 9 written on the front and back (there are no duplicate numbers). The total number on the front is "19", and the three from the left are turned over and the total number of visible numbers is "20". Similarly, it was "35" for the front and back front and back, "11" for the front and back front and back, and "31" for the front and back front and back. Answer the numbers on the front of the first card you saw in order!

`python`


hint = [[(1,1,1,1,1),19],
        [(0,0,0,1,1),20],
        [(0,0,1,0,0),35],
        [(1,1,0,1,0),11],
        [(1,0,1,0,1),31]]
r = range(10)
m = LpProblem()
v = addbinvars(10,10) # i:position,j:Numbers
for i in r:
    m += lpSum(v[i]) == 1
    m += lpSum(v[j][i] for j in r) == 1
for ptn,n in hint:
    m += lpSum(lpDot(r,v[i if j else i+5])
        for i,j in enumerate(ptn)) == n
m.solve()
print((np.vectorize(value)(v)@r).reshape(2,-1))
>>>
[[ 3.  1.  9.  2.  4.]
 [ 6.  8.  0.  7.  5.]]

The same answer.

Bill

Building It's a puzzle.

`python`


n = 5
hint = [([0,3,0,3,0], 1, 0,   0,   0,  0,  1),
        ([0,2,0,0,4], 1, 0,   0, n-1,  0, -1),
        ([0,0,0,4,0], 0, 1,   0,   0,  1,  0),
        ([0,4,0,0,0], 0, 1, n-1,   0, -1,  0)]
m = LpProblem()
v = np.array(addbinvars(n, n, n))
r = np.array(addvars(n, n))
def add(m, r, k, p, q, y, x):
    if k==0:
        return
    u = addbinvars(n-1)
    m += lpSum(u) == k - 1
    vmx = r[p,q]
    for i in range(1,n):
        vnx = r[p + y*i][q + x*i]
        m += vmx + n * u[i-1] >= vnx + 1
        m += vmx + 1 <= vnx + n - n * u[i-1]
        vtm = addvar()
        m += vmx <= vtm
        m += vnx <= vtm
        vmx = vtm
    m += vmx <= n
for i in range(n):
    for j in range(n):
        m += lpSum(v[i,j,:]) == 1
        m += lpDot(range(n), v[i,j]) + 1 == r[i,j]
        m += lpSum(v[i,:,j]) == 1
        m += lpSum(v[:,i,j]) == 1
    for h in hint:
        add(m, r, h[0][i], h[1]*i+h[3], h[2]*i+h[4], h[5], h[6])
m.solve()
np.vectorize(value)(r).astype(int).T.tolist()
>>>
[[4, 1, 3, 2, 5],
 [5, 4, 1, 3, 2],
 [3, 5, 2, 1, 4],
 [1, 2, 4, 5, 3],
 [2, 3, 5, 4, 1]]

This is the same answer.

that's all