消除矩阵中的舍入误差

Question

我的问题是，我有一个矩阵，其中所有行的总和和所有列的总和都为零。所有数字都舍入到x个小数位数。

然后，我将整个矩阵乘以介于0和1之间的数字（例如1/6），并将所有数字舍入到x个小数位数。现在，我不能确定行和列的总和是否仍然是零。我希望行和列的总和再次为零，并进行最小可能的调整（或者至少进行非常小的调整）。

是否存在一种算法可以解决这样的问题？

例如（非常简单的）矩阵：

    200  -200  0

    400  400  -800

   -600 -200  800

round2( (1/6)*matrix)

33.33  -33.33  0   

66.67  66.67   -133.33

-100   -33.33  133.33

- user1879280

---

我只需要将行和列相加，而不是测试它们是否等于零，而是测试总和的绝对值是否小于某个特定的容差值 - 在这种情况下，也许是 abs(sum) <= 0.01。 - Blazemonger

这不是一个“算法”问题。你通过四舍五入引入了一个问题，无论你如何“修复”，都会引入其他问题，例如打破矩阵中某些元素之间的对称性。难道你不能将四舍五入限制在“显示”的值上，同时保留数学处理的完整值吗？你仍然会有一些“噪音”，使得总和可能非零，但这个问题应该通过将“零”定义为“小于某个公差”来处理。 - Bert te Velde

Answer 1

这不是一个解决方案，只是对你所尝试实现的更数学化的描述（不判断此举是否正确）：

由于你将所有数字都舍入到x位小数，因此我们可以将这些数字视为整数（只需将它们乘以10^x）。

现在，你正在尝试解决以下问题：

给定矩阵。

A11+Adj11   A12+Adj12   ...   A1n+Adj1n
A21+Adj21   A22+Adj22   ...   A2n+Adj2n
A31+Adj31   A32+Adj32   ...   A3n+Adj3n
...         ...         ...   ...
Am1+Adjm1   Am2+Adjm2   ...   Amn+Adjmn

其中A11..Amn是常数整数，

查找整数Adj11...Adjmn

最小化sum(abs(Adjxy))

（或者您可能更喜欢：最小化sum((Adjxy)^2)）

遵循以下规定：

- for each row m: Adjm1+Adjm2+...+Adjmn = - (Am1+Am2+...+Amn)
- for each col n: Adj1n+Adj2n+...+Adjmn = - (A1n+A2n+...+Amn)

这是一个整数规划问题，有m*n个变量和m+n个约束条件。你要尝试最小化的函数不是线性的。

恐怕这个问题远非琐碎。我认为你最好在https://math.stackexchange.com/上发布它。

- Lior Kogan

Answer 2

您在这里经历的是精度误差。除非您根本不进行四舍五入，否则无法做任何事情。这类似于将照片保存为256色图像。您正在失去信息（基本上是由于离散化而导致的精度），您无能为力。对于图片，有算法可以使图像看起来更平滑/更接近原始图像（例如抖动），但是单个值数字没有这样的东西。

可能的解决方案（实际上只有一种，有两种不同的可视化方式）：

- 仅用于显示时进行四舍五入。用户应该能够理解数字被截断/四舍五入。在您的示例中，显然6.67实际上将是6.66666...。 - 根本不进行四舍五入，只需在固定小数位数后截断数字（如果需要，添加...；这实际上与另一种解决方案类似）。

一般来说，如果您想解决线性方程（或数学问题），请始终使用可用的（合理的；性能方面）具有最大精度的数据类型，通常是单精度或双精度值。否则，您会引入越来越严重的误差边界，随着计算次数的增加而变得越来越糟糕。

- Mario

---

总的来说是好建议，但“你无能为力”有些言过其实。例如，您可以寻找最小二乘解，即导致行和列总和为0并最小化平方差之和的调整组合。 - j_random_hacker

可能吧，我想你只是把问题转移了，如果你不幸的话，你的计算只会变得更加复杂和特定于情况。 - Mario

Answer 3

常见的确保小的舍入误差不会导致总和的大误差的方法是，检查每个部分和的误差是否不会变得太大。

对于一个一维向量 [a[1], a[2], ..., a[n]]，你可以计算部分和 [a[1], a[1]+a[2], ..., a[1]+a[2]+...+a[n]]，然后将其乘以一个倍数再通过将当前元素减去前一个元素来恢复原向量：[a[1]*b, (a[1]+a[2])*b-a[1]*b, ..., (a[1]+a[2]+...+a[n])*b-(a[1]+a[2]+...+a[n-1])*b]。通过使用这个技巧，任何部分和的误差都不会超过10^(-x)。

你可以根据以下3个步骤，将此方法适应于二维矩阵：

partial_sum(M) =
  for i = 0 to n-1 do
    for j = 1 to m-1 do
      M[i][j] += M[i][j-1]
    done
  done
  for i = 0 to n-1 do
    for j = 1 to m-1 do
      M[j][i] += M[j-1][i]
    done
  done

multiply(M, a) =
  for i = 0 to n-1 do
    for j = 0 to m-1 do
      M[i][j] *= a
    done
  done

restore(M) =
  for i = 0 to n-1 do
    for j = 1 to m-1 do
      M[i][j] -= M[i][j-1]
    done
  done
  for i = 0 to n-1 do
    for j = 1 to m-1 do
      M[j][i] -= M[j-1][i]
    done
  done

- Thomash

Answer 4

在使用浮点数时，无法消除舍入误差。您最好的解决方案可能是在矩阵中使用整数，然后将最终的1/6应用于结果。

- Mark Ransom