六面体の頂点をたどる(6)

六面体で同じように行列を作る。
ここで、行列の各成分を3倍すれば整数化できるので、その行列の固有値を考える。
この固有方程式は、

|A - λE| =
15120λ3+10080λ4-137484λ5-8916λ6+477972λ7
-192860λ8-753221λ9+591113λ10+497726λ11-643343λ12
-78020λ13+336853λ14-67406λ15-88995λ16+39366λ17
+9643λ18-8702λ19+419λ20+848λ21-181λ22
-22λ23+11λ2425 = 0

となった。
これを昔作った因数分解のプログラムにかけると、

3(1-λ)(2-λ)2(3-λ)(1+2λ-λ2)(2-λ2)(3-λ2)2
(2-5x24)(35+35x-47x2-12x3+13x456) = 0

となった。
最後の因子については、固有値にを具体的に表示できない。



Array.prototype.add = function(b) {
for(var i = 0; i < b.length; i++) {
if(this[i] == undefined)
this[i] = b[i];
else
this[i] += b[i];
}
return this;
}

Array.prototype.subtract = function(b) {
for(var i = 0; i < b.length; i++) {
if(this[i] == undefined)
this[i] = -b[i];
else
this[i] -= b[i];
}
return this;
}

Array.prototype.multiply = function(b) {
var a = [ ];
for(var i = 0; i < this.length; i++) {
for(var j = 0; j < b.length; j++) {
var k = i + j;
if(a[k] == undefined)
a[k] = this[i] * b[j];
else
a[k] += this[i] * b[j];
}
}
return a;
}

Array.prototype.divide = function(b) {
var a = [ ];
var r = this.slice(0);
var m = b[b.length-1];
for(var i = r.length - b.length; i >= 0; i--) {
var d = r[i+b.length-1] / m;
a[i] = d;
for(var j = 0; j < b.length; j++) {
var k = i + j;
r[k] -= d * b[j];
}
}
WScript.Echo("a : " + a);
return a;
}

Array.prototype.value = function(x) {
var v = this[this.length-1];
for(var i = this.length - 2; i >= 0; i--)
v = v * x + this[i];
return v;
}

Array.prototype.toStringEq = function() {
var str = "";
for(var i = 0; i < this.length; i++) {
var a = this[i];
if(a == 0)
continue;

if(i == 0) {
str += "" + a;
}
else {
if(a == 1) {
if(str == "")
str += "λ";
else
str += "+λ";
}
else if(a == -1) {
str += "-λ";
}
else if(a > 0) {
if(str == "")
str += a + "λ";
else
str += "+" + a + "λ";
}
else {
str += a + "λ";
}
if(i != 1) {
str += "" + i + "";
}
}
}

return str;
}

var neighbor = [ 1, 3, 4 ];

var purm = [ ]; // 置換群
make_purms();

var q = [ ]; // 遷移確率 * 3
var set_stat = [ ]; // 状態の集合
set_stat[1] = 1;
var queue_stat = [ 1 ]; // 遷移確率計算待ち行列
while(queue_stat.length != 0) {
var s = queue_stat.shift();
var set_next = [ ];
for(var i = 0; i < 3; i++) {
var s_next = proceed_stat(s, i);
increment(q, s, s_next);
set_next[s_next] = 1;
}

// 次の状態が現れていなかった状態ならキューに入れる
for(var s_next in set_next) {
if(set_stat[s_next] == undefined) {
set_stat[s_next] = 1;
queue_stat.push(s_next);
}
}
}

// 行列を作る
var A = [ ];
make_matrix(A);

// 固有方程式を作る
var eigen_function = make_eigen_function(A);
WScript.Echo(eigen_function.toStringEq());

function make_matrix(M) {
// 状態を昇順にソート
var ary_stat = [ ];
for(var s in q) {
if(!isNaN(s))
ary_stat.push(s - 0);
}
ary_stat.sort(function(a, b) { return a - b; });

var stat_to_index = [ ];
for(var i = 0; i < ary_stat.length; i++)
stat_to_index[ary_stat[i]] = i;

for(var s_next in q) {
var i = stat_to_index[s_next];
M[i] = [ ];
for(var j = 0; j < ary_stat.length; j++)
M[i][j] = 0;
var ary_terms = [ ];
var q_next = q[s_next];
for(var s in q_next) {
if(s - 1 == 0)
continue;
var j = stat_to_index[s];
var n = q_next[s];
M[i][j] = n;
}
WScript.Echo(M[i]);
}
}

function make_eigen_function(M) {
// L = M - λE
var L = [ ];
for(var i = 0; i < M.length; i++) {
L[i] = [ ];
for(var j = 0; j < M.length; j++) {
if(i == j)
L[i][j] = [ M[i][j], -1 ];
else if(M[i][j] == 0)
L[i][j] = [ ];
else
L[i][j] = [ M[i][j] ];
}
}

return calc_determinant(L);
}

function calc_determinant(M) {
if(M.length == 1)
return M[0][0];

var f = [ ];
for(var i = 0; i < M.length; i++) {
if(M[0][i].length == 0)
continue;
var minor_det = calc_determinant(make_minor_matrix(M, i));
if(i % 2 == 0)
f.add(M[0][i].multiply(minor_det));
else
f.subtract(M[0][i].multiply(minor_det));
}

return f;
}

// 0行と指定した列を除いた小行列を作る
function make_minor_matrix(M, c0) {
var L = [ ];
for(var i = 0; i < M.length - 1; i++) {
L[i] = [ ];
var c = 0;
for(var j = 0; j < M.length; j++) {
if(j == c0)
continue;
L[i][c] = M[i+1][j];
c++;
}
}
return L;
}

// 方程式の解を求める
function solve_equation(f) {
var s = [ ];
while(f.length > 1) {
// 剰余定理で根を一つ求める
var a = search_root(f);
WScript.Echo("aa : " + a);
if(a == null)
WScript.Quit();
f = f.divide([a, -1]);
s.push(a);
}

WScript.Echo("s : " + s);
return s;
}

// 剰余定理で根を一つ求める
function search_root(f) {
var limit = Math.abs(f[0]);
var n = 0;
WScript.Echo("f(1) : " + f.value(1));
WScript.Echo("f(2) : " + f.value(2));
WScript.Echo("f(3) : " + f.value(3));
WScript.Echo("f(4) : " + f.value(4));
do {
if(f.value(n) == 0)
return n;
if(n > 0)
n = -n;
else
n = 1 - n;
} while(n != -limit);

return null;
}

// 次の状態を正規化して返す
function proceed_stat(s, i) {
var pt_next = neighbor[i];
var next_stat = s | (1 << pt_next); // 置換前の次の状態
var min_stat = 256; // 置換後の次の状態
for(var k = 0; k < purm.length; k++) {
// 移動先の点が0に写る置換だけ選ぶ
if(purm[k][pt_next] == 0) {
var s_candidate = map(purm[k], next_stat);
if(s_candidate < min_stat)
min_stat = s_candidate;
}
}

return min_stat;
}

function increment(ary, s, s_next) {
if(ary[s_next] == undefined)
ary[s_next] = [ ];
if(ary[s_next][s] == undefined)
ary[s_next][s] = 1;
else
ary[s_next][s]++;
}

function map(sigma, s) {
var s_next = 0;
for(var i = 0; i < 8; i++) {
if(s & 1 << i)
s_next |= 1 << sigma[i];
}
return s_next;
}

function make_purms() {
var rot1 = [ 3, 0, 1, 2, 7, 4, 5, 6 ];
var sym1 = [ 4, 5, 6, 7, 0, 1, 2, 3 ];
var rot2 = [ 0, 4, 5, 1, 3, 7, 6, 2 ];
var sym2 = [ 0, 4, 7, 3, 1, 5, 6, 2 ];
var e = [ 0, 1, 2, 3, 4, 5, 6, 7 ];

for(var i = 0; i < 8; i++) {
var sigma = e;
for(var k = 0; k < (i & 3); k++)
sigma = synthesize(rot1, sigma);
if(i >> 2)
sigma = synthesize(sym1, sigma);
for(var j = 0; j < 6; j++) {
var sigma2 = sigma;
for(var k = 0; k < j % 3; k++)
sigma2 = synthesize(rot2, sigma2);
if(j >= 3)
sigma2 = synthesize(sym2, sigma2);
purm.push(sigma2);
}
}
}

function synthesize(a, b) {
var c = [ ];
for(var i = 0; i < a.length; i++)
c[i] = a[b[i]];
return c;
}