奇妙なことに、私は数学の学位を取得している学部生の夏を過ごして、この問題を研究し、それに対処するアルゴリズムを考え出しました。まず、他の回答についてコメントします。
Martin B: これをハミルトニアン パスの問題として正しく識別しました。ただし、グラフが通常のグリッドである場合 (コメントで説明したように)、ハミルトニアン パスを自明に見つけることができます (たとえば、行優先順序を蛇行させます)。
agnorest: ハミルトニアン パス問題が難しい問題のクラスにあることについて正しく説明しました。agnorest はまた、グラフ構造を悪用する可能性をほのめかしましたが、これは面白いことに、通常のグリッドの場合は些細なことです。
ここで、あなたが達成しようとしていることについて詳しく説明しながら、議論を続けます。コメントで述べたように、規則的な格子/グリッド上で交差しない「空間充填」の特定のクラスを見つけるのは簡単です。ただし、これらだけでは満足できず、グリッドをランダムにカバーする「興味深い」ウォークを見つけるアルゴリズムを見つけたいと考えているようです。しかし、その可能性を探る前に、これらのウォークの「交差しない」という特性が非常に重要であり、それらを列挙するのが難しい理由を指摘したいと思います。
実は、私が夏のインターンシップで学んだことは、Self Avoiding Walkと呼ばれるものです。驚くべきことに、SAW (自己回避歩行) の研究は、物理学におけるモデリングのいくつかのサブドメインにとって非常に重要です (そして、マンハッタン プロジェクトの重要な要素でした!) あなたが質問で与えたアルゴリズムは、実際にはアルゴリズムの変形です。 「成長」アルゴリズムまたはRosenbluthのアルゴリズムとして知られています( Marshal Rosenbluth以外にちなんで名付けられました)。このアルゴリズムの一般的なバージョン (SAW の統計を推定するために使用される) と物理学との関係の両方に関する詳細は、この論文のような文献で容易に見つけることができます。
2 次元の SAW は、研究が難しいことで有名です。2 次元の SAW に関する理論的な結果はほとんど知られていません。奇妙なことに、3 次元以上では、成長定数など、SAW のほとんどの特性が理論的に知られていると言えます。2 次元の SAW は非常に興味深いものです。
目の前の問題について話し合うためにこの議論に君臨すると、おそらく、成長アルゴリズムの実装が頻繁に「切断」され、ラティス全体を埋めるように拡張できないことに気付くでしょう。この場合、問題をハミルトニアン パス問題と見なす方が適切です。興味深いハミルトニアン パスを見つけるための私のアプローチは、問題を整数線形計画法として定式化し、ILP で使用される修正ランダム エッジを追加することです。ランダムなエッジを修正すると、生成プロセスにランダム性が与えられ、ILP 部分は、特定の構成が実現可能かどうかを効率的に計算し、可能であれば解を返します。
コード
次のコードは、任意の初期条件のハミルトニアン パスまたはサイクル問題を実装します。これは、4 連結性を備えた通常の 2D ラティスに実装されています。定式化は、標準のサブツアー除去 ILP ラグランジュに従います。サブツアーは遅延して削除されます。つまり、多くの反復が必要になる可能性があります。
これを拡張して、問題に対して「興味深い」と思われるランダムまたはその他の初期条件を満たすことができます。初期条件が実行不可能な場合は、早期に終了してこれを出力します。
このコードはNetworkXとPuLPに依存しています。
"""
Hamiltonian path formulated as ILP, solved using PuLP, adapted from
https://projects.coin-or.org/PuLP/browser/trunk/examples/Sudoku1.py
Authors: ldog
"""
# Import PuLP modeler functions
from pulp import *
# Solve for Hamiltonian path or cycle
solve_type = 'cycle'
# Define grid size
N = 10
# If solving for path a start and end must be specified
if solve_type == 'path':
start_vertex = (0,0)
end_vertex = (5,5)
# Assuming 4-connectivity (up, down, left, right)
Edges = ["up", "down", "left", "right"]
Sequence = [ i for i in range(N) ]
# The Rows and Cols sequences follow this form, Vals will be which edge is used
Vals = Edges
Rows = Sequence
Cols = Sequence
# The prob variable is created to contain the problem data
prob = LpProblem("Hamiltonian Path Problem",LpMinimize)
# The problem variables are created
choices = LpVariable.dicts("Choice",(Vals,Rows,Cols),0,1,LpInteger)
# The arbitrary objective function is added
prob += 0, "Arbitrary Objective Function"
# A constraint ensuring that exactly two edges per node are used
# (a requirement for the solution to be a walk or path.)
for r in Rows:
for c in Cols:
if solve_type == 'cycle':
prob += lpSum([choices[v][r][c] for v in Vals ]) == 2, ""
elif solve_type == 'path':
if (r,c) == end_vertex or (r,c) == start_vertex:
prob += lpSum([choices[v][r][c] for v in Vals]) == 1, ""
else:
prob += lpSum([choices[v][r][c] for v in Vals]) == 2, ""
# A constraint ensuring that edges between adjacent nodes agree
for r in Rows:
for c in Cols:
#The up direction
if r > 0:
prob += choices["up"][r][c] == choices["down"][r-1][c],""
#The down direction
if r < N-1:
prob += choices["down"][r][c] == choices["up"][r+1][c],""
#The left direction
if c > 0:
prob += choices["left"][r][c] == choices["right"][r][c-1],""
#The right direction
if c < N-1:
prob += choices["right"][r][c] == choices["left"][r][c+1],""
# Ensure boundaries are not used
for c in Cols:
prob += choices["up"][0][c] == 0,""
prob += choices["down"][N-1][c] == 0,""
for r in Rows:
prob += choices["left"][r][0] == 0,""
prob += choices["right"][r][N-1] == 0,""
# Random conditions can be specified to give "interesting" paths or cycles
# that have this condition.
# In the simplest case, just specify one node with fixed edges used.
prob += choices["down"][2][2] == 1,""
prob += choices["up"][2][2] == 1,""
# Keep solving while the smallest cycle is not the whole thing
while True:
# The problem is solved using PuLP's choice of Solver
prob.solve()
# The status of the solution is printed to the screen
status = LpStatus[prob.status]
print "Status:", status
if status == 'Infeasible':
print 'The set of conditions imposed are impossible to solve for. Change the conditions.'
break
import networkx as nx
g = nx.Graph()
g.add_nodes_from([i for i in range(N*N)])
for r in Rows:
for c in Cols:
if value(choices['up'][r][c]) == 1:
nr = r - 1
nc = c
g.add_edge(r*N+c,nr*N+nc)
if value(choices['down'][r][c]) == 1:
nr = r + 1
nc = c
g.add_edge(r*N+c,nr*N+nc)
if value(choices['left'][r][c]) == 1:
nr = r
nc = c - 1
g.add_edge(r*N+c,nr*N+nc)
if value(choices['right'][r][c]) == 1:
nr = r
nc = c + 1
g.add_edge(r*N+c,nr*N+nc)
# Get connected components sorted by length
cc = sorted(nx.connected_components(g), key = len)
# For the shortest, add the remove cycle condition
def ngb(idx,v):
r = idx/N
c = idx%N
if v == 'up':
nr = r - 1
nc = c
if v == 'down':
nr = r + 1
nc = c
if v == 'left':
nr = r
nc = c - 1
if v == 'right':
nr = r
nc = c + 1
return nr*N+c
prob += lpSum([choices[v][idx/N][idx%N] for v in Vals for idx in cc[0] if ngb(idx,v) in cc[0] ]) \
<= 2*(len(cc[0]))-1, ""
# Pretty print the solution
if len(cc[0]) == N*N:
print ''
print '***************************'
print ' This is the final solution'
print '***************************'
for r in Rows:
s = ""
for c in Cols:
if value(choices['up'][r][c]) == 1:
s += " | "
else:
s += " "
print s
s = ""
for c in Cols:
if value(choices['left'][r][c]) == 1:
s += "-"
else:
s += " "
s += "X"
if value(choices['right'][r][c]) == 1:
s += "-"
else:
s += " "
print s
s = ""
for c in Cols:
if value(choices['down'][r][c]) == 1:
s += " | "
else:
s += " "
print s
if len(cc[0]) != N*N:
print 'Press key to continue to next iteration (which eliminates a suboptimal subtour) ...'
elif len(cc[0]) == N*N:
print 'Press key to terminate'
raw_input()
if len(cc[0]) == N*N:
break