close all; clear; clc

pendulums = 100;

X = zeros(2 * pendulums, 1);
X(1:2:end) = linspace(-2.9, 2.9, pendulums);

A = diag(ones(1, 2 * pendulums - 1), 1);
A(2:2:end, :) = 0;

B = diag(ones(1, 2 * pendulums - 1), -1);
B = B(:, 1:2:end);

dt = 0.001;
t = 0:dt:10;

g = 9.8;
l = 1;

plot_data = zeros(pendulums, length(t));

tic
for i = 1:length(t)
    X = X + dt * (A * X - (g / l) * B * sin(X(1:2:end)));
    plot_data(:, i) = X(1:2:end);
end
toc

figure
plot(t, plot_data)