Revision history [back]
 | 1 | initial version |
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
| | 2 | No.2 Revision |
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}}
