Python | MATLAB | 粒子动力 | 数学 | 动能 | 矢量 | 势能 | Matplotlib
我们将假设一个由 N 个点粒子组成的系统,索引为 $i=1, \ldots, N$。每个粒子都有一个:
质量
位置 ri=[xi,yi,zi]r_{\boldsymbol{i}}=\left[x_i, y_i, z_i\right]ri=[xi,yi,zi]
速度 vi=[vxi,vyi,vzi]\mathbf{v}_{\boldsymbol{i}}=\left[\mathrm{vx}_i, \mathrm{v} y_i, \mathrm{v} z_i\right]vi=[vxi,vyi,vzi]
根据牛顿万有引力定律(著名的“平方反比定律”),每个粒子都会感受到所有其他粒子的引力。也就是说,每个粒子都会感受到加速度:
ai=G∑j≠imjrj−ri∣rj−ri∣3\mathbf{a}i=G \sum{j \neq i} m_j \frac{\mathbf{r}_j-\mathbf{r}_i}{\left|\mathbf{r}_j-\mathbf{r}_i\right|^3}ai=G∑j=imj∣rj−ri∣3rj−ri
位置和速度使用蛙跳方案进行更新。对于每个时间步 Δt,每个粒子都会受到半步“跳”:
vi=vi+Δt2×ai\mathbf{v}_i=\mathbf{v}_i+\frac{\Delta t}{2} \times \mathbf{a}_ivi=vi+2Δt×ai
上述代码的性能在Python中实际上可以通过向量化来提高。 也就是说,用向量和矩阵运算来表述问题。 它通常可以带来 100 倍的加速,还使代码更具可读性。 缺点是在矩阵内存储中间计算会占用大量内存。 这是计算所有粒子加速度的函数的矢量化版本:
我们的代码计算这些量并跟踪总能量,以确保通过数值方法近似守恒。
Last updated 6 months ago