Solving linear systems of equations is straightforward using the scipy command linalg.solve. This command expects an input matrix and a right-hand-side vector. The solution vector is then computed. An option for entering a symmetrix matrix is offered which can speed up the processing when applicable. As an example, suppose it is desired to solve the following simultaneous equations:
\begin{eqnarray*} x+3y+5z & = & 10\\ 2x+5y+z & = & 8\\ 2x+3y+8z & = & 3\end{eqnarray*}
We could find the solution vector using a matrix inverse:
\[ \left[\begin{array}{c} x\\ y\\ z\end{array}\right]=\left[\begin{array}{ccc} 1 & 3 & 5\\ 2 & 5 & 1\\ 2 & 3 & 8\end{array}\right]^{-1}\left[\begin{array}{c} 10\\ 8\\ 3\end{array}\right]=\frac{1}{25}\left[\begin{array}{c} -232\\ 129\\ 19\end{array}\right]=\left[\begin{array}{c} -9.28\\ 5.16\\ 0.76\end{array}\right].\]
However, it is better to use the linalg.solve command which can be faster and more numerically stable. In this case it however gives the same answer as shown in the following example:
>>> import numpy as np
>>> from scipy import linalg
>>> A = np.array([[1,2],[3,4]])
>>> A
array([[1, 2],
[3, 4]])
>>> b = np.array([[5],[6]])
>>> b
array([[5],
[6]])
>>> linalg.inv(A).dot(b) #slow
array([[-4. ],
[ 4.5]]
>>> A.dot(linalg.inv(A).dot(b))-b #check
array([[ 8.88178420e-16],
[ 2.66453526e-15]])
>>> np.linalg.solve(A,b) #fast
array([[-4. ],
[ 4.5]])
>>> A.dot(np.linalg.solve(A,b))-b #check
array([[ 0.],
[ 0.]])

Post a Comment