1.
function ode(t, y)
dydt = [y[2] * y[3]; - y[1] * y[3]; 0.51 * y[1] * y[2]]
return dydt
end
t,y, = ode45(ode, [0 5], [0; 1; 1]);
plot(t, y[:, 1], "-o", t, y[:, 2], "-o", t, y[:, 3], "-o"; markersize = 3.5)
legend(["y1", "y2", "y3"])
grid("on")
xlim([0, 5])
gcf().set_facecolor("none")
plt_print("某刚性体的运动", "-dpng")

2.
function ode_SIR(t, y)
global N
β = 0.4 # 感染率
γ = 0.1 # 康复率
S = y[1];
I = y[2];
R = y[3];
dSdt = -β * S * I / N
dIdt = β * S * I / N - γ * I
dRdt = γ * I
return [dSdt; dIdt; dRdt]
end
N = 10000 # 总人口
y0 = [N - 10; 10; 0]; # [S0; I0; R0]
t, y, = ode45(ode_SIR, [0 100], y0)
plot(t, y[:, 1], "-o", t, y[:, 2], "-o", t, y[:, 3], "-o"; linewidth = 0.8, markersize = 1.5)
xlabel("时间")
ylabel("人数")
legend(["易感染人数", "感染者人数", "康复者人数"])
grid("on")
xlim([0, 100])
gcf().set_facecolor("none")
plt_print("感染人数的预测", "-dpng")