这里的浮点精度是怎么回事？

Question

的观察

提交的 R 解决方案是一个可行的解决方案，但我们中的一些人（也许只有我）似乎对为什么需要初始 X=m 重新分配感到困惑。

代码被@Giuseppe 简化了一些，所以我会为 reader 写一些评论。

function(m){
            X=m
            # Re-assign input m as X

            while(any(X-(X=X%*%m))) 0
            # Instead of doing the meat of the calculation in the code block after `while`
            # OP exploited its infinite looping properties to perform the
            # calculations within the condition check.
            # `-` here is an abuse of inequality check and relies on `any` to coerce
            # the numeric to logical. See `as.logical(.Machine$double.xmin)`
            # The code basically multiplies the matrix `X` with the starting matrix `m`            
            # Until the condition is met: X == X%*%m

            X
            # Return result
           }

据我所知。乘 X%*%m 等价于 X%*%X，因为 X 只是 m 的迭代自乘版本。一旦矩阵收敛，乘以 m 或 X 的额外副本不会改变其值。将上述函数定义为v后，参见线性代数教科书或v(m)%*%v(m)%*%v(m)%*%v(m)%*%v(m)%*%m%*%m。好玩吧？

那么问题来了，为什么@CodesInChaos的实现这个想法行不通呢？

function(m){while(any(m!=(m=m%*%m)))0 m}

这是浮点精度问题造成的吗？或者这是由代码中的函数引起的，例如不等式检查或.Primitive("any")？我不相信这是由 as.logical 引起的，因为 R 似乎将小于 .Machine$double.xmin 的错误强制为 0.

以上为演示。我们只是循环并计算 m 和 m%*%m 之间的差异。当我们尝试收敛随机矩阵时，此误差变为 0。它似乎收敛然后最终根据输入吹到 0/INF。

mat = matrix(c(7/10, 4/10, 3/10, 6/10), 2, 2, byrow = T)

m = mat
for (i in 1:25) {
  m = m%*%m
  cat("Mean Error:", mean(m-(m=m%*%m)), 
      "\n Float to Logical:", as.logical(m-(m=m%*%m)),
      "\n iter", i, "\n")
}

关于为什么这是一个浮点数学问题的一些额外想法

1) 循环表明这可能不是 any 或任何逻辑 check/conversion 步骤的问题，而是与浮点矩阵数学有关。

2) @user202729 在原始线程中的评论说这个问题在 Jelly 中仍然存在，一种代码高尔夫语言使人们更加相信这可能是一个浮点问题。

Answer 1

这确实是一道浮点数学题。要查看它，请查看此函数的结果：

test2 <- function(m) {
  c <- 0
  res <- list()
  while (any(m!=(m=m%*%m))) {
    c <- c + 1
    res[[c]] <- m
  }
  print(c)
  res
}

要测试具有一定公差的相等性，您可以使用：

test3 <- function(m) {
  while (!isTRUE(all.equal(m, m <- m %*% m))) 0
  m
}

Answer 2

不同的方法迭代不同的函数，都从种子值 m 开始。如果给定的不动点稳定并且种子在该不动点的吸引力盆地内，函数迭代只会收敛到给定的不动点。

在原始代码中，您正在迭代函数

f <- function(X) X %*% m

极限矩阵是一个稳定的 fixed-point 假设（在 Code Gulf 问题中陈述）存在 well-defined 极限。由于函数定义取决于 m，不动点是 m.

的函数也就不足为奇了

另一方面，使用 m = m %*% m 的建议变化是通过迭代函数

获得的

g <- function(X) X %*% X

注意所有幂等矩阵都是这个函数的不动点，但显然它们不可能都是稳定的不动点。显然，原定函数中的极限矩阵并不是g的稳定不动点（虽然是不动点）

要真正解决这个问题，您需要深入研究函数迭代下的矩阵不动点理论，以显示 为什么 [=16= 情况下的不动点] 不稳定。

这里的浮点精度是怎么回事？

What is going on with floating point precision here?

r

floating-point

matrix-multiplication

stochastic

stochastic-process