2012-03-25 143 views
2

可以说我想模拟一个粒子状态,它可以是正常的(0)或兴奋的(1)在给定的帧。粒子处于激发状态f%的时间。如果粒子处于激发态,则它持续〜L帧(具有泊松分布)。我想模拟N个时间点的状态。因此,输入例如:matlab:praticle状态模拟

N = 1000; 
f = 0.3; 
L = 5; 

,结果会像

state(1:N) = [0 0 1 1 1 1 0 0 0 0 1 1 1 0 0 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 ... and so on] 

与总和(州)/ N接近0.3

如何做到这一点? 谢谢!

+0

粒子翻转状态的概率是多少? – 2012-03-25 08:39:53

+0

我真的不知道你是什么意思。我真正想要做的是模拟具有两种不同扩散系数的粒子的扩散行为,并且在一个状态或另一个状态中定义了更快和更慢的分量(f)和某种类型的寿命的一部分。我想首先模拟状态(在本例中是两个,但可能更多),然后根据状态(更快或更慢)模拟位移和坐标。我不知道这是否是最好的方法,但这是我心中的第一个方法:) – Art 2012-03-25 08:53:09

+0

@ NoamN.Kremen当f = 0.3,状态1的长度为5.状态0的长度为on平均值应该是17(5/0.3),所以翻转从0到1的变化是0.06。编辑:不确定这个陈述是否完全正确。 – Bernhard 2012-03-25 10:00:00

回答

2
%% parameters 
f = 0.3; % probability of state 1 
L1 = 5; % average time in state 1 
N = 1e4; 
s0 = 1; % init. state 
%% run simulation 
L0 = L1 * (1/f - 1); % average time state 0 lasts 
p01 = 1/L0; % probability to switch from 0 to 1 
p10 = 1/L1; % probability to switch from 1 to 0 
p00 = 1 - p01; 
p11 = 1 - p10; 
sm = [p00, p01; p10, p11]; % build stochastic matrix (state machine) 
bins = [0, 1]; % possible states 
states = zeros(N, 1); 
assert(all(sum(sm, 2) == 1), 'not a stochastic matrix'); 
smc = cumsum(sm, 2); % cummulative matrix 
xi = find(bins == s0); 
for k = 1 : N 
    yi = find(smc(xi, :) > rand, 1, 'first'); 
    states(k) = bins(yi); 
    xi = yi; 
end 
%% check result 
ds = [states(1); diff(states)]; 
idx_begin = find(ds == 1 & states == 1); 
idx_end = find(ds == -1 & states == 0); 
if idx_end(end) < idx_begin(end) 
    idx_end = [idx_end; N + 1]; 
end 
df = idx_end - idx_begin; 
fprintf('prob(state = 1) = %g; avg. time(state = 1) = %g\n', sum(states)/N, mean(df)); 
+0

这真棒,谢谢很多Serg! – Art 2012-03-26 03:28:58

2

激发态的平均长度为5.正常状态的平均长度,因此应该在12左右才能获得。

策略可以是这样的。

  • 开始在状态0
  • 从与平均L*(1-f)/f
  • 泊松分布绘制的随机数aa零填充状态阵列
  • 从与平均L一个泊松分布绘制的随机数b
  • 用状态数组填充b
  • 重复

另一种选择是考虑在切换概率,其中0-> 1和1-> 0的概率是不相等的条件。