# -*- 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("")