在 Julia 中模拟粒子的碰撞答案

【问题标题】：simulating the collision of particles in Julia在 Julia 中模拟粒子的碰撞
【发布时间】：2021-07-14 14:38:53
【问题描述】：

我想模拟一个盒子内粒子的碰撞。
更具体地说，我想创建一个函数（我们称之为collision!），在每次交互后更新粒子速度，如图所示。

我将粒子（半径等于 1）定义如下：

mutable struct Particle
    pos :: Vector{Float64}
    vel :: Vector{Float64}
end

p = Particle( rand(2) , rand(2) )
# example for the position
p.pos
> 2-element Vector{Float64}:
  0.49339012018408135
  0.11441734325871078

对于碰撞

function collision!(p1::Particle, p2::Particle)
     
    # ... #
    
    return nothing
end

主要思想是当两个粒子碰撞时，它们“交换”平行于粒子中心的速度矢量（矢量 n hat）。
为了做到这一点，需要将速度向量转换为碰撞法线（n hat）的正交基础。

然后交换平行分量，在原来的基础上旋转回来。

我认为我的数学是正确的，但我不确定如何在代码中实现它

【问题讨论】：

标签： julia linear-algebra collision particles

【解决方案1】：

请注意，我根本没有检查过数学，您提供的 2d 案例的一种实现可能是：

struct Particle
    pos :: Vector{Float64}
    vel :: Vector{Float64}
end

p1 = Particle( rand(2) , rand(2) )
p2 = Particle( rand(2) , rand(2) )

function collision!(p1::Particle, p2::Particle)
    # Find collision vector
    n = p1.pos - p2.pos
    # Normalize it, since you want an orthonormal basis
    n ./= sqrt(n[1]^2 + n[2]^2)
    # Construct M
    M = [n[1] n[2]; -n[2] n[1]]
    # Find transformed velocity vectors
    v1ₙ = M*p1.vel
    v2ₙ = M*p2.vel
    # Swap first component (or should it be second? Depends on how M was constructed)
    v1ₙ[1], v2ₙ[1] = v2ₙ[1], v1ₙ[1]
    # Calculate and store new velocity vectors
    p1.vel .= M'*v1ₙ
    p2.vel .= M'*v2ₙ
    return nothing
end

几点：

您不需要mutable struct；只是一个普通的struct 就可以正常工作，因为Vector 本身是可变的
这个实现有很多多余的分配，如果你可以就地工作或者更可行地在堆栈上工作（例如，使用某种静态数组而不是基本数组），你可以避免这些分配作为位置和速度矢量的基础）。如果您只是创建另一个结构（例如“CollisionEvent”）来保存 M、n、v1n 和 v2n 的预分配缓冲区，并将其传递给 collision! 函数，那么就地实际上可能不会太难。
虽然我没有深入了解，但也许可以在像 https://github.com/JuliaMolSim/Molly.jl 这样的分子动力学包中找到此类碰撞的有用参考实现

【讨论】：