Putting multiple ODE equations into one script

16 次查看（过去 30 天）
Zafina
Zafina 2020-10-4
评论： Alan Stevens 2020-10-5
I am trying to put multiple ODEs into one script. They model a CSTR in series, where each CSTR is it's own ODE. The first equation (dvdt) is supposed to model the water flow in the module and the second equation (dbdt) models the solute flow. There are 44 CSTRs in the module (ps.N)I , so I used a for loop.
function dkdt = ode_CSTR_series(t,x, ps)
dkdt = zeros(2,ps.N);
v = x(1);
b = x(2);
dvdt(1) = (ps.v_in - v(1)) - (ps.k .* ps.area .* ps.delta_P);
dbdt(1) = (ps.b_in - b(1))/ps.tau;
for i = 2:ps.N
dvdt(i) = (v(i-1) - v(i)) - (ps.k .* ps.area .* ps.delta_P);
dbdt(i) = (b(i-1) - b(i))/ps.tau;
end
dkdt = [dvdt; dbdt];
end
I am trying to get two solutions using ode15s. This is my separate script.
%Initial conditions and time span
initial = zeros(2,ps.N);
initial(1,1) = 0; % mol/s
initial(2,1) = 0;
tspan = t_water0;
[t,b] = ode15s(@ode_CSTR_series, tspan, initial, [], ps);
I get an error that says "Index exceeds the number of array elements (1)". As I understand, each equation should be kept in one row.
How do I solve this? Thanks.

请先登录，再回答此问题。

回答（1 个）

Alan Stevens
Alan Stevens 2020-10-4
编辑：Alan Stevens 2020-10-4
In
v = x(1);
b = x(2);
dvdt(1) = (ps.v_in - v(1)) - (ps.k .* ps.area .* ps.delta_P);
dbdt(1) = (ps.b_in - b(1))/ps.tau;
v and b are just scalars, so the latter two equations shoud be
dvdt(1) = (ps.v_in - v) - (ps.k .* ps.area .* ps.delta_P);
dbdt(1) = (ps.b_in - b)/ps.tau;
This explains why you get the error message, but doesn't completely solve your problem as you need v and b to be vectors, but they must be updated before being used in your for loop. You probably need to put the loop around the calling function ([t,b] = ode15s(...).
  4 个评论
显示 2更早的评论
Zafina
Zafina 2020-10-5
In the solution, k is not within any equation of the for loop, so how does it know to keep running until ps.N?
Alan Stevens
Alan Stevens 2020-10-5
When k is 1 all the calculations within the loop are carried out. Then k increments to 2 and all the calculations within the loop are carried out again. This repeats until k exceeds ps.N. This happens even if k isn't used explicitly within the loop. Incidentally, I just listed a skeleton code. I assume you will want to do other things within the loop (e.g. saving or plotting values of v an b).

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Mathematics 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by