[PYTHON] Try to solve the N Queens problem with SA of PyQUBO

I tried the N Queens problem to deepen my understanding of PyQUBO during the year-end and New Year holidays.

I think PyQUBO is a very convenient QUBO creation tool for using an annealing machine. "Solve the problem of setting only m to 1 out of n with an annealing machine using PyQUBO" https://qiita.com/morimorijap/items/196e7fc86ecff927bf40 I also use it in such cases, but I will try to solve the N queens problem again using PyQUBO.

What is the N queen problem? Wikipedia's 8 queens have 8x8 squares and 8 queens https://ja.wikipedia.org/wiki/Eight queen On the other hand, the mass says NxN problems.

Since the movement of the queen proceeds in eight directions diagonally up, down, left, and right as long as there is no obstruction, it becomes a problem to arrange the queen so that there is no queen in the column, row, or diagonally. In the form of a cost function to solve this on the annealing machine,

Hamiltonian of N-Queen problem, column column, row row, diagonal diagonal

H = \sum_{col}\left(\sum_{i \in col} x_i - 1\right)^2 + \sum_{row}\left(\sum_{i \in row} x_i - 1\right)^2 + \sum_{diag} \left((\sum_{i \in diag} x_i ) (\sum_{i \in diag} x_i - 1 )\right)

It will be.

The expression of the diagonal cost function is different, but it is an article of T-Wave of Tohoku University "Solving the N-Queen problem with a D-Wave machine" https://qard.is.tohoku.ac.jp/T- Wave /? P = 884 Will also be helpful.

In the case of 4x4,

Considering the expression of squares like, PyQUBO expresses an Array like Binary (q [0]) in QUBO format of 0,1.

from pyqubo import Array, Placeholder, solve_qubo, Sum
from pprint import pprint
import numpy as np

First, import what you need, such as PyQUBO. Suppose you want to make a 4x4 square.

# N-Number of Queen trout
N = 4

#Create 0 and 1 variables with pyQUBO for NxN squares
q = Array.create('q', (N*N), 'BINARY')

q_shape = np.reshape(q,(N,N))
print(q_shape)

Then

[[Binary(q[0]) Binary(q[1]) Binary(q[2]) Binary(q[3])] [Binary(q[4]) Binary(q[5]) Binary(q[6]) Binary(q[7])] [Binary(q[8]) Binary(q[9]) Binary(q[10]) Binary(q[11])] [Binary(q[12]) Binary(q[13]) Binary(q[14]) Binary(q[15])]]

In this way, you can create PyQUBO Binary in 4x4 array format. The line part is as follows in the formula, so

\sum_{row}\left(\sum_{i \in row} x_i - 1\right)^2

Expressing this in PyQUBO, you can slice each row, add them as they are, square them, and finally add them up.

#line
row_const = 0.0
for row in range(N):
    row_const += (sum(i for i in q_shape[row, :]) -1)**2

Also, the columns can be represented as follows.

\sum_{col}\left(\sum_{i \in col} x_i - 1\right)^2

#Column
col_const = 0.0
for col in range(N):
    col_const += (sum(i for i in q_shape[:, col]) -1)**2

It will be.

And the diagonal expression of the problem,

\sum_{diag} \left((\sum_{i \in diag} x_i ) (\sum_{i \in diag} x_i - 1 )\right)

Will be represented by PyQUBO, but it can be divided into two parts, "" and "/", by using Numpy conversion. (If you know a better way, please let me know)

#Diagonal
#\Who
diag_const = 0.0
xi = 0.0
xi_1 = 0.0
for k in range(-N+1,N):
    xi += (sum(i for i in np.diag(q_shape,k=k)))
    xi_1 += (sum(i for i in np.diag(q_shape,k=k))) -1 
    diag_const += xi * xi_1

When

#Diagonal
# /Replace for those who,\ do.
diag_const_f = 0.0
xi = 0.0
xi_1 = 0.0
for k in range(-N+1,N):
    xi += (sum(i for i in np.diag(np.fliplr(q_shape),k=k)))
    xi_1 += (sum(i for i in np.diag(np.fliplr(q_shape),k=k))) -1 
    diag_const_f += xi * xi_1

It will be.

By adding all of these, we can express the Hamiltonian of the energy function. Also, use PyQUBO's Placeholder to adjust the parameters and add alpha, beta, and gamma.

#energy(Hamiltonian)Build

alpha = Placeholder("alpha")
beta = Placeholder("beta")
gamma = Placeholder("gamma")

#Parameters
feed_dict = {'alpha': 1.0, 'beta': 1.0, 'gamma': 1.0}

H = alpha * col_const + beta * row_const  + gamma *(diag_const + diag_const_f)

Next, convert to QUBO.

#Compiling QUBO
model = H.compile()
qubo, offset = model.to_qubo(feed_dict=feed_dict)

pprint(qubo)

Just by itself, it will make QUBO.

{('q[0]', 'q[0]'): -11.0, ('q[0]', 'q[10]'): 6.0, ('q[0]', 'q[11]'): 4.0, ('q[0]', 'q[12]'): 2.0, ('q[0]', 'q[13]'): 6.0, ('q[0]', 'q[14]'): 6.0, ('q[0]', 'q[15]'): 6.0, ('q[0]', 'q[1]'): 6.0, ('q[0]', 'q[2]'): 4.0,

・・・ Because it is long, it is omitted on the way ...

('q[8]', 'q[8]'): -19.0, ('q[8]', 'q[9]'): 14.0, ('q[9]', 'q[9]'): -21.0}

is. This can be calculated as it is with SA of PyQUBO.

#Calculated by SA of PyQUBO
solution = solve_qubo(qubo)

#Decoding calculation results with SA of pyQUBO
decoded_solution, broken, energy = model.decode_solution(solution, vartype="BINARY", 
feed_dict=feed_dict)

pprint(solution)

Result is,

{'q[0]': 0, 'q[10]': 0, 'q[11]': 0, 'q[12]': 0, 'q[13]': 0, 'q[14]': 1, 'q[15]': 0, 'q[1]': 1, 'q[2]': 0, 'q[3]': 0, 'q[4]': 0, 'q[5]': 0, 'q[6]': 0, 'q[7]': 1, 'q[8]': 1, 'q[9]': 0}

It becomes ~~, and the conditions are met for the time being, but since Queen may be included in 13, I feel that it is necessary to adjust the parameters. ~~ Postscript: It seems that the optimum solution was not found because the diagonal corner was missing. This time, the solution came out properly.

Extra: Experiment with the chess board attached to Lego's Harry Potter Advent Calendar 2019 (laughs)

For PyQUBO, refer to the official document. https://pyqubo.readthedocs.io/

that's all.

Postscript: 2020/1/3 9:20pm I edited it because the diagonal condition was wrong.