Hello Jonathan,
I understand that you want to add 150 molecules exactly at 400 seconds and remove 150 molecules exactly at 800 seconds.
This can be done as follows. You are storing all the values of the time steps “t” in an array “t_nr”. This array can be used to access previous time step value. Use an if conditional to check if the current time step is greater than 400 and the previous time step from the array is less than 400. If this is satisfied, then it indicates that the time 400 seconds was just passed. This means that 150 molecules should be added at this instant. So, at this instance check whether the “is_perturbation” flag is still zero. Once this is satisfied add 150 to the current value of “x_nr” that is “x_nr(step_count)”. To prevent adding multiple times once the molecules are added the flag is set to 1.
if t>perturbation_time_400 && t_nr(step_counter-1)<400 && is_perturbation_400==0
x_nr(step_counter)=x_nr(step_counter-1)+150;
is_perturbation_400=1;
end
Similar method can be applied at 800 seconds as well.
if t>perturbation_time_800 && t_nr(step_counter-1)<800 && is_perturbation_800==0
x_nr(step_counter)=x_nr(step_counter-1)-150;
is_perturbation_800=1;
end
One issue with this is that once 150 molecules are added the “x_nr” values become big so thereby reducing the “wait_time” values to very small numbers this can increase the number of times the loop is run significantly. This will increase run time.
In the plot the values on the x axis correspond to the indexes of the plotted arrays. Since no x variable was mentioned in the “plot” the index values are being used and the 5000 is because the size of “x_nr” is nearly 5000. To change this, use the plot function with “t_nr” as the x variable. Using “hold on” makes them plot on the same graph but as y variables and hence the unwanted behaviour. To fix this plot the variables against “t_nr” and “t_ar” instead of along them. The three plots can be changed as follows to obtain the correct graphs.
For first figure:
figure,
plot(t_nr,x_nr,'r','LineWidth',2)
For the next two figures:
figure,
plot(t_ar,x_ar,'r','LineWidth',2)
figure,
plot(t_ar,x_ar,'r','LineWidth',2)
hold on
plot(t_nr,x_nr,'c','LineWidth',2)
legend('x_ar','x_nr')
You can refer to the below documentation to learn more about plot:
Hope this helps.
Thanks,
Charan.
