Hopfieldネットワーク( TSP )

/*********************************/
/* Hopfieldネットワーク( TSP ) */
/*      Coded by Y.Suganuma)     */
/*********************************/
import java.io.*;
import java.util.*;

public class Hopfield_c {

	static int n = 4;
	static double A = 10.0, B = 10.0, C = 10.0, D = 2.0, d[][] = new double [n][n];

	public static void main(String args[])
	{
					// 距離
		double r = Math.sqrt(2.0);
		for (int i1 = 0; i1 < n; i1++) {
			for (int i2 = 0; i2 < n; i2++) {
				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;
					// 初期状態
		double u[][] = new double [n][n];
		Random rn = new Random();
		for (int i1 = 0; i1 < n; i1++) {
			for (int i2 = 0; i2 < n; i2++)
				u[i1][i2] = 0.1 * rn.nextDouble() - 0.05;
		}

		Console con = System.console();
		con.printf("初期状態(出力):\n");
		for (int i1 = 0; i1 < n; i1++) {
			for (int i2 = 0; i2 < n; i2++)
				con.printf("%6.3f", out(u[i1][i2]));
			con.printf("\n");
		}
					// 更新
		double time = 0.0;
		double dx[] = new double[n*n];
		double uu[] = new double[n*n];
		double g[][] = new double [4][n*n];
		int k = 0;
		for (int i1 = 0; i1 < n; i1++) {
			for (int i2 = 0; i2 < n; i2++) {
				uu[k] = u[i1][i2];
				k++;
			}
		}
		for (int i1 = 0; i1 < 100; i1++)
			time = rkg(time, 1, uu, dx, g, n*n);

		con.printf("最終状態(出力):\n");
		for (int i1 = 0; i1 < n; i1++) {
			for (int i2 = 0; i2 < n; i2++)
				con.printf("%6.3f", out(uu[i1*n+i2]));
			con.printf("\n");
		}
	}

	/**************************/
	/* クロネッカのデルタ関数 */
	/**************************/
	static int delta(int i, int j)
	{
		if (i < 0)
			i = n - 1;
		else if (i >= n)
			i = 0;
		if (j < 0)
			j = n - 1;
		else if (j >= n)
			j = 0;
		int k = 0;
		if (i == j)
			k = 1;
		return k;
	}

	/******************/
	/* ユニットの出力 */
	/******************/
	static double out(double x)
	{
		return 0.5 * (1.0 + Math.tanh(x));
//		return 1.0 / (1.0 + Math.exp(-x));
	}

	/****************/
	/* 微係数の計算 */
	/****************/
	static void snx(double time, double u[], double du[])
	{
		for (int x = 0; x < n; x++) {
			for (int i = 0; i < n; i++) {
				int k = n * x + i;
				du[k] = 0.0;
				for (int j = 0; j < n; j++) {
					if (j != i)
						du[k] -= A * out(u[n*x+j]);
				}
				for (int y = 0; y < n; y++) {
					if (y != x)
						du[k] -= B * out(u[n*y+i]);
				}
				double N = 0.0;
				for (int xx = 0; xx < n; xx++) {
					for (int ii = 0; ii < n; ii++)
						N += out(u[n*xx+ii]);
				}
				du[k] -= C * (N - n);
				for (int y = 0; y < n; y++) {
					int m1 = (i + 1) % n;
					int 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]));
				}
			}
		}
	}

	/*************************************************/
	/* ルンゲ・クッタ法  dx/dt=f(t,x)                */
	/*      time : 現在の時間                        */
	/*      h : 時間刻み幅                           */
	/*      x : 現在の状態                           */
	/*      dx : 微係数(f(t,x):snxで計算する)      */
	/*      g : 作業域(g[4][n])                    */
	/*      n : 微分方程式の次数                     */
	/*      return : time+h                          */
	/*************************************************/
	static double rkg(double time, double h, double x[], double dx[], double g[][], int n)
	{
		int i1;
		double h2;

		h2 = 0.5 * h;

		snx(time, x, dx);
		for (i1 = 0; i1 < n; i1++)
			g[0][i1] = h * dx[i1];

		time += h2;
		for (i1 = 0; i1 < n; i1++)
			g[1][i1] = x[i1] + 0.5 * g[0][i1];
		snx(time, g[1], dx);
		for (i1 = 0; i1 < n; i1++)
			g[1][i1] = h * dx[i1];

		for (i1 = 0; i1 < n; i1++)
			g[2][i1] = x[i1] + 0.5 * g[1][i1];
		snx(time, g[2], dx);
		for (i1 = 0; i1 < n; i1++)
			g[2][i1] = h * dx[i1];

		time += h2;
		for (i1 = 0; i1 < n; i1++)
			g[3][i1] = x[i1] + g[2][i1];
		snx(time, g[3], dx);
		for (i1 = 0; i1 < n; i1++)
			g[3][i1] = h * dx[i1];

		for (i1 = 0; i1 < n; i1++)
			x[i1] = x[i1] + (g[0][i1] + 2.0 * g[1][i1] + 2.0 * g[2][i1] + g[3][i1]) / 6.0;

		return time;
	}
}