Set up a system of equations with 3 variables and 3 constraints and solve for x. Please write a function in R that will take two variables (matrix A & constraint vector b) and solve using elimination. Your function should produce the right answer for the system of equations for any 3-variable, 3-equation system. You don’t have to worry about degenerate cases and can safely assume that the function will only be tested with a system of equations that has a solution. Please note that you do have to worry about zero pivots, though. Please note that you should not use the built-in function solve to solve this system or use matrix inverses. The approach that you should employ is to construct an Upper Triangular Matrix and then back-substitute to get the solution. Alternatively, you can augment the matrix A with vector b and jointly apply the Gauss Jordan elimination procedure.
I’m employing the Upper triangular Matrix method and then using substition method to find the solution vectors. 3 for loops are being used.
First for loop is used to iterate the column except the last one. because we dont want to do transformation to the last column. Second for loop is used to iterate through the row( for 1st column, you need to transform row2, row3. Similarly for 2nd column, you need tranform row3). And the last ‘for loop’ is used to iterate through the row value( This is for doing row tranformation operation)
Also, i introduced a verbose flag to see how the matrix is getting transformed into the final upper triangular matrix.
Finally i’m using substitution technique to find the the value of solution vector.
UpperTrigSolve <- function(A, b, verbose=0) {
n <- nrow(A)
if(verbose==1){
print(A)
}
# Column iteration except the last column.
for(k in 1:(n-1)) {
# Row Iteration
for (i in (k+1):n) {
# Setting the multiplier
multiplier <- A[i,k]/A[k,k]
# Setting the element[i,k] to 0
A[i,k] <- 0
# Applying the row operation with pre-defined multiplier by iterating the rest of elements in the row
for (j in (k+1):n) {
A[i,j] <- A[i,j] - multiplier * A[k,j]
if(verbose==1){
print("=============")
print(A)
}
}
b[i] <- b[i] - multiplier * b[k]
}
}
# Substitution
x <- vector(mode="numeric",length=n)
x[n] <- b[n]/A[n,n]
x[n-1] <- (b[n-1] - (A[n-1,n] * x[n])) / A[n-1,n-1]
x[n-2] <- (b[n-2] - (A[n-2,n-1] * x[n-1]) - (A[n-2,n] * x[n])) / A[n-2,n-2]
return(x)
}
Create the test matrix with 3 X 3 size and with respective constraints. Pass both of them into the above defined function.
matrix_elements = c(1,2,-1,1,-1,-2,3,5,4)
MatrixA = matrix(matrix_elements,nrow=3, ncol=3)
b <- c(1,2,6)
solnvector<- UpperTrigSolve(MatrixA,b,verbose=1)
## [,1] [,2] [,3]
## [1,] 1 1 3
## [2,] 2 -1 5
## [3,] -1 -2 4
## [1] "============="
## [,1] [,2] [,3]
## [1,] 1 1 3
## [2,] 0 -3 5
## [3,] -1 -2 4
## [1] "============="
## [,1] [,2] [,3]
## [1,] 1 1 3
## [2,] 0 -3 -1
## [3,] -1 -2 4
## [1] "============="
## [,1] [,2] [,3]
## [1,] 1 1 3
## [2,] 0 -3 -1
## [3,] 0 -1 4
## [1] "============="
## [,1] [,2] [,3]
## [1,] 1 1 3
## [2,] 0 -3 -1
## [3,] 0 -1 7
## [1] "============="
## [,1] [,2] [,3]
## [1,] 1 1 3.000000
## [2,] 0 -3 -1.000000
## [3,] 0 0 7.333333
solnvector
## [1] -1.5454545 -0.3181818 0.9545455
As a side note and the reference, There is another method called “Reduced Echelon form” in which we skip the substition steps. In this approach, we try to make lower and upper side of the matrix be zero and all elements on the diagonal will be 1.
Reference link: https://www.youtube.com/watch?v=W95p0E-CB_s