caracas
vignettes/a11-linear-algebra.Rmd
a11-linear-algebra.Rmd
This vignette is based on caracas
version
r packageVersion("caracas")
. caracas
is
avavailable on CRAN at [https://cran.r-project.org/package=caracas] and on
github at [https://github.com/r-cas/caracas].
We now show different ways to create a symbolic matrix:
A <- matrix(c("a", "b", "0", "1"), 2, 2) %>% as_sym()
A
#> c: ⎡a 0⎤
#> ⎢ ⎥
#> ⎣b 1⎦
A <- matrix_(c("a", "b", "0", "1"), 2, 2) # note the '_' postfix
A
#> c: ⎡a 0⎤
#> ⎢ ⎥
#> ⎣b 1⎦
A <- as_sym("[[a, 0], [b, 1]]")
A
#> c: ⎡a 0⎤
#> ⎢ ⎥
#> ⎣b 1⎦
A2 <- matrix_(c("a", "b", "c", "1"), 2, 2)
A2
#> c: ⎡a c⎤
#> ⎢ ⎥
#> ⎣b 1⎦
B <- matrix_(c("a", "b", "0",
"c", "c", "a"), 2, 3)
B
#> c: ⎡a 0 c⎤
#> ⎢ ⎥
#> ⎣b c a⎦
b <- matrix_(c("b1", "b2"), nrow = 2)
D <- diag_(c("a", "b")) # note the '_' postfix
D
#> c: ⎡a 0⎤
#> ⎢ ⎥
#> ⎣0 b⎦
Note that a vector is a matrix in which one of the dimensions is one.
A + A2
#> c: ⎡2⋅a c⎤
#> ⎢ ⎥
#> ⎣2⋅b 2⎦
A %*% B
#> c: ⎡ 2 ⎤
#> ⎢ a 0 a⋅c ⎥
#> ⎢ ⎥
#> ⎣a⋅b + b c a + b⋅c⎦
reciprocal_matrix(A2)
#> c: ⎡1 1⎤
#> ⎢─ ─⎥
#> ⎢a c⎥
#> ⎢ ⎥
#> ⎢1 ⎥
#> ⎢─ 1⎥
#> ⎣b ⎦
reciprocal_matrix(A2, num = "2*a")
#> c: ⎡ 2⋅a⎤
#> ⎢ 2 ───⎥
#> ⎢ c ⎥
#> ⎢ ⎥
#> ⎢2⋅a ⎥
#> ⎢─── 2⋅a⎥
#> ⎣ b ⎦
M <- as_sym("[[1, 2, 3], [4, 5, 6]]")
pinv(M)
#> c: ⎡-17 ⎤
#> ⎢──── 4/9 ⎥
#> ⎢ 18 ⎥
#> ⎢ ⎥
#> ⎢-1/9 1/9 ⎥
#> ⎢ ⎥
#> ⎢ 13 ⎥
#> ⎢ ── -2/9⎥
#> ⎣ 18 ⎦
B <- as_sym("[[7], [8]]")
B
#> c: [7 8]ᵀ
z <- do_la(M, "pinv_solve", B)
print(z, rowvec = FALSE) # Do not print column vectors as transposed row vectors
#> c: ⎡ w₀ ₀ w₁ ₀ w₂ ₀ 55 ⎤
#> ⎢ ──── - ──── + ──── - ── ⎥
#> ⎢ 6 3 6 18 ⎥
#> ⎢ ⎥
#> ⎢ w₀ ₀ 2⋅w₁ ₀ w₂ ₀ 1⎥
#> ⎢- ──── + ────── - ──── + ─⎥
#> ⎢ 3 3 3 9⎥
#> ⎢ ⎥
#> ⎢ w₀ ₀ w₁ ₀ w₂ ₀ 59 ⎥
#> ⎢ ──── - ──── + ──── + ── ⎥
#> ⎣ 6 3 6 18 ⎦
Below we present convenient functions for performing linear algebra
operations. If you need a function in SymPy for which we have not
supplied a convinience function (see https://docs.sympy.org/latest/modules/matrices/matrices.html),
you can still call it with the do_la()
(short for “do
linear algebra”) function presented at the end of this section.
A <- matrix(c("a", "0", "0", "1"), 2, 2) %>% as_sym()
A
#> c: ⎡a 0⎤
#> ⎢ ⎥
#> ⎣0 1⎦
qr_res <- QRdecomposition(A)
qr_res$Q
#> c: ⎡ a ⎤
#> ⎢─── 0⎥
#> ⎢│a│ ⎥
#> ⎢ ⎥
#> ⎣ 0 1⎦
qr_res$R
#> c: ⎡│a│ 0⎤
#> ⎢ ⎥
#> ⎣ 0 1⎦
eigenval(A)
#> [[1]]
#> [[1]]$eigval
#> c: a
#>
#> [[1]]$eigmult
#> [1] 1
#>
#>
#> [[2]]
#> [[2]]$eigval
#> c: 1
#>
#> [[2]]$eigmult
#> [1] 1
evec <- eigenvec(A)
evec
#> [[1]]
#> [[1]]$eigval
#> c: 1
#>
#> [[1]]$eigmult
#> [1] 1
#>
#> [[1]]$eigvec
#> c: [0 1]ᵀ
#>
#>
#> [[2]]
#> [[2]]$eigval
#> c: a
#>
#> [[2]]$eigmult
#> [1] 1
#>
#> [[2]]$eigvec
#> c: [1 0]ᵀ
evec1 <- evec[[1]]$eigvec
evec1
#> c: [0 1]ᵀ
simplify(evec1)
#> c: [0 1]ᵀ
lapply(evec, function(l) simplify(l$eigvec))
#> [[1]]
#> c: [0 1]ᵀ
#>
#> [[2]]
#> c: [1 0]ᵀ
do_la
short for “do linear algebra”
args(do_la)
#> function (x, slot, ...)
#> NULL
The above functions can be called:
do_la(A, "QRdecomposition") # == QRdecomposition(A)
#> $Q
#> c: ⎡ a ⎤
#> ⎢─── 0⎥
#> ⎢│a│ ⎥
#> ⎢ ⎥
#> ⎣ 0 1⎦
#>
#> $R
#> c: ⎡│a│ 0⎤
#> ⎢ ⎥
#> ⎣ 0 1⎦
do_la(A, "inv") # == inv(A)
#> c: ⎡1 ⎤
#> ⎢─ 0⎥
#> ⎢a ⎥
#> ⎢ ⎥
#> ⎣0 1⎦
do_la(A, "eigenvec") # == eigenvec(A)
#> [[1]]
#> [[1]]$eigval
#> c: 1
#>
#> [[1]]$eigmult
#> [1] 1
#>
#> [[1]]$eigvec
#> c: [0 1]ᵀ
#>
#>
#> [[2]]
#> [[2]]$eigval
#> c: a
#>
#> [[2]]$eigmult
#> [1] 1
#>
#> [[2]]$eigvec
#> c: [1 0]ᵀ
do_la(A, "eigenvals") # == eigenval(A)
#> [[1]]
#> [[1]]$eigval
#> c: a
#>
#> [[1]]$eigmult
#> [1] 1
#>
#>
#> [[2]]
#> [[2]]$eigval
#> c: 1
#>
#> [[2]]$eigmult
#> [1] 1
do_la(A, "rank")
#> c: 2
B <- matrix(c(1, 3, 0, -2, -6, 0, 3, 9, 6), nrow = 3) %>% as_sym()
B
#> c: ⎡1 -2 3⎤
#> ⎢ ⎥
#> ⎢3 -6 9⎥
#> ⎢ ⎥
#> ⎣0 0 6⎦
columnspace(B)
#> c: ⎡1 3⎤
#> ⎢ ⎥
#> ⎢3 9⎥
#> ⎢ ⎥
#> ⎣0 6⎦
rowspace(B)
#> c: [1 -2 3 0 0 6]
x <- nullspace(B)
x
#> c: [2 1 0]ᵀ
rref(B)
#> $mat
#> c: ⎡1 -2 0⎤
#> ⎢ ⎥
#> ⎢0 0 1⎥
#> ⎢ ⎥
#> ⎣0 0 0⎦
#>
#> $pivot_vars
#> [1] 1 3
B %*% x
#> c: [0 0 0]ᵀ
B <- t(as_sym("[[ 2, 3, 5 ], [4, 6, 10], [8, 3, 6], [8, 3, 6] ]"))
B
#> c: ⎡2 4 8 8⎤
#> ⎢ ⎥
#> ⎢3 6 3 3⎥
#> ⎢ ⎥
#> ⎣5 10 6 6⎦
singular_values(B)
#> [[1]]
#> c: ______________
#> ╲╱ √30446 + 204
#>
#> [[2]]
#> c: ______________
#> ╲╱ 204 - √30446
#>
#> [[3]]
#> c: 0
#>
#> [[4]]
#> c: 0