Hopfieldネットワーク( TSP )

# -*- coding: UTF-8 -*-
from math import *
import numpy as np

############################################
# ルンゲ・クッタ法  dx/dt=f(t,x)
#      time : 現在の時間
#      h : 時間刻み幅
#      x : 現在の状態
#      dx : 微係数(f(t,x):subで計算する)
#      g : 作業域(g[4][n])
#      n : 微分方程式の次数
#      sub : 微係数を計算する関数の名前
#      return : time+h
#      coded by Y.Suganuma
############################################

def rkg(time, h, x, dx, g, n, sub) :

	h2 = 0.5 * h
	sub(time, x, dx)
	for i1 in range(0, n) :
		g[0][i1] = h * dx[i1]

	time += h2
	for i1 in range(0, n) :
		g[1][i1] = x[i1] + 0.5 * g[0][i1]
	sub(time, g[1], dx)
	for i1 in range(0, n) :
		g[1][i1] = h * dx[i1]

	for i1 in range(0, n) :
		g[2][i1] = x[i1] + 0.5 * g[1][i1]
	sub(time, g[2], dx)
	for i1 in range(0, n) :
		g[2][i1] = h * dx[i1]

	time += h2
	for i1 in range(0, n) :
		g[3][i1] = x[i1] + g[2][i1]
	sub(time, g[3], dx)
	for i1 in range(0, n) :
		g[3][i1] = h * dx[i1]

	for i1 in range(0, n) :
		x[i1] = x[i1] + (g[0][i1] + 2.0 * g[1][i1] + 2.0 * g[2][i1] + g[3][i1]) / 6.0

	return time

----------------------------------

# -*- coding: UTF-8 -*-
import numpy as np
import sys
from math import *
from random import *
from function import rkg

##################################
# Hopfieldネットワーク( TSP )
#      Coded by Y.Suganuma
##################################

##########################
# クロネッカのデルタ関数
##########################
def delta(i, j) :

	if i < 0 :
		i = n - 1
	elif i >= n :
		i = 0
	if j < 0 :
		j = n - 1
	elif j >= n :
		j = 0
	k = 0
	if i == j :
		k = 1
	return k

##################
# ユニットの出力
##################
def out(x) :
	return 0.5 * (1.0 + tanh(x))
#	return 1.0 / (1.0 + exp(-x))

################
# 微係数の計算
################
def snx(time, u, du) :

	for x in range(0, n) :
		for i in range(0, n) :
			k     = n * x + i
			du[k] = 0.0
			for j in range(0, n) :
				if j != i :
					du[k] -= A * out(u[n*x+j])
			for y in range(0, n) :
				if y != x :
					du[k] -= B * out(u[n*y+i])
			N = 0.0
			for xx in range(0, n) :
				for ii in range(0, n) :
					N += out(u[n*xx+ii])
			du[k] -= C * (N - n)
			for y in range(0, n) :
				m1 = (i + 1) % n
				m2 = i - 1
				if m2 < 0 :
					m2 = n - 1
				du[k] -= D * d[x][y] * (out(u[n*y+m1]) + out(u[n*y+m2]))

n = 4
A = 10.0
B = 10.0
C = 10.0
D = 2.0
d = np.empty((4, 4), np.float)
					# 距離
r = sqrt(2.0)
for i1 in range(0, n) :
	for i2 in range(0, n) :
		if i1 == i2 :
			d[i1][i2] = 0.0
		else :
			d[i1][i2] = 1.0
d[0][2] = r
d[1][3] = r
d[2][0] = r
d[3][1] = r
					# 初期状態
u = np.empty((4, 4), np.float)
for i1 in range(0, n) :
	for i2 in range(0, n) :
		u[i1][i2] = 0.1 * random() - 0.05

print("初期状態(出力):")
for i1 in range(0, n) :
	for i2 in range(0, n) :
		print(" " + str(out(u[i1][i2])), end="")
	print("")
					# 更新
time = 0.0
dx   = np.empty(16, np.float)
g    = np.empty((4, 16), np.float)
uu   = np.reshape(u, 16)
for i1 in range(0, 100) :
	time = rkg(time, 1, uu, dx, g, n*n, snx)

print("最終状態(出力):")
for i1 in range(0, n) :
	for i2 in range(0, n) :
		print(" " + str(out(u[i1][i2])), end="")
	print("")