Notebook 1 – Math 2121, Fall 2020

This is a Pluto notebook, using the Julia programming language.

Julia is at the cutting edge of programming languages for scientific computing.

It combines many of the best features of Matlab and Python but is also very fast. The language has a lot of built-in functions that will be useful for demonstrating concepts in Math 2121.

xxxxxxxxxx
md"# Notebook 1 -- Math 2121, Fall 2020

1.9 ms

Installing Pluto (optional)

Pluto is an interactive notebook for running Julia code.

You can access this notebook as a static HTML page on our course website. But if you install Julia and Pluto on your home computer, then you can also modify and run the code here yourself.

I can't improve upon the instructions at this link for installing both:

https://mitmath.github.io/18S191/Fall20/installation/

xxxxxxxxxx
md"## Installing Pluto (optional)

2.4 ms

Running this Pluto notebook

Once you have Pluto up and running, you can access the notebook you we are currently viewing by entering this link (right click -> Copy Link) in the Open from file menu in Pluto.

xxxxxxxxxx
md"## Running *this* Pluto notebook

17 μs

Working with matrices

The syntax for working with matrices in Julia is much like Matlab.

xxxxxxxxxx
md"## Working with matrices

3.8 μs

3×5 Array{Int64,2}:
 1  2  3   0    6
 4  1  6  -1    8
 7  8  1   2  200

xxxxxxxxxx
 
# constructs a 3-by-5 matrix called Y
Y = [1 2 3 0 6; 4 1 6 -1 8; 7 8 1 2 200]

26.9 ms

-1

xxxxxxxxxx
 
# access entries in a matrix
Y[2, 4]

5.2 ms

Make it interactive

row_index = 1

col_index = 1

xxxxxxxxxx
begin

117 ms

3×5 Array{Int64,2}:
 1  2  3   0    6
 4  1  6  -1    8
 7  8  1   2  200

xxxxxxxxxx
 
# here is Y again
Y

585 ns

xxxxxxxxxx
 
Y[row_index, col_index]

3.6 μs

xxxxxxxxxx
 
# we can also set entries in a matrix
Y[2, 4] = 60;

3.8 ms

3×5 Array{Int64,2}:
 1  2  3   0    6
 4  1  6  60    8
 7  8  1   2  200

xxxxxxxxxx
 
# this changes Y to be this matrix
Y

916 ns

2×2 Array{Int64,2}:
 3   0
 6  60

xxxxxxxxxx
 
# create a submatrix, the slice notation is inclusive
Y[1:2,3:4]

56.4 ms

1×2 Array{Int64,2}:
 3  0

xxxxxxxxxx
 
Y[1:1,3:4]

5.2 μs

3×1 Array{Int64,2}:
  0
 60
  2

xxxxxxxxxx
 
Y[1:3,4:4]

7.7 μs

xxxxxxxxxx
 
# we can get the dimensions of a matrix this way
(m, n) = size(Y)

15.3 μs

3×4 Array{Float64,2}:
 1.1                2.2               3.0   0.0
 3.141592653589793  2.23606797749979  6.0  -1.0
 7.0                8.0               1.0   2.0

xxxxxxxxxx
 
# we can also create matrices with non-integer entries
Z = [1.1 2.2 3 0; pi sqrt(5) 6 -1; 7 8 1 2]

77.2 ms

Later in the course we will see many other operations involving matrices. Most of these will also be easily available in the Julia programming language.

xxxxxxxxxx
md"Later in the course we will see many other operations involving matrices. Most of these will also be easily available in the Julia programming language."

2.3 μs

Augmented matrices

We can write our own functions to work with matrices.

(This requires some programming background, however.)

xxxxxxxxxx
md"## Augmented matrices

6.3 μs

The next cell implements a function called print_linear_system which takes a single input parameter. You can probably guess what it does.

(To see the code, click the Show/Hide icon on the left.)

xxxxxxxxxx
md"The next cell implements a function called `print_linear_system` which takes a single input parameter. You can probably guess what it does. 

3.9 μs

print_linear_system (generic function with 1 method)

xxxxxxxxxx
begin

68 μs


 x_1 + 2 x_2 + 3 x_3 = 6

 4 x_1 + x_2 + 6 x_3 + 60 x_4 = 8

 7 x_1 + 8 x_2 + x_3 + 2 x_4 = 200

xxxxxxxxxx
 
# prints something that looks like the linear system whose augmented matrix is Y
print_linear_system(Y)

100 ms

xxxxxxxxxx
 
# what happens if Y has only one column?
print_linear_system(Y[1:3,4:4])

22.1 μs


 4 x_1 + x_2 + 6 x_3 = 60

 7 x_1 + 8 x_2 + x_3 = 2

xxxxxxxxxx
 
# how about two rows
print_linear_system(Y[2:3,1:4])

22.9 μs

Here's another function (code hidden) that returns true or false depending on whether a given list s is a solution to the linear system with augemented matrix A

xxxxxxxxxx
md"Here's another function (code hidden) that returns **true** or **false** depending on whether a given list **s** is a solution to the linear system with augemented matrix **A**"

13.4 μs

is_solution (generic function with 2 methods)

xxxxxxxxxx
function is_solution(A, s, tolerance=10e-16)

37.2 μs

3×4 Array{Int64,2}:
 1  -2   1   0
 0   2  -8   8
 5   0  -5  10

xxxxxxxxxx
 
# we saw this matrix in lecture 
T = [1 -2 1 0; 0 2 -8 8; 5 0 -5 10]

9.4 μs


 x_1 - 2 x_2 + x_3 = 0

 2 x_2 - 8 x_3 = 8

 5 x_1 - 5 x_3 = 10

xxxxxxxxxx
 
print_linear_system(T)

18.7 μs

3×4 Array{Int64,2}:
 1  -2   1   0
 0   1  -4   4
 0   0   1  -1

xxxxxxxxxx
 
# T is row equivalent to this triangular matrix
U = [1 -2 1 0; 0 1 -4 4; 0 0 1 -1]

6.3 μs


 x_1 - 2 x_2 + x_3 = 0

 x_2 - 4 x_3 = 4

 x_3 = -1

xxxxxxxxxx
 
print_linear_system(U)

20.6 μs

true

xxxxxxxxxx
 
is_solution(U, [1; 0; -1])

224 ms

true

xxxxxxxxxx
 
is_solution(T, [1; 0; -1])

6.9 μs

Row operations

A handy feature of these notebooks is that we can very easily explore different ranges of inputs and parameters. Let's try this out with our elementary row operations.

xxxxxxxxxx
md"## Row operations

5.9 μs

Set some parameters

a = b = c = p =

d = e = f = q =

g = h = i = r =

xxxxxxxxxx
begin

11.5 ms

3×4 Array{Int64,2}:
 1  2   3  -1
 0  0   0  -2
 3  8  12   0

xxxxxxxxxx
 
# now define a matrix with these parameters
A = [a b c p; d e f q; g h i r]

11.3 μs


 x_1 + 2 x_2 + 3 x_3 = -1

 0 = -2

 3 x_1 + 8 x_2 + 12 x_3 = 0

xxxxxxxxxx
 
# here's the linear system
print_linear_system(A)

18.7 μs

Below are functions implementing our three row operations. Their syntax is:

rowop_replace(A, i, j, v): adds v times row i in matrix A to row j
rowop_scale(A, i, v): multiplies row i in matrix A by nonzero constant v
rowop_swap(A, i, j): swaps rows i and j in matrix A

Each function returns a new matrix without modifying the matrix A.

xxxxxxxxxx
md"Below are functions implementing our three row operations. Their syntax is:

6.2 μs

rowop_replace (generic function with 1 method)

xxxxxxxxxx
 
function rowop_replace(A, source_row, target_row, scalar_factor)
    @assert source_row != target_row
    A = copy(A)
    (m, n) = size(A)
    for j=1:n
        A[target_row,j] += A[source_row,j] * scalar_factor
    end
    A
end

57.4 μs

rowop_scale (generic function with 1 method)

xxxxxxxxxx
 
function rowop_scale(A, target_row, scalar_factor)
    @assert scalar_factor != 0
    A = copy(A)
    (m, n) = size(A)
    for j=1:n
        A[target_row,j] *= scalar_factor
    end
    A
end

49.4 μs

rowop_swap (generic function with 1 method)

xxxxxxxxxx
 
function rowop_swap(A, source_row, target_row)
    A = copy(A)
    (m, n) = size(A)
    for j=1:n
        A[target_row,j], A[source_row,j] = A[source_row,j], A[target_row,j]
    end
    A
end

35.7 μs

Another parameter

v = 1

xxxxxxxxxx
 
begin
    md"""**Another parameter**
    
    v = $(@bind v Slider(-30:30, default=1, show_value=true))"""
end

77.8 μs

3×4 Array{Int64,2}:
 1  2   3  -1
 0  0   0  -2
 3  8  12   0

xxxxxxxxxx
 
# our original matrix
A

1.4 μs

3×4 Array{Int64,2}:
 1   2   3  -1
 0   0   0  -2
 4  10  15  -1

xxxxxxxxxx
 
# add v * row 1 to row 3
rowop_replace(A, 1, 3, v)

20 ms

true

xxxxxxxxxx
 
# this row operation can be reversed
rowop_replace(rowop_replace(A, 1, 3, v), 1, 3, -v) == A

8.8 ms

3×4 Array{Int64,2}:
 1  2   3  -1
 0  0   0  -2
 3  8  12   0

xxxxxxxxxx
 
rowop_scale(A, 3, v)

16.8 ms

true

xxxxxxxxxx
 
# also can be reversed
rowop_scale(rowop_scale(A, 3, v), 3, 1/v) == A

18 ms

3×4 Array{Int64,2}:
 1  2   3  -1
 3  8  12   0
 0  0   0  -2

xxxxxxxxxx
 
rowop_swap(A, 2, 3)

18.4 ms

true

xxxxxxxxxx
 
# doing this one twice undoes everything
rowop_swap(rowop_swap(A, 2, 3), 2, 3) == A

8.5 μs

Solving a linear system by row operations

By applying a sequence of row operations, we can transform the augmented matrix of a linear system to a simpler matrix, whose associated system has the same solutions as the one we started with.

xxxxxxxxxx
md"## Solving a linear system by row operations

6.7 μs

3×4 Array{Int64,2}:
 1  2   3  -1
 0  9   0  -2
 3  8  12   0

xxxxxxxxxx
 
M = [1 2 3 -1; 0 9 0 -2; 3 8 12 0]

9.6 μs


 x_1 + 2 x_2 + 3 x_3 = -1

 9 x_2 = -2

 3 x_1 + 8 x_2 + 12 x_3 = 0

xxxxxxxxxx
 
print_linear_system(M)

34.2 μs

3×4 Array{Int64,2}:
 1  2  3  -1
 0  9  0  -2
 0  2  3   3

xxxxxxxxxx
 
B = rowop_replace(M, 1, 3, -3)

5.1 μs

3×4 Array{Int64,2}:
 1  0  0  -4
 0  9  0  -2
 0  2  3   3

xxxxxxxxxx
 
C = rowop_replace(B, 3, 1, -1)

8.1 μs

3×4 Array{Int64,2}:
 1  0    0   -4
 0  1  -12  -14
 0  2    3    3

xxxxxxxxxx
 
D = rowop_replace(C, 3, 2, -4)

4.9 μs

3×4 Array{Int64,2}:
 1  0    0   -4
 0  1  -12  -14
 0  0   27   31

xxxxxxxxxx
 
E = rowop_replace(D, 2, 3, -2)

4.2 μs

xxxxxxxxxx
 
# we would get an error here
# rowop_scale(E, 3, 1/27)

316 ns

3×4 Array{Float64,2}:
 1.0  0.0    0.0   -4.0
 0.0  1.0  -12.0  -14.0
 0.0  0.0   27.0   31.0

xxxxxxxxxx
 
# to avoid this we need to cast the matrix E to be a matrix of floating point numbers
F = float(E)

99 ms

3×4 Array{Float64,2}:
 1.0  0.0    0.0   -4.0
 0.0  1.0  -12.0  -14.0
 0.0  0.0    1.0    1.1481481481481481

xxxxxxxxxx
 
G = rowop_scale(F, 3, 1/27)

16.5 ms

3×4 Array{Float64,2}:
 1.0  0.0  0.0  -4.0
 0.0  1.0  0.0  -0.22222222222222143
 0.0  0.0  1.0   1.1481481481481481

xxxxxxxxxx
 
H = rowop_replace(G, 3, 2, 12)

17.8 ms


 x_1 = -4.0

 x_2 = -0.22222222222222143

 x_3 = 1.1481481481481481

xxxxxxxxxx
 
# our final matrix corresponds to a trivial linear system
print_linear_system(H)

149 ms

sol

Float641

-4.0

-0.222222

1.14815

xxxxxxxxxx
 
sol = H[1:3, 4]

41.8 ms

true

xxxxxxxxxx
 
# this definitely should be a solution to the linear system of H
is_solution(H, sol)

225 ms

true

xxxxxxxxxx
 
# as noted in class, it's also a solution to the linear system of M
is_solution(M, sol)

148 ms