如何在 MATLAB 中矢量化随机游走仿真

Question

我正在 MATLAB 中重写一个 Monte Carlo 仿真模型，重点是可读性。该模型涉及许多粒子，表示为 (x,y,z)，它们在具有特定终止概率的一小组状态上随机游走。与输出相关的信息是在给定状态终止的粒子数。

模拟需要足够的粒子，单独为每个粒子运行它的成本高得令人望而却步。矢量化似乎是从 MATLAB 中获得性能的方法，但是有没有任何惯用的方法可以在 MATLAB 中创建此仿真的矢量化版本？

我正在用头撞墙来完成这个 - 我什至尝试创建一个 (nStates x nParticles) 矩阵来表示每个粒子状态组合，但是这种方法在可读性方面很快就会失控，因为粒子从状态反弹相互独立地陈述。我应该硬着头皮改用更适合这个的语言吗？

Answer 1

只需像平常一样编写代码即可。几乎所有 Matlab 函数都可以接受并返回矢量化输入。例如，模拟 N 个粒子在一维中的布朗运动

position = zeros([N 1]); %start at origin
sigma = sqrt(D * dt); %D is diffusion coefficient, dt is time step
for j = 1:numSteps
    position = position + sigma*randn(size(position));
end

如果您想为每个位置使用不同的 sigma，则可以将 sigma 设为与位置大小相同的向量，并使用“点时间”表示法来指示逐个元素操作

position = position + sigma.*randn(size(position));

如果散射是位置和一些随机元素的任意函数，您只需编写一个矢量化函数，例如

function newstep = step(position)
%diffusion in a overdamped harmonic potential
newstep = -dt*k*position + D*randn(size(position));

for j = 1:numsteps; position = position + step(position);

等等