clear all; close all; clc;

% parameter for the run
N1 = 100;
om0 = 5;
omw = 2;
K = 2.5;

% Kc
display(2*sqrt(2)*omw/sqrt(pi))

% random frequencies
om = randn(N1,1)*omw + om0;
om = sort(om);

% solve model
opts = odeset('RelTol',1e-5,'AbsTol',1e-5,'MaxStep',0.1);
[t,th1] = ode45(@(t,th) Kuramoto_ode(t,th,om,K),[0 20],zeros(N1,1),opts);

% plot solutions
plot(t,sin(th1)+(1:N1)*2,'LineWidth',5)
figure;
plot(t,sum(sin(th1),2)/N1,'LineWidth',5)
xlim([5 20])
