1

我正在为 C++ 开发一个粒子引力相互作用模拟程序，我几乎可以让它工作，但我有一个愚蠢的无法解释的错误让我发疯。

该程序有 2048 个粒子，并以 200 个时间步长模拟它们；在每个状态下，它会根据所有其他粒子的合力更新每个粒子的新位置、速度和加速度。基本上，这是我编写的递归函数，它遍历 quad_tree 并更新所有粒子值：

    void iterateThroughQuadTree(double alpha, double delta_t)
    {
        std::cout << particles_count << std::endl;
        if (particles_count > 10)
        {
            //std::cout << particles_count << std::endl;
            children[0]->iterateThroughQuadTree(alpha,delta_t);
            children[1]->iterateThroughQuadTree(alpha,delta_t);
            children[2]->iterateThroughQuadTree(alpha,delta_t);
            children[3]->iterateThroughQuadTree(alpha,delta_t);
        }
        if (particles_count > 0)
        {

          std::cout << *(particles_x) << " " << *(particles_y) << " " << *(particles_vx) << " " << *(particles_vy) << std::endl;
          double acc_magnitude = 0, distance = 0, unit_vector_x = 0, unit_vector_y = 0, x_displacement = 0, y_displacement = 0, theta = 0, x_acc = 0, y_acc = 0;
          for (size_t i = 0; i < particles_count; i++) {
              double temp_x_acc = 0;
              double temp_y_acc = 0;
              for (size_t k = 0; k < particles_count; k++) {
                  if (i != k) {
                      x_displacement = particles_x[k] - particles_x[i];
                      y_displacement = particles_y[k] - particles_y[i];
                      distance = sqrt(pow(x_displacement,2) + pow(y_displacement,2));
                      unit_vector_x = x_displacement / distance;
                      unit_vector_y = y_displacement / distance;
                      acc_magnitude = particles_m[k] / pow(distance,2);
//                    //theta = atan(y_displacement/x_displacement);
                      x_acc = acc_magnitude * unit_vector_x;
                      y_acc = acc_magnitude * unit_vector_y;
//                    //x_acc = (x_displacement);
//                    //y_acc = (y_displacement);
                      temp_x_acc = temp_x_acc + x_acc;
                    //temp_y_acc = temp_y_acc + y_acc;
                  }
              }
              //std::cout << temp_x_acc << " " << temp_y_acc << std::endl;
              particles_x[i] += ((particles_vx[i]*delta_t) + alpha*temp_x_acc*pow(delta_t,2));
              particles_y[i] += ((particles_vy[i]*delta_t) + alpha*temp_y_acc*pow(delta_t,2));
              particles_vx[i] += 2*alpha*temp_x_acc*delta_t;
              particles_vy[i] += + 2*alpha*temp_y_acc*delta_t;
//            if (i < 1) {
//                //std::cout << *(particles_x + i) << " ";
//            }
          }
        }
    }

如果你注意两个函数中间的以下两行：

temp_x_acc = temp_x_acc + x_acc;
//temp_y_acc = temp_y_acc + y_acc;

错误是，如果其中一行被注释掉（不管哪一行），程序就可以正常工作。但是，当它们都未注释时，程序会在几次迭代后崩溃，从而导致分段错误。我使用了 ValGrind Memcheck，它说有堆栈溢出。我不明白为什么这两行不能一起出现在程序中；一个简单地将所有 x 分量加速度相加，一个将 y 分量加速度相加。当我有两个 += 运算符时，是否存在某种类型的内存问题？我每次都尝试重新初始化变量，但我仍然得到相同的结果。如果您对为什么会发生这种情况有任何想法，请您向我解释一下吗？谢谢，我很感激你的帮助。

如果需要，我会发布程序的附加代码，问题是程序很大，所以我不想发布超出需要的内容。

4

1 回答 1

0

没有看过你所有的代码，我有一个可能成立的理论。在您未显示的代码中的某处，您肯定必须修改particles_count 的值，它是一个全局变量。除非particles_count 低于10，否则您的递归函数不可能终止。

是否有可能当 temp_x_acc 和 temp_y_acc 都设置时，您的particles_count 值会以某种方式上升，导致循环的更多迭代，导致堆栈溢出？作为帮助验证的健全性检查，我建议您在函数顶部添加一个静态计数器，如下所示：

static int iteration=0;
iteration++;
cerr<<"Iteration: "<<iteration<<", "<<Particles count: "<<particles_count<<endl;

请注意 cerr 的使用，以便在崩溃之前获得所有输出。cout 崩溃时可能不会打印整个故事。

通常对于递归代码，像particles_count 这样的重要变量应该作为函数参数传递。我认为particles_count 的值可能会以您无法预料的方式发生变化。

于 2013-11-15T01:33:34.590 回答