C++ 2D N体引力模拟周期性性能波动问题排查求助及根因确认

阿华AIGC实验室

2026-4-2

最近我用C++开发一个2D N体引力模拟程序时，遇到了一个非常头疼的周期性性能波动问题——程序的更新时间会在“快状态”和“慢状态”之间规律性切换，导致帧率忽快忽慢。折腾了好一阵终于找到根因，这里把整个排查过程和结果分享出来，希望能帮到遇到类似问题的朋友。

问题现象

我的模拟需要计算5184个天体的引力相互作用，不管用std::execution::par_unseq并行执行还是std::execution::seq串行执行，都会出现明显的双峰性能分布：

快状态：单次更新耗时约60ms
慢状态：单次更新耗时在90-150ms之间

更糟的是，程序有70%以上的时间都处于慢状态，整体帧率表现很差。性能分析图能清晰看到这种周期性的跳变（曲线会在两个时间区间反复横跳）。

核心计算逻辑

下面是我从CpuAcceleratorPar.cpp里简化出来的核心计算代码，天体数据存在std::vector<Body>里，Body是包含位置、速度、质量的简单结构体：

#include <vector>
#include <ranges>
#include <execution>

template<typename T>
struct Vector2 {
    T x, y;
    T LengthSquared() const { return x*x + y*y; }
    Vector2 operator-(const Vector2& other) const {
        return {x - other.x, y - other.y};
    }
    Vector2& operator+=(const Vector2& other) {
        x += other.x; y += other.y;
        return *this;
    }
    Vector2 operator*(T scalar) const {
        return {x*scalar, y*scalar};
    }
};

struct Body {
    Vector2<double> Position;
    Vector2<double> Velocity;
    double Mass; // 简化为标量质量
};

void CpuAcceleratorPar::Update(std::vector<Body>& bodies, double dt) {
    const std::size_t n = bodies.size();
    if (n == 0) return;

    std::vector<std::size_t> idx(n);
    std::ranges::iota(idx, std::size_t{ 0 });

    // Pass 1: 计算每个天体的加速度
    std::vector<Vector2<double>> acc(n, Vector2<double>{});
    std::for_each(std::execution::par_unseq, idx.begin(), idx.end(), [&](std::size_t i) {
        Vector2<double> ai{};
        const auto pi = bodies[i].Position;
        for (std::size_t j = 0; j < n; ++j) {
            if (j == i) continue;
            const Vector2<double> r = bodies[j].Position - pi;
            const double dist2 = r.LengthSquared() + 1e-10; // 软化项避免除零
            const double inv = 1.0 / std::sqrt(dist2);
            const double inv3 = inv * inv * inv;
            ai += r * (6.67430e-11 * bodies[j].Mass * inv3);
        }
        acc[i] = ai;
    });

    // Pass 2: 积分更新速度和位置
    std::for_each(std::execution::par_unseq, idx.begin(), idx.end(), [&](std::size_t i) {
        bodies[i].Velocity += acc[i] * dt;
        bodies[i].Position += bodies[i].Velocity * dt;
    });
}

已尝试的排查方向

一开始我完全把注意力放在代码逻辑上，试了几个常见的优化方向：

快速平方根近似：怀疑1.0 / std::sqrt(dist2)会产生非归一化浮点数拖慢性能，特意换成了经典的快速平方根实现，但完全没效果——说明性能瓶颈根本不在平方根计算上。
排除多线程问题：把并行执行换成串行的std::execution::seq，结果波动依然存在，只是整体计算时间变长了，这说明问题不是多线程调度导致的，而是更深层的系统或硬件问题。