Processing math: 100%

First time here? Check out the FAQ!

Ask Your Question
1

How to use `solve()` for a system of equations using matrices?

asked 5 years ago

torito_verdejo gravatar image

I know I can solve a system of equations by inputing independently each equation in a same solve() expression using the syntax solve([[exp1],[exp2], ... [expn]], x1,x2, ... xn), but what should I do if, having defined a matrix A and two columns vectors x and y, I want to express the system of equations as A*x == y?

Each instance of A[i] * x returns a list, that is (of course) not treated as a valid argument for an equation in solve(). How could I "tranform" that list in order to use it with solve()?

Preview: (hide)

Comments

thank you for asking that. not knowing how to answer this question,(I did not dare to ask the question) I did like that but I'm sure it must be easier way ! see this post and the cocalc link in it

ortollj gravatar imageortollj ( 5 years ago )

2 Answers

Sort by » oldest newest most voted
3

answered 5 years ago

Emmanuel Charpentier gravatar image

updated 5 years ago

Try this :

sage: A=matrix([[var("a_{}_{}".format(u,v), latex_name="a_{{{},{}}}".format(u,v)) for v in (1..2)] for u in (1..2)]);A
[a_1_1 a_1_2]
[a_2_1 a_2_2]
sage: X=vector([var("x_{}".format(u)) for u in (1..2)]);X
(x_1, x_2)
sage: Y=vector([var("y_{}".format(u)) for u in (1..2)]);Y
(y_1, y_2)
sage: Sol=solve(map(lambda u,v:u==v, A*X, Y),[u for u in X]);Sol
[[x_1 == -(a_2_2*y_1 - a_1_2*y_2)/(a_1_2*a_2_1 - a_1_1*a_2_2), x_2 == (a_2_1*y_1 - a_1_1*y_2)/(a_1_2*a_2_1 - a_1_1*a_2_2)]]

i. e. [[x1=a2,2y1a1,2y2a1,2a2,1a1,1a2,2,x2=a2,1y1a1,1y2a1,2a2,1a1,1a2,2]].

which looks entirely reasonable...

HTH,

Preview: (hide)
link

Comments

Actually I don't understand how this actually works. A*x returns a list, as map() does, but A*x won't work while you're construction does. Why?

torito_verdejo gravatar imagetorito_verdejo ( 5 years ago )
1

A*x returns a list

Nope :

sage: (A*X).parent()
Vector space of dimension 2 over Symbolic Ring

But A*x is an iterable, as well as Y, and map iterates through those, producing at each step one equation, and returns the list of those equations. Similarly, I do not use X as the second argument of solve, but the list of its components.

?map should enlighten you :-).

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 5 years ago )
1

POST-SCRIPTUM : It turns out that vectors have a list` method, that you can use to lighten the notations :

sage: solve((A*X-Y).list(),X.list())
[[x_1 == -(a_2_2*y_1 - a_1_2*y_2)/(a_1_2*a_2_1 - a_1_1*a_2_2), x_2 == (a_2_1*y_1 - a_1_1*y_2)/(a_1_2*a_2_1 - a_1_1*a_2_2)]]

Lighter and cleaner, indeed...

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 5 years ago )
1

answered 5 years ago

Juanjo gravatar image

updated 5 years ago

I suppose you are aware that, in fact, you don't need to use solve: once defined A and Y, the solution of the linear system AX=Y is given by A\Y. For example,

sage: A = matrix([[1,2,3], [4,5,6], [-7,8,9]])
sage: Y = vector([3,9,1])
sage: print "Solution:", A\Y
Solution: (1, 1, 0)

Anyway, if you really need to use solve and write the system in the form AX=Y, you can continue the above example as follows:

sage: X = vector([var("x_{}".format(i)) for i in [1..3]])
sage: system = [A[i]*X==Y[i] for i in range(len(Y))]
sage: solve(system, *X)
[[x_1 == 1, x_2 == 1, x_3 == 0]]

Post-scriptum. The following function covers the cases considered in the comments below:

def my_solve(A, Y):
    m, n = A.dimensions()
    p, q = Y.dimensions()
    if m!=p:
        raise RuntimeError("The matrices have different numbers of rows")
    X = vector([var("x_{}".format(i)) for i in [1..n]])
    sols = []
    for j in range(q):
        system = [A[i]*X==Y[i,j] for i in range(m)]
        sols += solve(system, *X)
    return sols

Of course, A and Y should be two matrices with the same number of rows. Let us check the examples:

sage: A = matrix(QQ, [[1, 2, 2, 2], [2, 4, 6, 8], [3, 6, 8, 10], [0, 0, 0, 0]])
sage: Y = matrix(QQ, 4, 1, [0, 0, 0, 0])
sage: my_solve(A, Y)
[[x_1 == 2*r1 - 2*r2, x_2 == r2, x_3 == -2*r1, x_4 == r1]]

Please, note the way Y is given.

sage: A = matrix(QQ, [[2, 4, -3], [-3, -2, 4], [2, -5, -2]])
sage: Y = matrix(QQ, [[-5, -2], [3, 1], [0, 3]])
sage: my_solve(A, Y)
[[x_1 == 51, x_2 == 4, x_3 == 41], [x_1 == -5, x_2 == -1, x_3 == -4]]

Update. The code had an error in the way the number of unknowns was determined. I have modified the code to correct this bug. Let us check it with this new example:

sage: A = matrix([[1, 2, 3], [5, 7, 11]])
sage: Y = matrix([[13 ,17], [19, 23]])
sage: my_solve(A, Y)
[[x_1 == -1/3*r11 - 53/3, x_2 == -4/3*r11 + 46/3, x_3 == r11],
 [x_1 == -1/3*r12 - 73/3, x_2 == -4/3*r12 + 62/3, x_3 == r12]]

You can see these examples in this Sagemath Cell.

Preview: (hide)
link

Comments

Juanjo, the code i put herein cocalc works for this 2 types of Matrix: not your solution ?

N=matrix(QQ,[[ 1 , 2 ,  2 , 2 ],[ 2 , 4 , 6 , 8 ],[ 3, 6 , 8 , 10 ],[ 0, 0 , 0 , 0 ]])
M=matrix(QQ,[0,0,0,0]).transpose()

or

N=matrix(QQ,[[ 2  2 -2],[ 0  2 -4],[-5  4 -2]]

M=matrix(QQ,[[-4 -5],[-5 -4],[ 1  0]]

Y=\left(\begin{array}{rr}
-4 & 1 \\
0 & 2 \\
2 & 2
\end{array}\right)

X=\left(\begin{array}{rrr}
2 & -2 & 0 \\
-1 & 4 & -1 \\
-5 & 4 & 4
\end{array}\right)


or 

Y=\left(\begin{array}{r}
0 \\
0 \\
0 \\
0
\end{array}\right)

X=\left(\begin{array}{rrrr}
1 & 2 & 2 & 2 \\
2 & 4 & 6 & 8 \\
3 & 6 & 8 & 10 \\
0 & 0 & 0 & 0
\end{array}\right)
ortollj gravatar imageortollj ( 5 years ago )

what is headers for Latex code in this forum ? in the post above A=N and Y=M (now I cant edit my comment !!, too late)

ortollj gravatar imageortollj ( 5 years ago )

Emmanuel this does not give a solution too.

nEC=4
A=matrix(QQ,[[ 1 , 2 ,  2 , 2 ],[ 2 , 4 , 6 , 8 ],[ 3, 6 , 8 , 10 ],[ 0, 0 , 0 , 0 ]])
Y=matrix(QQ,[0,0,0,0]).transpose()
#A=matrix([[var("a_{}_{}".format(u,v), latex_name="a_{{{},{}}}".format(u,v)) for v in (1..nEC)] for u in (1..nEC)])
X=vector([var("x_{}".format(u)) for u in (1..nEC)])
#Y=vector([var("y_{}".format(u)) for u in (1..nEC)])
show(" A : ",A," X : ",X," Y : ",Y)
Sol=solve(map(lambda u,v:u==v, A*X, Y),[u for u in X]);Sol
ortollj gravatar imageortollj ( 5 years ago )

Emmanuel, your code does not work either with these two matrices:

A=matrix(QQ,[[ 2 , 4, -3],[-3 ,-2 , 4],[ 2 ,-5 ,-2]])
Y=matrix(QQ,[[-5, -2],[ 3 , 1],[ 0 , 3]])
ortollj gravatar imageortollj ( 5 years ago )

I forgot to precise, that my code is OK with these matrices

ortollj gravatar imageortollj ( 5 years ago )

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

2 followers

Stats

Asked: 5 years ago

Seen: 4,682 times

Last updated: Jun 20 '19