% Oscillations of a discrete number of masses.
% Normal modes
om0=1;
N=5;om(1:N)=0;
% Define the normal mode frequencies
for i=1:N
    om(i)=2*om0*sin(i*pi/2/(N+1));
end

nstep=1000;dt=0.1;
t=0:dt:nstep*dt;

x(1:N)=1:N;
fi0(1:N)=0;
y0(1:N,1:N)=0;
A(1:N)=[0.8 0.0 -0.1 0 0.08];    % amplitudes of each mode

sum=0*y0(:,1);
for i=1:N             % initial conditions - modes
y0(1:N,i)=A(i)*sin(i*x*pi/(N+1));
sum=sum+y0(:,i);
end
%plot(x,y0(:,i),'-o');axis ([0 N+1 -2 2])
plot(x,sum,'-o', 'MarkerSize',20);axis ([0 N+1 -2 2])
    line([0 x(1)],[0 sum(1)]);
    line([x(N) N+1],[sum(N) 0]);
%end
%%
% Simulate oscillation in time
y(1:N,1:N)=y0;
for j=1:nstep
    sum=0*y(:,1);
    for i=1:5   % sum over modes
        t=j*dt;
        y(:,i)=y0(:,i)*cos(om(i)*t+fi0(i));
        sum=sum+y(:,i);
    end
    plot(x,sum,'-o','MarkerSize',20);axis ([0 N+1 -3 3]);
    line([0 x(1)],[0 sum(1)]);
    line([x(N) N+1],[sum(N) 0]);
    F = getframe;
end