Help with Elastic Collision

Hi! I am trying to write code to model elastic collisions. Right now the position vs time graphs are completely off (it doesn't even start at t=0). Second it appears to be an inelastic collision right now. Third, there seems to be in acurrate boundary conditions because the object on the right, the red, moves past x=10. I think that is because it is drawn not from a center of a circle but from the left most point of the 4-point "rectangle" with curvature (aka a circle). Any tips?

I have also tried plot (t, xx1) and (t, yy1) , etc within the loop but that also doesn't work.

Interestingly if you take out the line:

    hold on
    plot (xx1, yy1, '-r');
    plot (xx2, yy2, '-b');
    hold off

it has much trouble ploting anything for the Position vs Time graphs.

       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% Modeling Elastic Collisions  %%%
%%% Final Project Matlab 225     %%%
%%% K. Gifford August 2020        %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

t=0;
dt=0.5;

x1=1;
y1=0;
v1x=0;
v1y=0;
m1=1;
% centers = [x1, y1];

x2=-8;
y2=0;
v2x=2;
v2y=0;
m2=1;
% centers2 = [x2, y2];

% t=0, plot position of particles
subplot (3, 1, 1)
plot ([-10, 10], [-10, 10], 'w')
% obj1 = viscircles(centers,1)
% obj2 = viscircles (centers2, 1)
r1 = rectangle ('Position', [x1,y1,1,1], 'FaceColor', 'r', 'Curvature', [1,1]);
r2 = rectangle ('Position', [x2,y2,1,1], 'FaceColor', 'b', 'Curvature', [1,1]);
axis ([-10 10 -10 10]);

subplot (3, 1, 2)
plot (t, x1)
plot (t, x2)
xlabel ('Time')
ylabel ('Position-X')

subplot (3, 1, 3)
plot  (t, y1)
plot  (t, y2)
xlabel ('Time')
ylabel ('Position-Y')


xx1 = [ ];
xx2 = [ ];
yy1 = [ ];
yy2 = [ ];

% Start moving :)
while t < 100
    if 1 < sqrt((x2-x1)^2+(y2-y1)^2)

        x1 =  x1 + (v1x * dt);
        x2 =  x2 + (v2x * dt);  

        y1 =  y1 + (v1y * dt);
        y2 = y2 + (v2y * dt);

        xx1 = [xx1, x1];
        xx2 = [xx2, x2];
        yy1 = [yy1, y1];
        yy2 = [yy2, y2];
       
    else 
        v1x = ((m1-m2)/(m1+m2))*v1x + ((2*m2)/(m1+m2)) * v2x;
        v1y = ((m1-m2)/(m1+m2))*v1y + ((2*m2)/(m1+m2)) * v2y;
        
        v2x = ((m2-m1)/(m1+m2))*v2x + ((2*m1)/(m1+m2)) * v1x;
        v2y = ((m2-m1)/(m1+m2))*v2y + ((2*m1)/(m1+m2)) * v1y;

        x1 =  x1 + (v1x * dt);
        x2 =  x2 + (v2x * dt);  
    
        y1 =  y1 + (v1y * dt);
        y2 =  y2 + (v2y * dt);
        
        
        xx1 = [xx1, x1];
        xx2 = [xx2, x2];
        yy1 = [yy1, y1];
        yy2 = [yy2, y2];
    end

    % At Boundary - Turn Direction
            if x1 ==10
                v1x= -v1x;
            end

            if x2 == 10
                v2x=-v2x;
            end

            if x1 == -10
                v1x=abs(v1x);
            end

            if x2 == -10
                v2x=abs(v2x);
            end

            if y1 == 10
                v1y=-v1y;
            end

            if y2 == 10
                v2y=-v2y;
            end

            if y1 == -10
                v1y=abs(v1y);
            end

            if y2 == -10
                v2y=abs(v2y);
            end
 
    %update the position of the rectangles 
    set(r1,'Position', [x1,y1,1,1]);
    set(r2,'Position', [x2,y2,1,1]);
    
    subplot (3, 1, 2)
    plot (t, x1)
    plot (t, x2)
    xlabel ('Time')
    ylabel ('Position-X')
    
    subplot (3, 1, 3)
    plot (t, y1)
    plot (t, y2)
    xlabel ('Time')
    ylabel ('Position-Y')
    

    hold on
    plot (xx1, yy1, '-r');
    plot (xx2, yy2, '-b');
    hold off

    drawnow

      t= t+ dt;
 end

% 
   

Answer 1

[Alan Stevens]

You need to be careful about pre and post collision velocities, and the boundary checks:

t=0;
dt=0.5;
x1=1;
y1=0;
v1x=0;
v1y=0;
m1=1;
% centers = [x1, y1];
x2=-8;
y2=0;
v2x=2;
v2y=0;
m2=1;
% centers2 = [x2, y2];
% t=0, plot position of particles
subplot (3, 1, 1)
plot ([-10, 10], [-10, 10], 'w')
% obj1 = viscircles(centers,1)
% obj2 = viscircles (centers2, 1)
r1 = rectangle ('Position', [x1,y1,1,1], 'FaceColor', 'r', 'Curvature', [1,1]);
r2 = rectangle ('Position', [x2,y2,1,1], 'FaceColor', 'b', 'Curvature', [1,1]);
axis ([-10 10 -10 10]);
subplot (3, 1, 2)
plot (t, x1)
plot (t, x2)
xlabel ('Time')
ylabel ('Position-X')
subplot (3, 1, 3)
plot  (t, y1)
plot  (t, y2)
xlabel ('Time')
ylabel ('Position-Y')
xx1 = [ ];
xx2 = [ ];
yy1 = [ ];
yy2 = [ ];

% Start moving :)
while t < 100
    % Update "old" values
    v1xold = v1x; v2xold = v2x; v1yold = v1y; v2yold = v2y;
    
    if 1 < sqrt((x2-x1)^2+(y2-y1)^2)
        x1 =  x1 + (v1x * dt);
        x2 =  x2 + (v2x * dt);  
        y1 =  y1 + (v1y * dt);
        y2 = y2 + (v2y * dt);
        xx1 = [xx1, x1];
        xx2 = [xx2, x2];
        yy1 = [yy1, y1];
        yy2 = [yy2, y2];
       
    else   % Use "old" values on RHS
        v1x = ((m1-m2)/(m1+m2))*v1xold + ((2*m2)/(m1+m2)) * v2xold; 
        v1y = ((m1-m2)/(m1+m2))*v1yold + ((2*m2)/(m1+m2)) * v2yold;
        
        v2x = ((m2-m1)/(m1+m2))*v2xold + ((2*m1)/(m1+m2)) * v1xold;
        v2y = ((m2-m1)/(m1+m2))*v2yold + ((2*m1)/(m1+m2)) * v1yold;
        x1 =  x1 + (v1x * dt);
        x2 =  x2 + (v2x * dt);  
    
        y1 =  y1 + (v1y * dt);
        y2 =  y2 + (v2y * dt);
        
        
        xx1 = [xx1, x1];
        xx2 = [xx2, x2];
        yy1 = [yy1, y1];
        yy2 = [yy2, y2];
    end
    % At Boundary - Turn Direction
            if x1+dt*v1x>9
                v1x= -v1x;
            end
            if x2+dt*v2x>9
                v2x=-v2x;
            end
            if x1+dt*v1x<-10
                v1x=abs(v1x);
            end
            if x2+dt*v2x<-10
                v2x=abs(v2x);
            end
            if y1+dt*v1y>9
                v1y=-v1y;
            end
            if y2+dt*v2y>9
                v2y=-v2y;
            end
            if y1+dt*v1y<-10
                v1y=abs(v1y);
            end
            if y2+dt*v2y<-10
                v2y=abs(v2y);
            end
            
 
    %update the position of the rectangles 
    set(r1,'Position', [x1,y1,1,1]);
    set(r2,'Position', [x2,y2,1,1]);
    
    subplot (3, 1, 2)
    plot (t, x1,'ro',t,x2,'bs')
    axis([0 100 -10 10])
    xlabel ('Time')
    ylabel ('Position-X')
    
    subplot (3, 1, 3)
    plot (t, y1,'ro',t,y2,'bs')
    axis([0 100 -10 10])
    xlabel ('Time')
    ylabel ('Position-Y')

    drawnow
      t= t+ dt;
 end