20

在我的编码过程中,我遇到了以下问题:在二维空间中找到具有最高粒子密度的固定大小的区域。可以认为粒子通常随机分布在整个空间中,但理论上应该有一些区域具有更高的密度。

例如,100 个粒子随机放置在 500x500 的 2D 网格中,我需要找到粒子最多(密度最高)的 50x50 区域。

除了蛮力测试每个可能的区域(在这种情况下大约超过 200000 个区域)之外,还有其他方法来计算最佳区域吗?对于 n 长度的轴,这将在 O(n^2) 处放大。

4

3 回答 3

10

算法 1

创建一个 500x500 二维数组,其中每个单元格包含该单元格中粒子数的计数。然后将该数组与 50x50 内核卷积,生成的数组将包含每个单元格中 50x50 区域中的粒子数。然后找到具有最大值的单元格。

如果您使用 50x50 的框作为区域,则可以将内核分解为两个单独的卷积,每个轴一个。生成的算法是 O(n^2) 空间和时间,其中 n 是您正在搜索的 2D 空间的宽度和高度。

提醒一下,带有 boxcar 函数的一维卷积可以在 O(n) 时间和空间内完成,并且可以就地完成。设 x(t) 为 t=1..n 的输入,设 y(t) 为输出。为 t<1 和 t>n 定义 x(t)=0 和 y(t)=0。将内核 f(t) 定义为 0..d-1 的 1 和其他地方的 0。卷积的定义为我们提供了以下公式:

y(t) = sum ix(ti) * f(i) = sum i=0..d-1 x(ti)

这看起来需要时间 O(n*d),但我们可以将其重写为递归:

y(t) = y(t-1) + x(t) - x(td)

这表明一维卷积是 O(n),与 d 无关。要执行二维卷积,您只需对每个轴执行一维卷积。之所以可行,是因为可以分解 boxcar 内核:通常,大多数内核无法分解。高斯核是另一个可以分解的核,这就是为什么图像编辑程序中的高斯模糊如此之快。

对于您指定的数字类型,这将非常快。500x500 是一个极小的数据集,你的计算机最多可以在几毫秒内检查 202,500 个区域。您将不得不问自己是否值得额外的数小时、数天或数周的时间来进一步优化。

这与justhalf的解决方案相同,除了由于分解卷积,区域大小不会影响算法的速度。

算法 2

假设至少有一个点。不失一般性,将二维空间视为整个平面。令d为区域的宽度和高度。设 N 为点数。

引理:存在一个密度最大的区域,其左边缘有一个点。

证明:令 R 为最大密度区域。设 R' 是同一个区域,向右平移 R 的左边缘和 R 中最左边的点之间的距离。R 中的所有点也必须位于 R' 中,因此 R' 也是最大密度的区域。

算法

  1. 将所有点插入到 KD 树中。这可以在 O(N log 2 N) 时间内完成。

  2. 对于每个点,请考虑宽度为d且高度为 2 d的区域,其中该点以该区域的左边缘为中心。将此区域称为 R。

  3. 在 KD 树中查询区域 R 中的点。将此集合称为 S。这可以在 O(N 1/2 +|S|) 时间内完成。

  4. 找到 R 的 dxd 子区域,其中包含 S 中最大数量的点。这可以在 O(|S| log |S|) 时间内完成,方法是按 y 坐标对 S 进行排序,然后执行线性扫描。

结果算法的时间为 O(N 3/2 + N |S​​| log |S|)。

比较

当密度高时,算法#1 优于算法#2。算法#2 仅在粒子密度非常低的情况下更为优越,并且算法#2 优越的密度随着总板尺寸的增加而降低。

请注意,可以认为连续情况的密度为零,此时只有算法#2 有效。

于 2013-10-23T05:40:44.010 回答
1

我不知道您使用什么蛮力方法,但最蛮力的方法是O(n^2 d^2),通过及时迭代每个区域O(n^2),然后及时计算该区域中的粒子数,即您O(d^2)所在d区域的大小。

这个问题和这个问题完全一样:Rat Attack,因为区域面积是固定的,所以密度和计数一样,对于它的解是O(n^2 + k*d^2),其中

  1. n是整个区域的大小(边长)
  2. k是粒子数
  3. d是每个区域的大小(边长)

通过这个算法:

  1. O(d^2)对于每个粒子,更新受此粒子影响的区域的计数
  2. 遍历所有O(n^2)可能的区域,找到最大值

如此代码所示,我将相关部分复制到这里供您参考:

using namespace std;

int mat [1024 + 3] [1024 + 3]; // Here n is assumed to be 1024

int main ()
{
    int testCases; scanf ("%d", &testCases);

    while ( testCases-- ) {

        Set(mat, 0);

        int d; scanf ("%d", &d); // d is the size of the region
        int k; scanf ("%d", &k); // k is the number of particles

        int x, y, cost;

        for ( int i = 0; i < k; i++ ) {
            scanf ("%d %d %d", &x, &y, &cost); // Read each particle position

            // Update the count of the d^2 region affected by this particle
            for ( int j = max (0, x - d); j <= min (x + d, 1024); j++ ) {
                for ( int k = max (0, y - d); k <= min (y + d, 1024); k++ ) mat [j] [k] += cost;
            }
        }

        int resX, resY, maxi = -1;

        // Find the maximum count over all regions
        for ( int i = 0; i < 1025; i++ ) {
            for ( int j = 0; j < 1025; j++ ) {
                if ( maxi < mat [i] [j] ) {
                    maxi = mat [i] [j];
                    resX = i;
                    resY = j;
                }
            }
        }

        printf ("%d %d %d\n", resX, resY, maxi);

    }
    return 0;
}

我已将我的评论放在代码中向您解释。

于 2013-10-23T06:40:12.380 回答
0

将区域划分为 1000x1000 并计算每个(重叠)2x2 中的粒子数。您可以简单地通过规范化 0..1、缩放 0..999 和转换为整数来对它们进行分区。计数可以很容易地存储为整数的二维数组(ushort、uint 或 ulong...mmmm tea)。这相当于宽相碰撞检测中使用的简单二维空间分割。

于 2013-10-23T05:49:27.843 回答