快速模拟矩阵操作的方法

Question

快速模拟矩阵操作的方法

4

这是一道老的奥林匹克练习题：

假设你有一个1000x1000的网格，在其中单元格(i,j)包含数字i*j。（行和列从1开始编号。）

每一步，我们从上一步的网格构建一个新的网格，在该网格中，每个单元格(i,j)包含上一个网格中(i,j)及其最多8个邻居的平均值的“邻域平均值”。例如，如果网格角落的4个数字为1、2、5、7，则在下一步中，角落将被计算为(1+2+5+7)/4=3。

最终，我们将达到一个所有数字都相同且网格不再改变的点。目标是找出需要多少步才能达到此点。

我尝试了简单的模拟，但这并不起作用，因为答案似乎是O(n^2)步，而每个模拟步骤需要O(n^2)的处理时间，导致O(n^4)，对于n=1000来说太慢了。

有没有更快的方法？

- CaptainCodeman

你试过从较小的网格开始吗（比如5到20）？很有可能会出现一种模式，让你能够使用公式计算任何网格大小所需的步骤数。这听起来更像是需要数学优化而不是编程的问题。 - Nuclearman

我编写了从n=1到200的所有步数（我们称之为f(n)），但没有找到任何规律。我甚至尝试编写f(n)-f(n-1)，但除了它通常以O（n）的速度增长外，我没有看到任何模式。看起来相当随机。 - CaptainCodeman

你能大概定义一下“太慢”的意思吗？一万亿次迭代是很多的，特别是取决于每次迭代的速度...但这并不是“宇宙寿命”级别的。这需要有多快？ - Nemo

@Nemo 这个问题集的设计是在4小时的时间范围内解决，而且只需要输出结果，所以基本上任何能在一个小时左右内完成的东西都是合理的。根据我的估计，仅运行模拟就需要超过4个小时。 - CaptainCodeman

@CaptainCodeman：这听起来像是微观优化任务...我认为你可以通过重复使用前一行的中间结果来大大减少计算量。让我写一些想法作为答案。 - Nemo

@Nemo，这可能是一个优化任务，但需要比仅仅减少操作次数更好的方法。例如，如果我们可以消除加邻居所需的8倍因子，那可能就足够了。 - CaptainCodeman

3个回答

1

“floor” 步骤让我怀疑解析解不太可能，并且这实际上是一个微优化练习。这是我的想法。

暂时忽略角落和边缘，只有 3996 个需要特殊处理的单元格。

对于内部单元格，您需要添加 9 个元素来获得其下一个状态。但反过来说：每个内部单元格都必须成为 8 个加法的一部分。

或者吗？从三个连续的行 A[i]，B[i] 和 C[i] 开始，计算三个新行：

A'[i] = A[i-1] + A[i] + A[i+1]
B'[i] = B[i-1] + B[i] + B[i+1]
C'[i] = C[i-1] + C[i] + C[i+1]

（注意，你可以使用“滑动窗口”稍微更快地计算每个值，因为 A'[i+1] = A'[i] - A[i-1] + A[i+1]。运算次数相同，但加载次数较少。）

现在，要获取位置 B[j] 的新值，只需计算 A'[j] + B'[j] + C'[j]。

到目前为止，我们还没有节省任何工作量；我们只是重新排列了加法操作。

但是现在，在计算更新的行 B 后，您可以丢弃 A' 并计算下一行：

D'[i] = D[i-1] + D[i] + D[i+1]

您可以使用数组B'和C'来计算行C的新值，而无需重新计算B'或C'。（当然，您会通过将行B'和C'移动成为A'和B'来实现这一点...但是这种方式更容易解释。也许。我想。）

对于每一行，比如B，我们扫描它一次以产生B'，执行2n个算术运算，并再次扫描以计算更新后的B，这也需要2n个操作，因此总共每个元素执行四个加减法，而不是八个。

当然，在实践中，您会在更新B时计算C'，但操作数相同，但局部性更好。

这是我唯一的结构性想法。 SIMD优化专家可能会有其他建议...

- Nemo

谢谢，这是目前为止最好的答案。 - CaptainCodeman

0

如果您查看初始矩阵，您会注意到它是symmetric即m[i][j] = m[j][i]。因此，m[i][j]的邻居将具有与m[j][i]的邻居相同的值，因此您只需要为每个步骤计算略多于一半的矩阵值。

这种优化将每个网格的计算次数从N^2减少到((N^2)+N)/2.

- FuzzyTree

在我的测试中，中心元素确实发生了变化。 - CaptainCodeman

@CaptainCodeman，你看到中心元素改变的奇数值是多少？ - FuzzyTree

几乎所有的情况都是如此。例如，n=5时，中心从9开始，到4结束。这种情况不会发生在开头，但通常你会得到一个低数向右下方从左上角开始传播，直到填满整个网格。 - CaptainCodeman

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Abhishek Bansal · Accepted Answer

有一种略微更快的方法如下：

如果您注意到任何一个不在矩阵边界（x，y）上的单元格，其原始值应为x * y。

此外，第一次迭代后单元格的值应为：

V1 = (  xy    +    x(y+1)   +   x(y-1)
    +(x+1)y  + (x+1)(y+1)  + (x+1)(y-1)
    +(x-1)y  + (x-1)(y+1)  + (x-1)(y-1)
     ) / 9
   = xy

对于左侧竖边上的元素（不包括角落）

v2 = ( xy  + (x-1)y + (x+1)y + x(y+1) + (x-1)(y+1) + (x+1)(y+1) ) / 6
   = xy + x/2.

对于右侧垂直边缘上的元素（不包括角落）

v3 = ( xy  + (x-1)y + (x+1)y + x(y-1) + (x-1)(y-1) + (x+1)(y-1) ) / 6
   = xy - x/2.

同样的，对于顶部和底部的水平边缘和角落同样适用。

因此，在第一次迭代后，只有边框元素会改变其值，非边框元素将保持不变。

对于随后的迭代，这种变化将从矩阵的边框向内传播。

因此，您可以通过仅更改预计在前N/2次迭代中发生更改的那些元素来稍微减少计算量。注意：通过这样做，复杂度不会改变，但常数因子将减少。

另一个可能的方法如下：

您知道，中心最元素在N/2次迭代之前都不会改变。

因此，您可以考虑一种从中心最元素向外开始迭代的方法。

也就是说，如果您能找到一个关于N/2次迭代后元素变化的递增数学公式，则可以将算法的复杂度降低N倍。