% simuloidaan usean massan ja jousen yhdistelmää.

clear;
close all;

x = [0 1.5/6 2/6 3/6 4/6 5/6 1;
    0  0.5    0   0   0   0  0]; % massojen sijainnit x,y-tasossa
% ensimmäinen ja viimeinen piste ovat kiinnityspisteet 'seiniin'
v = zeros( size( x)); % alkunopeudet
a = zeros( size( x)); % alustetaan kiihtyvyysmatriisi sopivan kokoiseksi
F = zeros( size( x)); % ...sama juttu voimille

t = 0; % aika
dt = 0.01; % ajanhetkien erotus

K = 250; % jousivakio, N/m
m = 0.5; % massa, kg

iter = 0;
while iter < 1000,
    % lasketaan voimat ja niiden avulla nopeudet
    for i = 2:length(x)-1,
        F(:,i) = K*(x(:,i-1)-x(:,i))+K*(x(:,i+1)-x(:,i));
        a(:,i) = F(:,i)/m;
        v(:,i) = v(:,i)+dt*a(:,i);                      
    end
    % päivitetään sijainnit nopeuksien avulla
    for i=2:length(x)-1,
        x(:,i) = x(:,i) + v(:,i)*dt;
    end
    
    % piirretään kuva tilanteesta
    figure(1);
    cla;
       
    for i=2:length(x)-1,
        draw_spring( x(:,i-1), x(:,i), 7, 0.03);
        hold on
    end
    draw_spring( x(:,end-1), x(:,end), 7, 0.03);
    for i=2:length(x)-1,
        h = plot( x(1,i), x(2,i), 'ro', 'MarkerSize', 10, ...
            'MarkerFaceColor' ,'r');
    end

    % seinät
    patch( [0 0 -0.05 -0.05], [-1 1 1 -1], 'k');
    patch( [1 1 1.05 1.05], [-1 1 1 -1], 'k');
    

    
    axis( [0 1 -1 1]);
    axis equal;
    disp('please press a key');
    title(['t = ' num2str(t) ' s']);
    pause;
    x
    iter = iter + 1;
end