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.
# This function creates an augmented matrix from the 2 input matrices with the co-efficients and the constraints
augmentMatrix = function(matrixA,matrixb) {
r = dim(matrixA)[1]
c = dim(matrixb)[2]+dim(matrixA)[2]
matrixB = matrix(c(matrixA,matrixb), ncol=c, nrow=r)
return(matrixB)
}
# This function modifies the matrix and brings all zero rows to the bottom
modifyMatrix=function(matrixA) {
numberOfRows=as.numeric(dim(matrixA)[1])
n=0
z=0
zeroRoWIndex<-c()
nonZeroRoWIndex<-c()
matrixZ<-c()
matrixN<-c()
rowIndex<-as.numeric(0)
r=0
# matrixB<-matrixA[order(matrixA[,1], decreasing = TRUE),]
for(rowIndex in 1:numberOfRows) {
if(all(matrixA[rowIndex,]== 0)) {
zeroRoWIndex<-c(zeroRoWIndex,rowIndex)
matrixZ<-rbind(matrixZ,matrixA[rowIndex,])
z=z+1
}
else {
nonZeroRoWIndex<-c(nonZeroRoWIndex,rowIndex)
matrixN<-rbind(matrixN,matrixA[rowIndex,])
n=n+1
}
}
matrixC<-rbind(matrixN,matrixZ)
return(matrixC)
}
# This is the main function. It converts the input matrices to Reduced Row Echelon Form
reduceMatrix<-function(matrixA, matrixb) {
matrixB<-augmentMatrix(A,b)
matrixB<-modifyMatrix(matrixB)
numberOfRows=(dim(matrixB)[1])
numberOfColumns=(dim(matrixB)[2])
swap_row=c()
colIndex=0
rowIndex=0
r=0
#for(colIndex in seq_along(numberOfColumns))
while(colIndex<=numberOfColumns) {
r=colIndex+1
rowIndex=r
colIndex=colIndex+1
# cat("Processing Col index: ",colIndex,"\n")
k=0
# cat("r: ",r,"\n")
# cat("Row index: ",rowIndex,"\n")
while(rowIndex<=numberOfRows) {
if (all(matrixB[(1:numberOfRows),colIndex]== 0)) {
# cat("All rows in this column :",colIndex,"are 0 \n")
rowIndex=numberOfRows
rowIndex=rowIndex+1
}
else if (matrixB[rowIndex,colIndex]!=0) {
# rowIndex=min(which(matrixB[(r+1:numberOfRows),colIndex]!=0))
# cat("Row index: ",rowIndex,"\n")
# cat("r: ",r,"\n")
swap_row<-matrixB[r,]
matrixB[r,]=matrixB[rowIndex,]
matrixB[rowIndex,]=swap_row
# print(matrixB)
matrixB[r,]=(matrixB[r,]/matrixB[r,colIndex])
for(k in 1:numberOfRows) {
if(k==r){
next
}
else {
# cat("Processing row ",k,"for col ",colIndex,"\n")
matrixB[k,]=matrixB[k,]-(matrixB[k,colIndex]/matrixB[r,colIndex]*matrixB[r,])
}
# print(matrixB)
}
rowIndex=numberOfRows
rowIndex=rowIndex+1
}
else if (matrixB[rowIndex,colIndex]==0) {
# cat("cell in row :",rowIndex ,"and col :",colIndex, "has 0. Nothing to do \n")
rowIndex=rowIndex+1
}
}
}
# print(matrixB)
return(matrixB)
}