26、N-Body算法的CUDA实现与优化

N-Body算法的CUDA实现与优化

1. 引言

在处理N个物体的系统时，每个物体的位置为$x_i$，速度为$v_i$（$1 \leq i \leq N$），物体$i$受到物体$j$的力向量$f_{ij}$由以下公式给出：
$f_{ij} = G \frac{m_i m_j}{d_{ij}^2} \cdot \frac{\vec{d} {ij}}{d {ij}}$
其中，$m_i$和$m_j$分别是物体$i$和物体$j$的质量，$\vec{d} {ij}$是从物体$i$到物体$j$的差异向量，$G$是引力常数。由于除法溢出问题，当$d {ij}$的模很小时，该表达式会发散。为了补偿这一点，通常会应用一个软化因子，用于模拟两个Plummer质量（表现得像球形星系的质量）之间的相互作用。对于软化因子$\epsilon$，得到的表达式为：
$f_{ij} = G \frac{m_i m_j \vec{d} {ij}}{(d {ij}^2 + \epsilon^2)^{\frac{3}{2}}}$
物体$i$由于与其他$N - 1$个物体的相互作用而受到的总力$F_i$，是通过对所有相互作用求和得到的：
$F_i = \sum_{j=1,j\neq i}^{N} f_{ij} = \sum_{j=1,j\neq i}^{N} G \frac{m_i m_j \vec{d} {ij}}{(d {ij}^2 + \epsilon^2)^{\frac{3}{2}}}$
为了更新每个物体的位置和速度，作用在物体$i$上的力（加速度）为$a_i = \frac{F_i}{m_i}$，因此可以