Sage wont solve simultaneous equation

edit

It is hard to understand what happens behind, but the following works:

vars = var('a b c d e f g h u')

X = matrix(3, 3, [[0, 1, 0], [0, 0, 1], [1, 0, 0]])
P = matrix(3, 3, [[a, b, c], [d, e, f], [g, h, u]])

Pdagger = P.transpose()
Xdagger = X.transpose()    # not needed

Q = X * Pdagger
B = Q.solve_left(X)

R     = range(0,3)
eqns  = [ B[i][j] for i in R for j in R ]
xlist = [ X[i][j] for i in R for j in R ]
final = [ B[i][j] == X[i][j] for i in R for j in R ]

sols = solve(final, vars) 
# show(sol)
for sol in sols:
    print sol

and i have some comments:

• when programming, please use a decent spacing, else nobody will see an error when an obvious one exists. "pep8".
• if i is a variable... ok, no problem we overwrite that $\sqrt {-1}$, then please do not redefine it in a loop as a number. (This would have been a second alarm while reusing i...) Nobody does this, so this atypical error cannot be seen at a first glance. In our case this brings the error as follows. After the last i-loop this "global variable" has the value $2$, and we try to solve a system of equations, where one of the equation variables is the constant $2$. The error message is also misleading.
• give a minimal code showing the error, only needed lines of code.
• use list comprehension.
• do not call a list eqns, when it does not collect equations, but rather entries of a matrix. The solve of these "equations" will be hardly detected with bare eyes as an error.
• try to give pointed prints. The show command delivers a mess, is not a good control of the operations.
• the two definitions of eqns and xlist are not needed. The wording for xlist is also not proper. Something like X_entries would reflect the provenience, but we do not need this new variable, the matrix X already has the information.

Please understand that these "critics" are just kind advices, at any rate not offensive, in my case, knowing these as i started programming in python would have been a real benefit.

Let's be lazy (and avoid identifier clashes)

sage: V=var("v", n=9)
sage: P=matrix(V,nrows=3)
sage: X=matrix(SR,3,3,[[0,1,0],[0,0,1],[1,0,0]])
sage: Pdagger=P.transpose()
sage: Xdagger=X.transpose()
sage: Q=X*Pdagger
sage: B=Q.solve_left(X)

Let's be lazy again : don't build a list of equations representing B==X, but a list representing the elements of B-X

sage: E=[B[r][c]-X[r][c] for c in range(3) for r in range(3)]

Look (lazily) for solutions

sage: S=solve(E,V)

A loook at the solutions show that they are rational expressions. Are they acceptable ? This entails that their denominators must not be null. But :

sage: [[e.denominator().subs(s) for s in S] for e in E]
[[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1],
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1]]

Check the only acceptable one :

sage: [e.subs(S[16]) for e in E]
[0, 0, 0, 0, 0, 0, 0, 0, 0]

And, indeed :

sage: matrix([v.subs(S[16]) for v in V],nrows=3)==X
True

Note that Sympy does not seem to be able to solve this system :

sage: sympy.solve(map(lambda e:e._sympy_(), E),map(lambda v:v._sympy_(), V))
[]

On the other hand, Mathematica is able to solve the system and does not gives "spurious" solutions. But expressing the exuations is harder in Mathematica than in Sage : that's why I keep the first part in Sage :

sage: V=var("v",n=9)
sage: P=matrix(SR,V,nrows=3);P
[v0 v1 v2]
[v3 v4 v5]
[v6 v7 v8]
sage: X=matrix(SR,[[0,1,0],[0,0,1],[1,0,0]]);X
[0 1 0]
[0 0 1]
[1 0 0]
sage: Q=X*(P.transpose())
sage: B=Q.solve_left(X)
sage: E=[B[l][c]==X[l][c] for c in range(3) for l in range(3)];E
[(v2*v7/v1 - v8)*(v0/v1 + v2*(v3 - v0*v4/v1)/(v1*(v2*v4/v1 - v5)))/((v2*v4/v1 - v5)*(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) + v2/(v1*(v2*v4/v1 - v5)) == 0,
 -(v0/v1 + v2*(v3 - v0*v4/v1)/(v1*(v2*v4/v1 - v5)))/(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1) == 0,
 -(v0/v1 + v2*(v3 - v0*v4/v1)/(v1*(v2*v4/v1 - v5)))*(v4*(v2*v7/v1 - v8)/(v1*(v2*v4/v1 - v5)) - v7/v1)/(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1) + 1/v1 - v2*v4/(v1^2*(v2*v4/v1 - v5)) == 1,
 -1/(v2*v4/v1 - v5) - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/((v2*v4/v1 - v5)^2*(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) == 1,
 (v3 - v0*v4/v1)/((v2*v4/v1 - v5)*(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) == 0,
 (v3 - v0*v4/v1)*(v4*(v2*v7/v1 - v8)/(v1*(v2*v4/v1 - v5)) - v7/v1)/((v2*v4/v1 - v5)*(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) + v4/(v1*(v2*v4/v1 - v5)) == 0,
 -(v2*v7/v1 - v8)/((v2*v4/v1 - v5)*(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) == 0,
 (1/(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1)) == 1,
 (v4*(v2*v7/v1 - v8)/(v1*(v2*v4/v1 - v5)) - v7/v1)/(v6 - (v3 - v0*v4/v1)*(v2*v7/v1 - v8)/(v2*v4/v1 - v5) - v0*v7/v1) == 0]
sage: mathematica.Solve(mathematica.And(*E), V)
{{v0 -> 0, v1 -> 1, v2 -> 0, v3 -> 0, v4 -> 0, v5 -> 1, v6 -> 1, v7 -> 0, 
  v8 -> 0}}