gusucode.com > 多个小球碰撞有背景颜色支持更换小球个数项目 > Billiards.m

    function[]=Billiards(NumberOfBalls)
% Call: Billiards(Number of Balls), e.g.:
% Billiards(10)
if nargin == 0
    NumberOfBalls = 10;%规定球的个数    重要
end
close all;
hold on;
drawflag=1;
factor=82; %Adjust this factor when scaling  重要
DT=5e-3;
Bound=[-4 4 -2 2];
BallColour=[[1 0 0];[1 0 0.5];[1 0.5 0];[0 1 0];[0 0 1];[1 1 0];[1 1 1];...
[0 0.3 0];[0 0 0];[0.65 0.65 0.65];[0 0.75 0.75];[0.3 0 0.6];[0.95 0.65 0.75];...
[0.5 0.25 0];[0 0.2 0.4];[0.9 0.4 0.7];[0.4 0.2 0.3];[0.65 0.55 0.15];[0.25 0.35 0.25];[0.5 0 0]];
TableColour=[.4 .5 .8];
%Plothandle
axis(Bound);
set(gca,'Color',TableColour,'xcolor',TableColour,'ycolor',...
TableColour,'PlotBoxAspectRatio',[1 abs((Bound(3)-Bound(4))/(Bound(2)-Bound(1))) 1],'xtick',[],'ytick',[])
%===============================================================
r=0.15+0.25*rand(NumberOfBalls,1);
%Radii==========================================================
for count1=1:NumberOfBalls;
    for count2=1:count1,
        rmatrix(count2,count1)=r(count1)+r(count2);
    end;
end;
rmatrix=rmatrix.*triu(abs(-1+eye(size(rmatrix))));
%Mass ==========================================================
massa=pi*r.^2;
X=[(Bound(2)-Bound(1)-2*max(r))*rand(NumberOfBalls,1)+Bound(1)+max(r),...
(Bound(4)-Bound(3)-2*max(r))*rand(NumberOfBalls,1)+Bound(3)+max(r)];
for j=1:NumberOfBalls;
    for i=1:j;
        distmatrix(i,j)=sqrt((X(j,1)-X(i,1))^2+(X(j,2)-X(i,2))^2);
    end;
end;
%Initial Edgedetectionmatrix====================================
Botsmatrix=(distmatrix-rmatrix)+tril(abs(-1+eye(size(distmatrix))))+eye(size(distmatrix));
while find(Botsmatrix<=0);
    X=[(Bound(2)-Bound(1)-2*max(r))*rand(NumberOfBalls,1)+Bound(1)+max(r),...
    (Bound(4)-Bound(3)-2*max(r))*rand(NumberOfBalls,1)+Bound(3)+max(r)];
    for j=1:NumberOfBalls;
        for i=1:j;
            distmatrix(i,j)=sqrt((X(j,1)-X(i,1))^2+(X(j,2)-X(i,2))^2);
        end;
    end;
    Botsmatrix=(distmatrix-rmatrix)+tril(abs(-1+eye(size(distmatrix))))+eye(size(distmatrix));
end
%Initial velocities
V=5*(-1+2*rand(NumberOfBalls,2));
%Plot startingpositions
for k=1:NumberOfBalls;
    h(k) = plot(X(k,1),X(k,2),'o', 'MarkerEdgeColor',BallColour(mod(k-1,length(BallColour))+1,:),...
    'MarkerFaceColor',BallColour(mod(k-1,length(BallColour))+1,:),'MarkerSize',factor*r(k));
end
pause(1e-1);
%Loop
while drawflag==1;
    for k=1:NumberOfBalls;
        delete(h(k));
    end
    %Edgedetecton positive======================================
    for j=1:NumberOfBalls;
        for i=1:2;
            d=X(j,i)+r(j)-Bound(2*i);
            if d>=0
                dt=d/V(j,i);
                X(j,i)=X(j,i)-V(j,i)*dt;
                V(j,i)=-V(j,i);
            end
        end
    end
    %Edgedetecton negative======================================
    for j=1:NumberOfBalls;
        for i=1:2;
            d=X(j,i)-r(j)-Bound(2*i-1);
            if d<=0
                dt=d/V(j,i);
                X(j,i)=X(j,i)-V(j,i)*dt;
                V(j,i)=-V(j,i);
            end
        end
    end
    %Distancematrix============================================
    for j=1:NumberOfBalls;
        for i=1:j;
            distmatrix(i,j)=sqrt((X(j,1)-X(i,1))^2+(X(j,2)-X(i,2))^2);
        end;
    end;
    %Collisiondetectionmatrix=========================================================
    Botsmatrix=(distmatrix-rmatrix)+tril(abs(-1+eye(size(distmatrix))))+eye(size(distmatrix));
    %===========================================================================
    if find(Botsmatrix<0);
        [I,J]=find(Botsmatrix<0);
        for i=1:length(I)
            normdist=normr([X(I(i),1)-X(J(i),1) X(I(i),2)-X(J(i),2)]);
            vaA=(V(I(i),1)*normdist(1)+V(I(i),2)*normdist(2));
            vaB=(V(J(i),1)*normdist(1)+V(J(i),2)*normdist(2));
            dt=abs(r(I(i))+r(J(i))-distmatrix(I(i),J(i)))/(abs(vaA)+abs(vaB));
            X(I(i):J(i),:)=X(I(i):J(i),:)-V(I(i):J(i),:)*dt;
        end
        for i=1:length(I)
            normdist=normr([X(I(i),1)-X(J(i),1) X(I(i),2)-X(J(i),2)]);
            vaA=(V(I(i),1)*normdist(1)+V(I(i),2)*normdist(2));
            vbA=(-V(I(i),1)*normdist(2)+V(I(i),2)*normdist(1));
            vaB=(V(J(i),1)*normdist(1)+V(J(i),2)*normdist(2));
            vbB=(-V(J(i),1)*normdist(2)+V(J(i),2)*normdist(1));
            va_A=vaA+2*(vaB-vaA)/(1+massa(I(i))/massa(J(i)));
            va_B=vaB+2*(vaA-vaB)/(1+massa(J(i))/massa(I(i)));
            V(I(i),1)=va_A*normdist(1)-vbA*normdist(2);
            V(I(i),2)=va_A*normdist(2)+vbA*normdist(1);
            V(J(i),1)=va_B*normdist(1)-vbB*normdist(2);
            V(J(i),2)=va_B*normdist(2)+vbB*normdist(1);
            X(I(i):J(i),:)=X(I(i):J(i),:)+V(I(i):J(i),:)*dt;
        end
    end
    %Propagation=================================================
    X=X+V*DT;
    %Plotting====================================================
    for k=1:NumberOfBalls;
        h(k) = plot(X(k,1),X(k,2),'o', 'MarkerEdgeColor',BallColour(mod(k-1,length(BallColour))+1,:),...
        'MarkerFaceColor',BallColour(mod(k-1,length(BallColour))+1,:),'MarkerSize',factor*r(k));
    end
    drawnow;
end