也許你可以用這段代碼幫助我。我只得到了 4 個正確且需要的值，但我還需要制作一個包含所有迭代的表，直到我得到這些正確的值，但我不知道如何從我擁有的這個函式中獲取這些迭代的數量：此處提供的示例 在此處輸入影像描述
E <- matrix(c(1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1), byrow = TRUE, nrow = 4)
B <- matrix(c(1.2, 2, 3, 1.5), nrow = 4, byrow = TRUE)
D <- matrix(c(0.18, 0.57, 0.38, 0.42,
0.57, 0.95, 0.70, 0.44,
0.38, 0.70, 0.37, 0.18,
0.42, 0.44, 0.18, 0.40), byrow = TRUE, nrow = 4)
#my matrix
A <- D 0.1*(8 3)*E
A
# Define a convenience function matching `numpy.inner`
inner <- function(a, b) as.numeric(as.numeric(a) %*% as.numeric(b))
conjugate_grad <- function(A, b) {
n <- dim(A)[1]
x <- numeric(n)
z <- b - A %*% x
p <- z
z_old <- inner(z, z)
for (i in 1:60) {
teta <- z_old / inner(p, (A %*% p))
x <- x teta * p
z <- z - teta * A %*% p
z_new <- inner(z, z)
if (sqrt(z_new) < 0.001)
break
beta <- z_new / z_old
p <- z beta * p
z_old <- z_new
}
return(x)
}
conjugate_grad(A,B)
先感謝您！
uj5u.com熱心網友回復：
我認為您的代碼的這種修改為您提供了您想要的：
conjugate_grad <- function(A, b) {
n <- dim(A)[1]
x <- numeric(n)
history <- NULL # Insert this line
z <- b - A %*% x
p <- z
z_old <- inner(z, z)
for (i in 1:60) {
teta <- z_old / inner(p, (A %*% p))
x <- x teta * p
history <- cbind(history, x) # Insert this line
z <- z - teta * A %*% p
z_new <- inner(z, z)
if (sqrt(z_new) < 0.001)
break
beta <- z_new / z_old
p <- z beta * p
z_old <- z_new
}
return(history) # return history instead of x
}
conjugate_grad(A,B)
# [,1] [,2] [,3] [,4]
# [1,] 0.4326431 0.1288703 0.09609509 0.08061088
# [2,] 0.7210718 0.1695202 0.15988971 0.16900513
# [3,] 1.0816077 1.8689058 1.85211220 1.85311415
# [4,] 0.5408039 0.6284172 0.70703113 0.70548042
當您在一個物件中計算結果時，您必須累積結果，這里稱為history.
轉載請註明出處，本文鏈接：https://www.uj5u.com/gongcheng/517644.html
標籤：r功能for循环迭代