#THIS PROGRAM RUNS IN SAGEMATH 8.1

#MARKOV PARTITION FOR (2,1,1,1) WITH A refinement of SQUARE S OF SIDE l WITH VERTEX ON 0

M=matrix(2,2,[2,1,1,1])
S=M.eigenvectors_right()
evs=S[0][0]
evu=S[1][0]
EVs=vector(S[0][1][0])
EVu=vector(S[1][1][0])
Mi=M.inverse()
EVu=EVu/EVu.norm()
EVs=EVs/EVs.norm()

P=arrow((0,0), evu*EVu, legend_label='unstable direction', legend_color='red',color='red')+arrow( evs*EVs,(0,0), legend_label='stable direction', legend_color='green',color='green')


l=1/sqrt(1+(1-evs)^2)
A1=EVu*l
B1=-EVs*l
C1=A1+B1
D1=vector([0,0])

A2=D1+vector([1,0])
B2=C1+vector([0,-1])
C2=A1
D2=B1+vector([1,-1])


def Itrans(M,v):
	return [M[0]+v,M[1]+v,M[2]+v,M[3]+v]

R1=polygon([D1, A1, C1,B1], fill = True, color='blue',alpha=0.3)
R1L=polygon(Itrans([D1, A1, C1,B1],vector([1,0])), fill = True, color='blue',alpha=0.3)
R1B=polygon(Itrans([D1, A1, C1,B1],vector([0,-1])), fill = True, color='blue',alpha=0.3)
R1LB=polygon(Itrans([D1, A1, C1,B1],vector([1,-1])), fill = True, color='blue',alpha=0.3)

R2=polygon([D2, A2, C2,B2], fill = True, color='orange',alpha=0.5)
R2L=polygon(Itrans([D2, A2, C2,B2],vector([1,0])), fill = True, color='orange',alpha=0.5)
R2B=polygon(Itrans([D2, A2, C2,B2],vector([0,-1])), fill = True, color='orange',alpha=0.5)
R2LB=polygon(Itrans([D2, A2, C2,B2],vector([1,-1])), fill = True, color='orange',alpha=0.5)


R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB

MR1=polygon([M*D1, M*A1, M*C1,M*B1], fill = True, color='red',alpha=0.2)
MR1L=polygon(Itrans([M*D1, M*A1, M*C1,M*B1],vector([1,0])), fill = True, color='red',alpha=0.2)
MR1B=polygon(Itrans([M*D1, M*A1, M*C1,M*B1],vector([0,-1])), fill = True, color='red',alpha=0.2)
MR1LB=polygon(Itrans([M*D1, M*A1, M*C1,M*B1],vector([1,-1])), fill = True, color='red',alpha=0.2)

R1+R1L+R1B+R1LB+MR1+MR1L+MR1B+MR1LB

MiR1=polygon([M^(-1)*D1, M^(-1)*A1, M^(-1)*C1,M^(-1)*B1], fill = True, color='red',alpha=0.2)
MiR1L=polygon(Itrans([M^(-1)*D1, M^(-1)*A1, M^(-1)*C1,M^(-1)*B1],vector([1,0])), fill = True, color='red',alpha=0.2)
MiR1B=polygon(Itrans([M^(-1)*D1, M^(-1)*A1, M^(-1)*C1,M^(-1)*B1],vector([0,-1])), fill = True, color='red',alpha=0.2)
MiR1LB=polygon(Itrans([M^(-1)*D1, M^(-1)*A1, M^(-1)*C1,M^(-1)*B1],vector([1,-1])), fill = True, color='red',alpha=0.2)

R1+R1L+R1B+R1LB+MR1+MR1L+MR1B+MR1LB+MiR1LB+MiR1+MiR1L+MiR1LB


################################################

MR2=polygon([M*D2, M*A2, M*C2,M*B2], fill = True, color='red',alpha=0.4)
MR2L=polygon(Itrans([M*D2, M*A2, M*C2,M*B2],vector([1,0])), fill = True, color='red',alpha=0.4)
MR2B=polygon(Itrans([M*D2, M*A2, M*C2,M*B2],vector([0,-1])), fill = True, color='red',alpha=0.4)
MR2LB=polygon(Itrans([M*D2, M*A2, M*C2,M*B2],vector([1,-1])), fill = True, color='red',alpha=0.4)

R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB+MR2+MR2L+MR2B+MR2LB

MiR2=polygon([M^(-1)*D2, M^(-1)*A2, M^(-1)*C2,M^(-1)*B2], fill = True, color='red',alpha=0.4)
MiR2L=polygon(Itrans([M^(-1)*D2, M^(-1)*A2, M^(-1)*C2,M^(-1)*B2],vector([1,0])), fill = True, color='red',alpha=0.4)
MiR2B=polygon(Itrans([M^(-1)*D2, M^(-1)*A2, M^(-1)*C2,M^(-1)*B2],vector([0,-1])), fill = True, color='red',alpha=0.4)
MiR2LB=polygon(Itrans([M^(-1)*D2, M^(-1)*A2, M^(-1)*C2,M^(-1)*B2],vector([1,-1])), fill = True, color='red',alpha=0.4)

R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB+MiR2+MiR2L+MiR2B+MiR2LB


MMiR2=polygon([M^(-2)*D2, M^(-2)*A2, M^(-2)*C2,M^(-2)*B2], fill = True, color='red',alpha=0.4)
MMiR2L=polygon(Itrans([M^(-2)*D2, M^(-2)*A2, M^(-2)*C2,M^(-2)*B2],vector([1,0])), fill = True, color='red',alpha=0.4)
MMiR2B=polygon(Itrans([M^(-2)*D2, M^(-2)*A2, M^(-2)*C2,M^(-2)*B2],vector([0,-1])), fill = True, color='red',alpha=0.4)
MMiR2LB=polygon(Itrans([M^(-2)*D2, M^(-2)*A2, M^(-2)*C2,M^(-2)*B2],vector([1,-1])), fill = True, color='red',alpha=0.4)

R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB+MMiR2+MMiR2L+MMiR2B+MMiR2LB

MMMiR2=polygon([M^(-3)*D2, M^(-3)*A2, M^(-3)*C2,M^(-3)*B2], fill = True, color='red',alpha=0.4)
MMMiR2L=polygon(Itrans([M^(-3)*D2, M^(-3)*A2, M^(-3)*C2,M^(-3)*B2],vector([1,0])), fill = True, color='red',alpha=0.4)
MMMiR2B=polygon(Itrans([M^(-3)*D2, M^(-3)*A2, M^(-3)*C2,M^(-3)*B2],vector([0,-1])), fill = True, color='red',alpha=0.4)
MMMiR2LB=polygon(Itrans([M^(-3)*D2, M^(-3)*A2, M^(-3)*C2,M^(-3)*B2],vector([1,-1])), fill = True, color='red',alpha=0.4)

R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB+MMMiR2+MMMiR2L+MMMiR2B+MMMiR2LB


MMMMiR2=polygon([M^(-4)*D2, M^(-4)*A2, M^(-4)*C2,M^(-4)*B2], fill = True, color='red',alpha=0.4)
MMMMiR2L=polygon(Itrans([M^(-4)*D2, M^(-4)*A2, M^(-4)*C2,M^(-4)*B2],vector([1,0])), fill = True, color='red',alpha=0.4)
MMMMiR2B=polygon(Itrans([M^(-4)*D2, M^(-4)*A2, M^(-4)*C2,M^(-4)*B2],vector([0,-1])), fill = True, color='red',alpha=0.4)
MMMMiR2LB=polygon(Itrans([M^(-4)*D2, M^(-4)*A2, M^(-4)*C2,M^(-4)*B2],vector([1,-1])), fill = True, color='red',alpha=0.4)

R1+R1L+R1B+R1LB+R2+R2L+R2B+R2LB+MMMMiR2+MMMMiR2L+MMMMiR2B+MMMMiR2LB