Okay, there’s all sorts of stuff — NumPy, linear algebra, matrix arithmetic — I haven’t looked at. Consider, masks, boolean operations, fancy indexing, sorting, determinants, eigenvectors, inversion, etc. And, I am not going to do so at this time. But, I do want to take a quick, very quick, look at solving equations using NumPy.

Systems of Linear Equations

A linear system of equations is a collection of linear equations

$$\begin{align} a_{0,0}x_0 + a_{0,1}x_2 + \cdots + a_{0,n}x_n & = b_0 \\ a_{1,0}x_0 + a_{1,1}x_2 + \cdots + a_{1,n}x_n & = b_1 \\ & \vdots \\ a_{m,0}x_0 + a_{m,1}x_2 + \cdots + a_{m,n}x_n & = b_m \\ \end{align}$$

In matrix notation, a linear system is \(A \mathbf{x}= \mathbf{b}\) where

$$A = \begin{bmatrix} a_{0,0} & a_{0,1} & \cdots & a_{0,n} \\ a_{1,0} & a_{1,1} & \cdots & a_{1,n} \\ \vdots & & & \vdots \\ a_{m,0} & a_{m,1} & \cdots & a_{m,n} \\ \end{bmatrix} \ \ , \ \ \mathbf{x} = \begin{bmatrix} x_0 \\ x_1 \\ \vdots \\ x_n \end{bmatrix} \ \ , \ \ \mathbf{b} = \begin{bmatrix} b_0 \\ b_1 \\ \vdots \\ b_m \end{bmatrix}$$

Solving Linear Systems from Mathematical Python by Patrick Walls, an instructor in the Department of Mathematics at the University of British Columbia

Our First Problem: 2 Equations

Let’s assume we have a small bakery. Earlier this week we sold 31 loaves of rye and 22 loaves of pumpernickel for $167. The following day we sold 21 rye and 32 pumpernickel for $161. If the prices remained unchanged both days what was the price for a loaf of each type of bread. Let’s build a system of equations.

31r + 22p = 167
21r + 32p = 161

So our A matrix and b vector look like:

$$A = \begin{bmatrix} 31 & 22 \\ 21 & 32 \\ \end{bmatrix} \ \ , \ \ \mathbf{b} = \begin{bmatrix} 167 \\ 161 \end{bmatrix}$$

In [2]:
# okay let's build our matrices/vectors
A = np.array([[31, 22], [21, 32]])
b = np.array([167, 161])
print('A =', A, '\nb =', b)

A = [[31 22]
 [21 32]] 
b = [167 161]
In [3]:
x = np.linalg.solve(A,b)
print(x)

[3.4 2.8]

And, the rye was going for $3.40 a loaf and the pumpernickel for $2.80 a loaf. Let’s make sure.

In [5]:
print(31 * 3.4 + 22 * 2.8)
print(21 * 3.4 + 32 * 2.8)

167.0
161.0

Our Second Problem: 3 Equations

Let’s mess with a system of 3 linear equations.

Our friend Sharon sells hand made teddy bears/bunnys at craft fairs. She prices the bears according to size: small critters cost $15, medium size critters cost $25, and huggable crittes cost $55. (Well the materials used also affect the price, but for this exercise we can simplify things a wee bit.) She usually sells as many small critters as medium and large critters combined. She also sells twice as many medium softies as huggable ones.

This is the most expensive craft fair Sharon attends, by a fair margin. A booth at the craft fair costs $300.

If her sales go as usual, how many of each size bear/bunny must she sell to pay for the booth?

Let’s sort out our system of linear equations.

15s + 25m + 55h = 300
m = 2h
s = m + h

Which becomes:

15s + 25m + 60h = 300
0s + m - 2h = 0
s - m - h = 0

So our A matrix and b vector look like:

$$A = \begin{bmatrix} 15 & 25 & 55 \\ 1 & -1 & -1 \\ 0 & 1 & -2 \\ \end{bmatrix} \ \ , \ \ \mathbf{b} = \begin{bmatrix} 300 \\ 0 \\ 0 \end{bmatrix}$$

In [6]:
# let's build our arrays
A2 = np.array([[15, 25, 55], [1, -1, -1], [0, 1, -2]])
b2 = np.array([300, 0, 0])
x3 = np.linalg.solve(A2, b2)
print(x3)

[6. 4. 2.]

So, we need to sell 6 small, 4 medium and 2 huggable to cover the cost of the booth. Let’s make sure our result is valid.

In [7]:
s, m, h = x3
print(f"15s + 25m + 55h = {round(15 * s + 25 * m + 55 * h)}")
print(f"1st row of A dot result vector = {round(np.dot(A2[0, :], x3))}")
print(f"s - m - h = {round(s - m - h)}")
print(f"2nd row of A dot result vector = {round(np.dot(A2[1, :], x3))}")
print(f"m - 2h = {round(m - 2*h)}")
print(f"3rd row of A dot result vector = {round(np.dot(A2[2, :], x3))}")

15s + 25m + 55h = 300
1st row of A 'dot' result vector = 300
s - m - h = 0
2nd row of A 'dot' result vector = 0
m - 2h = 0
3rd row of A 'dot' result vector = 0

You did notice that the dot product of each row of coefficients with the solution vector will give us the value for the equation with that set of coefficients.

matrix vs ndarray

You may have someone tell you that you should be using numpy.matrix() to generate your matrices for linear algebra uses rather than numpy.array(). This, at one time, was probably true. But, according to the NumPy site, that is no longer the case.

It is strongly advised not to use the matrix subclass. As described below, it makes writing functions that deal consistently with matrices and regular arrays very difficult. Currently, they are mainly used for interacting with scipy.sparse. We hope to provide an alternative for this use, however, and eventually remove the matrix subclass.

numpy.matrix

You have been warned!

And, We’re Done, m’thinks

I believe we now have a decent introduction to NumPy.

Feel free to download my notebook covering the above and play around.

Next post I plan to start looking at pandas.

Resources