Plotting outputs from for loops in a 3D space

So, I managed to get a simple function to work approximately the way that I wanted, which was to take a point in 3 dimensional space, described by a vectors (x,y,z) and ascribe forces to it, based on its proximity to other points ( in this case, fixed points) described by vectors. The code is as follows.

    %particle  vector [x,y,z] position
pm= [1;1;1];
%particle mass, arbitrary metric
mm1=2.5;
%inital velocity of particle
vm0=[.00001;-.0000003;.000006];
%initial acceleration of particle
am0=[-.00001;.00005;-.000002];
%two forces, respectively defined by proximity to two points [5,5,5] and [-3,-3,-3]
f1= .00000006*(pm - ([5;5;5])).^2;
f2= (.0000006*(1./(pm-[-3;-3;-3]).^2));
%acceleration of particle based on the previous two forces 
am= am0+(f1+f2)./mm1;

 %a for loop that iteratively calculates new positions based on the previous
%values, with respect to the time variable
for t=1:.01:100
pm= ((am*(t.^2))./2)+(vm0*t)+pm
plot(pm)
end

I want to display the progression of my particle, its trajectory, in a 3D plot. How would I go about this?

Also, I would like to be able to generalize this to 'n' particles, all starting at various positions (only interacting with [5,5,5] and [-3,-3,-3], not each other). Does anyone have any thoughts on how that could be achieved elegantly? Again, I'd like to be able to visualize the trajectories of these particles as well.

I appreciate any and all help.

Answer 1

[John Chilleri]

Hello,

This should be somewhat like you wanted:

%particle  vector [x,y,z] position
pm= [1;1;1];
%particle mass, arbitrary metric
mm1=2.5;
%inital velocity of particle
vm0=[.00001;-.0000003;.000006];
%initial acceleration of particle
am0=[-.00001;.00005;-.000002];
%two forces, respectively defined by proximity to two points [5,5,5] and [-3,-3,-3]
f1= .00000006*(pm - ([5;5;5])).^2;
f2= (.0000006*(1./(pm-[-3;-3;-3]).^2));
%acceleration of particle based on the previous two forces
am= am0+(f1+f2)./mm1;
%a for loop that iteratively calculates new positions based on the previous
%values, with respect to the time variable
figure(1)
hold on
%axis([-200 50 0 900 -30 10]) % uncomment this if you dont want the axes to ever change
for t=1:.01:100
pm= ((am*(t.^2))./2)+(vm0*t)+pm;
%if (mod(t,1) == 0)           % uncomment this if you want to plot only 100 equally spaced points
plot3(pm(1),pm(2),pm(3),'*');
drawnow
%end                          % uncomment this if you want to plot only 100 equally spaced points
end

If you want other specifications, such as color, one developing line, size of points, or anything, I can help you with that as well. It's also important to note that this runs somewhat slowly because you're plotting 9901 times. There are ways to make this faster, such as plotting every 10 points which is also simple to implement.

You may also notice that I included three lines to plot only 100 points rather than all 9901 and an axis call that will maintain the view. I would suggest running it without these once as it does an intriguing curve at the beginning which is hard to see once it's fully plotted.

Hope this helps!

[John Chilleri]

A simple solution to plotting n particles would be to compute all of their "pm's" and store those as columns in a matrix (call it pm). Then, use a for loop like the one below to plot them all simultaneously:

% Assume pm is an (9901 x (n*3)) matrix.
figure
hold on
for i = 1:9901
   for j = 1:3:n-2
      plot3(pm(i,j),pm(i,j+1),pm(i,j+2),'*') %specify color if desired
      drawnow
   end
end

An issue you may encounter is the axis rapidly switching to accommodate the separate plots, try to determine axis limits that will contain the entirety of all the particles and you can set it before the loop as was done with your original code.

[Erik Anderson]

Thank you very much! This is precisely what I was looking for.

I'd like to learn how to add a trailing line in the 3D plot, as I think it would make the progression, and the sense of movement in 3 space more apparent.

You are correct, the only issue with the plot is the rapid switching of axes after the first several hundred iterations. Setting axis limits should take care of that.

I will try the technique of storing position data in a (t x n*3) matrix.

Is there any good way of generating a whole set of positions ('pm's), without rewriting my existing code n-times?

I appreciate you leaving the code for the first 100 iterations. It is an interesting curve, and really becomes indistinct in the larger plot. In fact, after the first 80 units of time, the trajectory is fairly stable, and probably not necessary for me to calculate. This wouldn't be as readily apparent without that visual.

Thanks again.

[John Chilleri]

Hello,

I'm glad you found my answer helpful!

I think the simplest way to generate a whole set of positions without rewriting your existing code would be to define pm as a cell.

pm = cell(n,1);
for i = 1:n
   pm{i} = %whatever;
   % Run your code
end

As an important note, when you access pm in the code, you'll need to access the elements with:

pm{i}(element numbers)

% If pm{i} = [3;5;7];
% pm{i}(1) = 3
% pm{i}(2) = 5
% pm{i}(3) = 7

I hope this helps!

[John Chilleri]

Also, checkout addpoints for improved animation!