Constructing matrix A:
A <- matrix(
c(2, 4, -4, 1,
3, 6, 1, -2,
-1, 1, 2, 3,
1, 1, -4, 1), byrow=T, nrow=4, ncol=4)
A
## [,1] [,2] [,3] [,4]
## [1,] 2 4 -4 1
## [2,] 3 6 1 -2
## [3,] -1 1 2 3
## [4,] 1 1 -4 1
Constructing vector b:
b <- matrix(c(0,-7,4,2), nrow=4, ncol=1)
b
## [,1]
## [1,] 0
## [2,] -7
## [3,] 4
## [4,] 2
Implementing GEM method.
p <- nrow(A)
U.pls <- cbind(A, b)
U.pls[1,] <- U.pls[1,]/U.pls[1,1]
i <- 2
while (i < p+1) {
j <- i
while (j < p+1) {
U.pls[j, ] <- U.pls[j, ] - U.pls[i-1, ] * U.pls[j, i-1]
j <- j+1
}
while (U.pls[i,i] == 0) {
U.pls <- rbind(U.pls[-i,], U.pls[i,])
}
U.pls[i,] <- U.pls[i,]/U.pls[i,i]
i <- i+1
}
for (i in p:2){
for (j in i:2-1) {
U.pls[j, ] <- U.pls[j, ] - U.pls[i, ] * U.pls[j, i]
}
}
U.pls
## [,1] [,2] [,3] [,4] [,5]
## [1,] 1 0 0 0 1.000000e+00
## [2,] 0 1 0 0 -1.000000e+00
## [3,] 0 0 1 0 -4.440892e-16
## [4,] 0 0 0 1 2.000000e+00
The final result of the system is:
| Variable | Value |
|---|---|
| \(x_1\) | 1 |
| \(x_2\) | -1 |
| \(x_3\) | -4.4408921^{-16} |
| \(x_4\) | 2 |
Before everything, lets declare a function for GEM:
GEM <- function(A, b){
p <- nrow(A)
U.pls <- cbind(A, b)
U.pls[1,] <- U.pls[1,]/U.pls[1,1]
i <- 2
while (i < p+1) {
j <- i
while (j < p+1) {
U.pls[j, ] <- U.pls[j, ] - U.pls[i-1, ] * U.pls[j, i-1]
j <- j+1
}
while (U.pls[i,i] == 0) {
U.pls <- rbind(U.pls[-i,], U.pls[i,])
}
U.pls[i,] <- U.pls[i,]/U.pls[i,i]
i <- i+1
}
for (i in p:2){
for (j in i:2-1) {
U.pls[j, ] <- U.pls[j, ] - U.pls[i, ] * U.pls[j, i]
}
}
return(U.pls)
}
Constructing matrix A:
A <- matrix(
c(1, 1,
-10, 1e5), byrow=T, nrow=2, ncol=2)
A
## [,1] [,2]
## [1,] 1 1e+00
## [2,] -10 1e+05
Constructing vector b:
b <- matrix(c(3, 1e5), nrow=2, ncol=1)
b
## [,1]
## [1,] 3e+00
## [2,] 1e+05
Implementing GEM method.
U.pls <- GEM(A, b)
U.pls
## [,1] [,2] [,3]
## [1,] 1 0 1.9998
## [2,] 0 1 1.0002
The final result of the system is:
| Variable | Value |
|---|---|
| \(x_1\) | 1.9998 |
| \(x_2\) | 1.0002 |
Constructing matrix A:
A <- matrix(
c(3, 2,
7, 9), byrow=T, nrow=2, ncol=2)
A
## [,1] [,2]
## [1,] 3 2
## [2,] 7 9
Constructing vector b:
b <- matrix(c(13, 713), nrow=2, ncol=1)
b
## [,1]
## [1,] 13
## [2,] 713
Implementing GEM method.
U.pls <- GEM(A, b)
U.pls
## [,1] [,2] [,3]
## [1,] 1 0 -100.6923
## [2,] 0 1 157.5385
The final result of the system is:
| Variable | Value |
|---|---|
| \(x_1\) | -100.6923077 |
| \(x_2\) | 157.5384615 |
Constructing matrix \(A\):
A <- matrix(
c(7, 90,
3, 2), byrow=T, nrow=2, ncol=2)
A
## [,1] [,2]
## [1,] 7 90
## [2,] 3 2
Constructing vector \(b\):
b <- matrix(c(713, 13), nrow=2, ncol=1)
b
## [,1]
## [1,] 713
## [2,] 13
Implementing GEM method.
U.pls <- GEM(A, b)
U.pls
## [,1] [,2] [,3]
## [1,] 1 0 -1
## [2,] 0 1 8
The final result of the system is:
| Variable | Value |
|---|---|
| \(x_1\) | -1 |
| \(x_2\) | 8 |
Constructing matrix \(A\):
A <- matrix(
c(1, 2, 4, 17,
3, 6, -12, 3,
2, 3, -3, 2,
0, 2, -2, 6), byrow=T, nrow=4, ncol=4)
A
## [,1] [,2] [,3] [,4]
## [1,] 1 2 4 17
## [2,] 3 6 -12 3
## [3,] 2 3 -3 2
## [4,] 0 2 -2 6
Constructing vector \(b\):
b <- matrix(c(17, 3, 3, 4), nrow=4, ncol=1)
b
## [,1]
## [1,] 17
## [2,] 3
## [3,] 3
## [4,] 4
Using Matrix library and \(lu()\) function, we create a model. We can access \(L\) and \(U\) matrix by commands below:
Matrix::lu(A) -> lum
Matrix::expand(lum) -> model
L <- model$L
U <- model$U
P <- model$P
L
## 4 x 4 Matrix of class "dtrMatrix" (unitriangular)
## [,1] [,2] [,3] [,4]
## [1,] 1.0000000 . . .
## [2,] 0.0000000 1.0000000 . .
## [3,] 0.3333333 0.0000000 1.0000000 .
## [4,] 0.6666667 -0.5000000 0.5000000 1.0000000
U
## 4 x 4 Matrix of class "dtrMatrix"
## [,1] [,2] [,3] [,4]
## [1,] 3 6 -12 3
## [2,] . 2 -2 6
## [3,] . . 8 16
## [4,] . . . -5
P
## 4 x 4 sparse Matrix of class "pMatrix"
##
## [1,] . . | .
## [2,] | . . .
## [3,] . . . |
## [4,] . | . .
Lets confirm the result is valid by calculating residuals:
P %*% L %*% U
## 4 x 4 Matrix of class "dgeMatrix"
## [,1] [,2] [,3] [,4]
## [1,] 1 2 4 17
## [2,] 3 6 -12 3
## [3,] 2 3 -3 2
## [4,] 0 2 -2 6
Lets finally solve the equation \(Ax = b\).
solve(A, b)
## [,1]
## [1,] 2
## [2,] -1
## [3,] 0
## [4,] 1
Constructing matrix \(A\):
A <- matrix(
c(1, 1, 1, 1,
1, 2, 3, 4,
1, 3, 6, 10,
1, 4, 10, 20), byrow=T, nrow=4, ncol=4)
A
## [,1] [,2] [,3] [,4]
## [1,] 1 1 1 1
## [2,] 1 2 3 4
## [3,] 1 3 6 10
## [4,] 1 4 10 20
The Choleski factorization simply can be done by R-base function \(chol()\)!
chol(A)
## [,1] [,2] [,3] [,4]
## [1,] 1 1 1 1
## [2,] 0 1 2 3
## [3,] 0 0 1 3
## [4,] 0 0 0 1
A
A <- matrix(
c(-0.4, 1, -0.8,
1.2, -2, 1.4,
-0.6, 1, -0.2), byrow=T, nrow=3, ncol=3)
A
## [,1] [,2] [,3]
## [1,] -0.4 1 -0.8
## [2,] 1.2 -2 1.4
## [3,] -0.6 1 -0.2
The condition number of a matrix simply can be calculated by R-base function \(kappa()\)!
kappa(A)
## [1] 20.90613