15

So I was given a homework assignment that requires solving the coefficients of cubic splines. Now I clearly understand how to do the math on paper as well as with MatLab, I want to solve the problem with Python. Given an equation Ax = b where I know the values of A and b, I want to be able to solve for x with Python and I am having trouble finding a good resource to do such a thing.

Ex.

A = |1 0 0|
    |1 4 1|
    |0 0 1|

x = Unknown 3x1 matrix

b = |0 |
    |24| 
    |0 |

Solve for x

Bach
  • 6,145
  • 7
  • 36
  • 61
Scalahansolo
  • 2,615
  • 6
  • 26
  • 43
  • 1
    @MattDMo: OP already tagged `numpy`. – Amadan Mar 04 '14 at 04:43
  • I have looked at NumPy a little but its a lot to get through, do you know which NumPy function(s), or at least which area of NumPy would handle this best? – Scalahansolo Mar 04 '14 at 04:44
  • 2
    If `Ax = B`, `x = (A^-1)B`. Take a look at `inv` and `dot` functions. – Amadan Mar 04 '14 at 04:45
  • @Amadan - true, didn't see that... – MattDMo Mar 04 '14 at 04:46
  • 2
    as a general reference, take a look at the [NumPy for Matlab Users page](http://wiki.scipy.org/NumPy_for_Matlab_Users) if you haven't come across it already. Scrolling down, there's a big list of linear algebra equivalents that may be helpful, as well as a variety of other comparisons to help you from Matlab to the exciting world of Python :) – MattDMo Mar 04 '14 at 04:48
  • @Amadan, that was a big help. Thanks! Exactly what I was looking for. – Scalahansolo Mar 04 '14 at 04:53
  • 1
    Extra marks for disclosing this is a homework assignment. That lets us treat it accordingly, giving you the benefit of the work you have to do to understand the solutions. – holdenweb Mar 04 '14 at 05:32
  • Solving linear equations by calculating the inverse like @Amadan suggested is a nice way to write the solution mathematically but it should be avoided for real code as there are far better algorithms for it. – Dirk Mar 04 '14 at 08:24
  • @holdenweb Only extra marks for disclosing that it's a homework assignment if it actually is. ;-) (I came here to learn Python as a teacher.) – jvriesem Mar 18 '21 at 23:05

3 Answers3

20

In a general case, use solve:

>>> import numpy as np
>>> from scipy.linalg import solve
>>> 
>>> A = np.random.random((3, 3))
>>> b = np.random.random(3)
>>> 
>>> x = solve(A, b)
>>> x
array([ 0.98323512,  0.0205734 ,  0.06424613])
>>> 
>>> np.dot(A, x) - b
array([ 0.,  0.,  0.])

If your problem is banded (which cubic splines it often is), then there's http://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.solve_banded.html

To comment on some of the comments to the question: better not use inv for solving linear systems. numpy.lstsq is a bit different, it's more useful for fitting.

As this is homework, you're really better off at least reading up on ways of solving tridiagonal linear systems.

ev-br
  • 24,968
  • 9
  • 65
  • 78
8

Numpy is the main package for scientific computing in Python. If you are windows user then download it here: http://www.lfd.uci.edu/~gohlke/pythonlibs/#numpy else follow these instructions: http://www.scipy.org/install.html.

import numpy
A = [[1,0,0],[1,4,1],[0,0,1]]
b = [0,24,0]
x = numpy.linalg.lstsq(A,b)
BushMinusZero
  • 1,202
  • 16
  • 21
2

In addition to the code of Zhenya, you might also find it intuitive to use the np.dot function:

import numpy as np
A = [[1,0,0],
    [1,1,1],
    [6,7,0]]
b = [0,24,0]
# Now simply solve for x
x = np.dot(np.linalg.inv(A), b) 
#np.linalg.inv(A)  is simply the inverse of A, np.dot is the dot product
print x

Out[27]: array([  0.,   0.,  24.])
moldovean
  • 3,132
  • 33
  • 36