9

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

模拟需要足够多的粒子,因此为每个粒子单独运行它的成本太高了。向量化似乎是从 MATLAB 中获得性能的方法,但是有没有什么惯用的方法可以在 MATLAB 中创建这个仿真的向量化版本?

我正在努力实现这一点——我什至尝试创建一个 (nStates x nParticles) 矩阵来表示每个粒子状态组合,但是这种方法在可读性方面很快就失控了,因为粒子从状态反弹彼此独立地陈述。我应该硬着头皮改用更适合这个的语言吗?

4

1 回答 1

3

只需像往常一样编写代码。几乎所有的 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);

等等

于 2010-09-24T13:07:46.267 回答