为这些耦合ode编写一个函数。在这个函数中，您可以使用您能想到的所有常规编程工具和技巧。像这样的东西应该工作:
函数dOISdt = your_chemistry_ode(t,OIS,kP,kD)
O = OIS(1:end-2);
I = OIS(end-1);
S = OIS(结束);
dSdt = kD*I;
dIdt = kP*sum(O(1:end-1)) - kD*I;
dOdt = 0 (size(O(:)))
dOdt(1) = -kP*O(1);不太清楚你想如何对待O的第一个品种
为i1 = 2:数字(O)
dOdt(i1) = kP*(O(i1-1)-kP*O(i1));
结束
dOISdt = [dOdt;dIdt;dSdt];
结束
然后，将其与任何odenn函数集成应该相当简单:
初始条件
我只是模仿一些东西，完全不知道是什么链
% 100种物种向I和S you模型衰减
O0 = [100;rand(99,1)];
I0 = pi;
S0 = 0;
Ois0 = [o0; i0; s0];
kP = 1e6;%反应率，也不知道时间常数
kD = 1e4;
T_span = linspace(0,1,1e6+1);%相同，ms的步长为1 s…
[tOut,OISout] = ode45(@(t,IOS) your_chemistry_ode(t,OIS,kP,kD)，t_span,OIS0);
HTH
