How to Use the Linear Algebra Package in Julia
Overview
Julia supports [linear algebra](../../categories/Linear Algebra) as well as MATLAB does, if not better. The intuitive and elegant syntax of Julia gives a feeling that it has been well-designed since its inception1.
Code
julia> A = [
1 0 3
0 5 1
3 1 9
]
3×3 Matrix{Int64}:
1 0 3
0 5 1
3 1 9
As you can see, defining matrices is intuitive and easy from the get-go. Now let’s look at some common functions that should exist. Links to related posts are provided in the subsection titles, and no further explanation is given.
Trace tr()
julia> tr(A)
15
Determinant det()
julia> det(A)
-1.000000000000003
Inverse inv()
julia> inv(A)
3×3 Matrix{Float64}:
-44.0 -3.0 15.0
-3.0 6.10623e-16 1.0
15.0 1.0 -5.0
julia> round.(Int64, inv(A))
3×3 Matrix{Int64}:
-44 -3 15
-3 0 1
15 1 -5
Diagonal matrix and diagonal elements diag(), diagm()
julia> diag(A)
3-element Vector{Int64}:
1
5
9
julia> diagm([1,5,9])
3×3 Matrix{Int64}:
1 0 0
0 5 0
0 0 9
Norm norm()
julia> norm(A, 1)
23.0
Eigenvalues eigvals()
julia> eigvals(A)
3-element Vector{Float64}:
-0.020282065792505244
4.846013411157458
10.174268654635046
julia> eigvecs(A)
3×3 Matrix{Float64}:
-0.944804 0.117887 0.305692
-0.0640048 -0.981459 0.180669
0.321322 0.151132 0.934832
julia> eigmax(A)
10.174268654635046
Matrix Decomposition factorize()
julia> factorize(A)
BunchKaufman{Float64, Matrix{Float64}}
D factor:
3×3 Tridiagonal{Float64, Vector{Float64}}:
-0.0227273 0.0 ⋅
0.0 4.88889 0.0
⋅ 0.0 9.0
U factor:
3×3 UnitUpperTriangular{Float64, Matrix{Float64}}:
1.0 -0.0681818 0.333333
⋅ 1.0 0.111111
⋅ ⋅ 1.0
permutation:
3-element Vector{Int64}:
1
2
3
julia> svd(A)
SVD{Float64, Float64, Matrix{Float64}}
U factor:
3×3 Matrix{Float64}:
-0.305692 0.117887 -0.944804
-0.180669 -0.981459 -0.0640048
-0.934832 0.151132 0.321322
singular values:
3-element Vector{Float64}:
10.174268654635044
4.846013411157461
0.02028206579250516
Vt factor:
3×3 Matrix{Float64}:
-0.305692 -0.180669 -0.934832
0.117887 -0.981459 0.151132
0.944804 0.0640048 -0.321322
Refer to the [Matrix Algebra](../../categories/Matrix Algebra) category’s [Matrix Decomposition](../../categories/Matrix Algebra#Matrix-Decomposition-and-Least-Squares) section. Depending on the form of the matrix, it automatically chooses the appropriate decomposition. Of course, if desired and the conditions are met, one can directly use a specific decomposition function.
Matrix Operations
julia> B = [
1 0 1
1 1 0
2 1 1
]
3×3 Matrix{Int64}:
1 0 1
1 1 0
2 1 1
julia> A + B
3×3 Matrix{Int64}:
2 0 4
1 6 1
5 2 10
julia> A - B
3×3 Matrix{Int64}:
0 0 2
-1 4 1
1 0 8
julia> A * B
3×3 Matrix{Int64}:
7 3 4
7 6 1
22 10 12
julia> A .* B
3×3 Matrix{Int64}:
1 0 3
0 5 0
6 1 9
julia> B / A
3×3 Matrix{Float64}:
-29.0 -2.0 10.0
-47.0 -3.0 16.0
-76.0 -5.0 26.0
julia> B * inv(A)
3×3 Matrix{Float64}:
-29.0 -2.0 10.0
-47.0 -3.0 16.0
-76.0 -5.0 26.0
julia> A / B
ERROR: SingularException(3)
All the operations we consider common sense apply. Division is naturally the same as multiplying by the inverse of multiplication, and if the inverse matrix doesn’t exist like B, it raises a singular exception.
Block Matrix []
Compared to other languages, creating block matrices is incredibly convenient.
julia> [A B]
3×6 Matrix{Int64}:
1 0 3 1 0 1
0 5 1 1 1 0
3 1 9 2 1 1
julia> [A;B]
6×3 Matrix{Int64}:
1 0 3
0 5 1
3 1 9
1 0 1
1 1 0
2 1 1
julia> [A,B]
2-element Vector{Matrix{Int64}}:
[1 0 3; 0 5 1; 3 1 9]
[1 0 1; 1 1 0; 2 1 1]
Placing a space between two matrices stacks them horizontally, while a semicolon stacks them vertically. The comma doesn’t stack matrices but arranges them as an array, following the same syntax used for arrays in general.
Full Code
Content related to complex matrices and inner products is omitted, but included in the full code.
using LinearAlgebra
A = [
1 0 3
0 5 1
3 1 9
]
tr(A)
det(A)
inv(A)
round.(Int64, inv(A))
diag(A)
diagm([1,5,9])
norm(A, 1)
eigvals(A)
eigvecs(A)
eigmax(A)
factorize(A)
svd(A)
B = [
1 0 1
1 1 0
2 1 1
]
det(B)
rank(B)
eigvals(B)
Symmetric(B) # |> issymmetric
transpose(B)
B'
C = [
im im 1
2 im 0
im 1 2
]
C'
B'B
x = [1,2,3]
y = [0,1,2]
x'y
A + B
A - B
A * B
A .* B
B / A
B * inv(A)
[A B]
[A;B]
[A,B]
x' * y
y * x'
Environment
- OS: Windows
- julia: v1.7.0