RCCOUT_20060922_0827.mw

> Delta:=evalm((R-B)/2):A:=evalm((R+B)/2): `R`=matrans(R);`B`=matrans(B);

R = vector([4, 3, 1, 2])

B = vector([3, 4, 4, 3])

> print("pi=",pi);

        AA is (numerical) Asup2. Use A2 for two-out subgraph corresponding to the coloring

> eau:=evalm(AA&*uu);idx:={}:for i to NN do if (eau[i]<1) then idx:=idx union {i} fi od:print("not crossed",idx);

eau := vector([1/2, 3/4, 3/4, 3/4, 3/4, 1])

> unassign('t'):map(evalf,eigenvalues(A));map(abs,map(evalf,[eigenvalues(A)]));   map(abs,map(evalf,[eigenvalues(AA)]));

0., 1., -.5000000000, -.5000000000

[0., 1., .5000000000, .5000000000]

[0., 0., 0., .2500000000, .8090169942, .3090169942]

> A2:=evalm((1/2)*(R2+B2)):a2:=evalm(AA+Del2):iszero(A2-a2),"det",det(evalm(J2-AA)),"synced?",det(J2-A2);

true,

> map(abs,map(evalf,[eigenvalues(A2)]));

[1., .5000000000, .5000000000, .5000000000, .5000000000, .5000000000]

> print(P,Q);

matrix([[1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0]]), matrix([[1, 1, 0, 0]])

> rdd:=read2(Del2):for i to nops(rdd[1]) do print(convert(choose(n,2)[i],set),i,rdd[1][i],"               ",rdd[2][i]) od;

{1, 2}, 1, [[6, 6], [

{1, 3}, 2, [[3, 6], [2,

{1, 4}, 3, [[5,

{2, 3}, 4, [[2,

{2, 4}, 5, [[4, 6], [5,

{3, 4}, 6, [[1, 6], [2, 5]],

> unassign('x'):pi2:=evalm(1/x[1]*(linsolve(transpose(J2-A2),vector(NN,0),'r',x)));u2:=evalm((1/x[1]*linsolve(J2-A2,vector(NN,0),'r',x)));

pi2 := vector([1, 0, 0, 0, 0, 2])

u2 := vector([3, 2, 1, 1, 2, 3])

> Omega2:=abel(a2);readcycles(Omega2);

Omega2 := matrix([[1/3, 0, 0, 0, 0, 2/3], [2/9, 0, 0, 0, 0, 4/9], [1/9, 0, 0, 0, 0, 2/9], [1/9, 0, 0, 0, 0, 2/9], [2/9, 0, 0, 0, 0, 4/9], [1/3, 0, 0, 0, 0, 2/3]])

{{1, 6}}

> k:=2:print(seq({i,choose(n,k)[i]},i=1..binomial(n,k)));

{1, [1, 2]}, {2, [1, 3]}, {3, [1, 4]}, {4, [2, 3]}, {5, [2, 4]}, {6, [3, 4]}

> NN:=binomial(n,2):J2:=IdentityMatrix(binomial(n,2)):unassign('e','x'):e:=vector(n):ee:=vector(n,1):ee:=evalm(e):phi:=diag(seq(ee[q],q=1..n)):deta:=multiply(phi,Delta):psi:=sympow(phi,2):

> APART:=ecliffe(evalm(symult(A,J)+symult(deta,J)));

APART := matrix([[1, 1/2-1/2*e[1]*e[2], 1/4-1/4*e[1]*e[3], 1/4+1/4*e[1]*e[4]], [1/2-1/2*e[1]*e[2], 1, 1/4+1/4*e[2]*e[3], 1/4-1/4*e[2]*e[4]], [1/4-1/4*e[1]*e[3], 1/4+1/4*e[2]*e[3], 1, 0], [1/4+1/4*e[1]...

> A2:=evalm(AA+D2):#"xA2",xtend(A2);

> J2A2:=evalm(J2-A2):evalm(J2A2);

matrix([[1, 0, 0, 0, 0, -1/2-1/2*e[1]*e[2]], [0, 3/4+1/4*e[1]*e[3], -1/4-1/4*e[1]*e[3], 0, 0, -1/4-1/4*e[1]*e[3]], [0, 0, 1, -1/4+1/4*e[1]*e[4], -1/4-1/4*e[1]*e[4], -1/4+1/4*e[1]*e[4]], [0, -1/4-1/4*e...

        DJ2 is the determinant with e's

> dt2:=det(J2-AA):"FROM det Asup2",dt2;DJ2:=ecliffit(det(J2A2));gamma:=zcliffsubsmat(DJ2,vector(n,0)):"LOOKING FOR",gamma;

DJ2 := -3/32*e[2]*e[3]*e[1]*e[4]-1/32*e[1]*e[2]-1/32*e[1]*e[3]+1/32*e[1]*e[4]+1/32*e[2]*e[3]-1/32*e[2]*e[4]-1/32*e[3]*e[4]+5/32

> evalf(dt2),"gamma",evalf(gamma),evalf(dt2-gamma);

.1875000000,

> "PI-SIDE",matvec(ecliffe(map(simplify,evalm((linsolve(transpose(J2-AA-Del2),vector(NN,0),'r',x))))));

> "U-SIDE",matvec(ecliffe(map(simplify,evalm((linsolve((J2-AA-Del2),vector(NN,0),'r',x))))));

> V2:=transpose(linsolve(transpose(J2-AA),transpose(Del2),'r',a));

V2 := matrix([[1/3, 1/2, 1/6, 1/6, 1/2, 4/3], [1/6, (-1)/16, 13/48, 7/48, 5/16, 2/3], [(-1)/6, (-5)/16, (-7)/48, (-13)/48, 1/16, (-2)/3], [(-1)/6, 1/16, (-13)/48, (-7)/48, (-5)/16, (-2)/3], [1/6, 5/16...

>                 #####################################################################################

>                    #################### !!!!!!!!!!!!!!!!!! START RUNS HERE !!!!!!!!!!!!! ######################################

>                 #####################################################################################

> print(matrans(R),matrans(B));

vector([4, 3, 1, 2]), vector([3, 4, 4, 3])

> "rank Delta"=rank(Delta);"rank A"=rank(A);"min-poly",factor(minpoly(Delta,x));unassign('phi','ee','e');print("pi",pi);

> d:=rank(Delta);n:=rowdim(A):zeta:=vector(n,0):ND:=nullspace(Delta);P:=concat(seq(ND[k],k=1..nops(ND)),seq(zeta,j=1..d)):nullspace(J-P);

d := 3

ND := {vector([1, 1, 1, 1])}

{vector([1, 1, 1, 1])}

        DETTA DEFINED HERE

> ee:=[1,1,1,1,1,1,1,1,1,1,1, 1,1, 1]:phi:=diag(seq(ee[q],q=1..n)):Detta:=evalm(phi&*Delta):R:=evalm(A+Detta):B:=evalm(A-Detta):matrans(R),matrans(B);factor(minpoly(R,s)),factor(minpoly(B,s));nullspace(transpose(Detta));minpoly(Detta,x),factor(minpoly(Detta,x));D2:=sympow(Detta,2):psi:=sympow(phi,2):"MORPHISM CHECK",iszero(evalm(psi&*Del2-D2));A2:=evalm(AA+D2);

vector([4, 3, 1, 2]), vector([3, 4, 4, 3])

(s-1)*(s+1)*(s^2+1), s*(s-1)*(s+1)

{vector([1, 1, 0, 0])}

1/4*x+1/4*x^2+x^4, 1/4*x*(2*x+1)*(2*x^2-x+1)

A2 := matrix([[0, 0, 0, 0, 0, 1], [0, 0, 1/2, 0, 0, 1/2], [0, 0, 0, 0, 1/2, 0], [0, 1/2, 0, 0, 0, 0], [0, 0, 0, 1/2, 0, 1/2], [1/2, 0, 0, 0, 0, 1/2]])

> map(abs,map(evalf,[eigenvalues(AA)]));map(abs,map(evalf,[eigenvalues(A2)]));

[0., 0., 0., .2500000000, .8090169942, .3090169942]

[1., .5000000000, .5000000000, .5000000000, .5000000000, .5000000000]

                PLOT OF det(I-Asup2-t times Del2)

> unassign('t'):dp:=simplify(det(J2-AA-t*D2));subs(t=0,dp),subs(t=1,dp);plot(dp,t=-0.1..1.1);

dp := -1/128*t^4-3/64*t^3-11/128*t^2-3/64*t+3/16

3/16, 0

[Plot]

> trace(D2);V2:=multiply(D2,inverse(J2-AA)):trace(V2)*(-1);

(-1)/4

(-1)/4

        FULLY SYMMETRIC ARE COMING HERE

> AS2:=sym2(A):RS2:=sym2(R):BS2:=sym2(B):as2:=evalm(1/2*(RS2+BS2)):

> print("A-sup-2",AA,"A",A,"A-SYM2",AS2);

> JX:=IdentityMatrix(NN+n):unassign('t'):factor(det(JX-t*AS2));

> "pi for ASYM2",nullspace(transpose(JX-AS2));

1/256*(t-1)*(t+2)^2*(t-4)^3

> print("RS2",RS2,"BS2",BS2,"as2 combines RS2 and BS2",as2);

> unassign('x'):psi2:=evalm((linsolve(transpose(JX-as2),vector(n+NN,0),'r',x)));us2:=evalm((linsolve(JX-as2,vector(n+NN,0),'r',x)));

psi2 := vector([x[1], 0, 0, 0, 0, 2*x[1], x[2], x[2], 2*x[2], 2*x[2]])

us2 := vector([x[2], x[1], 2*x[1]-x[2], 2*x[1]-x[2], x[1], x[2], 3*x[1]-2*x[2], 3*x[1]-2*x[2], 3*x[1]-2*x[2], 3*x[1]-2*x[2]])

        HERE THE ABSORBING STATE HAS BEEN EXTENDED TO n STATES

> abas2:=abel(as2);

abas2 := matrix([[1/3, 0, 0, 0, 0, 2/3, 0, 0, 0, 0], [2/9, 0, 0, 0, 0, 4/9, 1/18, 1/18, 1/9, 1/9], [1/9, 0, 0, 0, 0, 2/9, 1/9, 1/9, 2/9, 2/9], [1/9, 0, 0, 0, 0, 2/9, 1/9, 1/9, 2/9, 2/9], [2/9, 0, 0, 0...


                LEVEL TWO HERE

> R2:=sympow(R,2):B2:=sympow(B,2):

> N2:=binomial(n,2):J2:=evalm(IdentityMatrix(N2)):UJ:=matrix(n,n,1):

> nx:=evalm(linsolve(Omega,vector(n,0),'r',x));

nx := vector([-x[1]-2*x[2]-2*x[3], x[1], x[2], x[3]])

> unassign('x'):pi2:=evalm(1/x[1]*linsolve(J2-transpose(A2),vector(N2,0),'r',x)):u2:=evalm(1/x[1]*linsolve(J2-A2,vector(N2,0),'r',x)):

> print("pi2",pi2);print("u2",evalm(u2));uu:=vector(N2,1):

> pssi2:=subs({x[1]=1,x[2]=0},evalm(psi2));uss2:=subs({x[1]=2,x[2]=3},evalm(us2));

pssi2 := vector([1, 0, 0, 0, 0, 2, 0, 0, 0, 0])

uss2 := vector([3, 2, 1, 1, 2, 3, 0, 0, 0, 0])

> rank(abas2);

2