a <- matrix(1:12, ncol = 3, nrow = 4)
a
## [,1] [,2] [,3]
## [1,] 1 5 9
## [2,] 2 6 10
## [3,] 3 7 11
## [4,] 4 8 12
at<-t(a)
aat<-a%*%at
aat
## [,1] [,2] [,3] [,4]
## [1,] 107 122 137 152
## [2,] 122 140 158 176
## [3,] 137 158 179 200
## [4,] 152 176 200 224
ata<-at%*%a
ata
## [,1] [,2] [,3]
## [1,] 30 70 110
## [2,] 70 174 278
## [3,] 110 278 446
a <- matrix(c(1,2,3,2,4,5,3,5,6), ncol = 3, nrow = 3)
a
## [,1] [,2] [,3]
## [1,] 1 2 3
## [2,] 2 4 5
## [3,] 3 5 6
at<-t(a)
att<-a%*%at
att
## [,1] [,2] [,3]
## [1,] 14 25 31
## [2,] 25 45 56
## [3,] 31 56 70
ata<-at%*%a
ata
## [,1] [,2] [,3]
## [1,] 14 25 31
## [2,] 25 45 56
## [3,] 31 56 70
Write an R function to factorize a square matrix A into LU or LDU, whichever you prefer.
x<-matrix(c(11,4,9,8,88,7,1,77,32,3,99,67,67,45,3,5),ncol=4,nrow=4)
x
## [,1] [,2] [,3] [,4]
## [1,] 11 88 32 67
## [2,] 4 7 3 45
## [3,] 9 1 99 3
## [4,] 8 77 67 5
l<-function(b){
A<-diag(nrow(b))
for (i in 2:(nrow(b)))
{for (j in i:(nrow(b)))
{
temp<-b[j,i-1]/b[i-1,i-1]
b[j,]<-b[j,]-(b[j,i-1]/b[i-1,i-1])*b[i-1,]
A[j,i-1]<-temp
}
}
print(b)
print(A)
return(A%*%b)
}
A<-l(x)
## [,1] [,2] [,3] [,4]
## [1,] 11 88 32.000000 67.00000
## [2,] 0 -25 -8.636364 20.63636
## [3,] 0 0 97.345455 -110.42545
## [4,] 0 0 0.000000 11.51207
## [,1] [,2] [,3] [,4]
## [1,] 1.0000000 0.00 0.0000000 0
## [2,] 0.3636364 1.00 0.0000000 0
## [3,] 0.8181818 2.84 1.0000000 0
## [4,] 0.7272727 -0.52 0.4030631 1
A
## [,1] [,2] [,3] [,4]
## [1,] 11 88 32 67
## [2,] 4 7 3 45
## [3,] 9 1 99 3
## [4,] 8 77 67 5