Parallel-izing, speeding up simulation
2 次查看(过去 30 天)
显示 更早的评论
Hello,
I am trying to simulate a stochastic equation and empirically determine it's pdf. I have the following code,
NN=100;
ind = @(x) (x>0);
x_vect=zeros(NN,1001);
for j=1:NN
ts = 0.01;
x0=0;
x_vect(j,1)=x0;
if (mod(j,1000)==0)
rng('shuffle')
end
for i=1:100000
w = rand(1,1);
if ~(w < ts)
if (x0-2*ind(x0)*ts>0)
x_vect(j,i+1)=x0-2*ind(x0)*ts;
else
x_vect(j,i+1)=0;
end
else
x_vect(j,i+1)=x0+1;
end
x0=x_vect(j,i+1);
end
end
Each sample path should last for 1000 seconds (so the inner loop runs till 10000) and I need to do a total of NN=100000 sample paths. Currently this would take me about 7 hours to run the simulation. Does anyone see an obvious way to speed up the simulation and/or parallelize this code?
0 个评论
采纳的回答
Peter O
2017-2-21
Taking a quick look. There's probably more you can do.
1. Eliminating the anonymous function call on the conditional inside the inner loop should save you some time since it has a healthy amount of overhead (unless Mathworks did some magic in 2016b's JIT). I'd see what inlining it right under the ~(w < ts) does for speed.
2. You're going against MATLAB's memory indexing, which is column major order. Since your inner loop i is striding columns you're paging RAM like nobody's business. Can you rearrange the data so that the inner loop is the one looping on the rows? You want to go down rows on your inner loop. I think you can just do a transpose on how you're storing your data in this case.
3. If you're only altering the seed every 1000 rows, why not grab your rand() call then rather than call it inside the inner loop 100k times? If you've got the RAM to spare. rand is fairly quick, so I'm not positive this will outweigh the memory hit. Depends on the size of the vector, I'd wager.
4. It doesn't seem that your j's are interacting with previous ones (which makes sense since you're doing this stochastically). So, MATLAB should no problem treating this as a parfor loop. The outer loop should be the parfor one.
更多回答(1 个)
Jan
2017-2-21
编辑:Jan
2017-2-21
- x_vect is pre-allocated to the final size.
- anonymous function ind is replaced by direct comparison.
- x0 - 2 * (x0 > 0) * ts is calculated once only.
- x_vect accessed in columnwise order.
NN = 100;
ni = 100000;
x_vect = zeros(ni+1, NN);
for j=1:NN
ts = 0.01;
x0 = 0;
x_vect(j,1) = x0;
if mod(j,1000) == 0
rng('shuffle')
end
for i = 2:ni+1
w = rand;
if w >= ts
v = x0 - 2 * (x0 > 0) * ts;
if v > 0
x_vect(i, j) = v;
end
else
x_vect(i, j) = x0 + 1;
end
x0 = x_vect(i, j);
end
end
This give a boost from 35.1 sec to 0.60 sec on my machine. Now start experiments with parfor.
2 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Parallel for-Loops (parfor) 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!