I am trying to do this problem, but I'm having some issues with steps 2 3 and 4. I thought I did it right, but the values stored in the "velocites" vector all share the same value and they shouldn't. Any help would be appreciated

1 view (last 30 days)
InitialVelocity = 2 * rand(1) - 1;
F = 1; tau = 1; N = T/dt; T = 10; dt = 0.01; n = 0;
velocities = zeros(1,N); for i = 1:N n = n + 1; Pc = dt/tau; X = rand; if X < Pc vel = ((n-1)*dt) + F * dt; %Collision occurs else vel = F * dt; %No collision occurs end velocities(i) = vel; end
velocities
  1 Comment
Jan
Jan on 31 Aug 2017
Edited: Jan on 31 Aug 2017
Please use the "{} code" button for posting readable code. Increase the size of the screenshot, it is a pain to read the small characters.

Sign in to comment.

Accepted Answer

Jan
Jan on 31 Aug 2017
Edited: Jan on 31 Aug 2017
Your assumption is wrong: the values stored in the "velocites" vector do not share the same value. About 1 of 100 is different. Simply try it:
all(velocities == velocities(1))
You get FALSE in the vast majority of cases.
But at least most of the elements of velocities contain the same value , because this is the expected result of the code. In 99% of the cases the elements are set to:
vel = F * dt;
Velocity is force times time. Hm. Perhaps you mean:
a = F / m; % With m is the mass and a the acceleration
vel = vel + a * dt

More Answers (1)

John BG
John BG on 5 Sep 2017
Edited: John BG on 6 Sep 2017
Hi Caleb
I have reproduced the steps listed in your question, not to optimise code or speed, but to precisely reproduce the experiment.
1.
start variables
clear all;clc;close all
dt=.01;tau=1;m=1;N=1e3;T=N*dt
F=1; % N = kg*m/s^2
pcol=dt/tau % collision probability
n=[1:1:N]; % we can use it but k comes handy
t=[0:dt:T-dt]; % not used
2.
Single run
v=zeros(1,N);
v(1)=randi([-100 100],1,1)/100; % step 1)
for k=2:1:N
X=randi([0 1],1,1);
if X v(k)=F*dt; end
if ~X v(k)=v(k-1)+F*dt; end
end
3.
100 runs
N2=100;
Vn=zeros(N2,N);
for s=1:1:N2
v=zeros(1,N);
v(1)=randi([-100 100],1,1)/100;
for k=2:1:N
X=randi([0 1],1,1);
if X v(k)=F*dt; end
if ~X v(k)=v(k-1)+F*dt; end
end
Vn(s,:)=v;
end
.
zoomed image, full runs are N=1e3 samples long.
.
.
if you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?
To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link
thanks in advance
John BG
.
For the readers
.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!