Error using odearguments (line 21) When the first argument to ode23s is a function handle, the tspan argument must have at least two elements.

13 次查看（过去 30 天）
Dursman Mchabe
Dursman Mchabe 2021-3-26
评论： Bjorn Gustavsson 2021-3-26
Hi everyone,
I am trying to solve a set of odes using ode23s, yet I get the error message that states:
Error using odearguments (line 21)
When the first argument to ode23s is a function handle, the tspan argument must have at least two elements.
Error in ode23s (line 121)
= odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in funcKinetics1 (line 3)
[~, y] = ode23s(@(t,y) funcODE(t, y, p), t, [cs1c0,cs2c0,cs3c0,cs4c0,cs5c0,cs6c0,ce1c0,ce2c0,ce3c0,ce4c0]);
Error in RunFile1 (line 76)
[cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4,pH] = funcKinetics1(t, cs1c0,cs2c0,cs3c0,cs4c0,cs5c0,cs6c0,ce1c0,ce2c0,ce3c0,ce4c0, p);
>>
I am not sure which two elements are reffered to in this error message.
my code is:
function [cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4,pH] = funcKinetics1(t, cs1c0,cs2c0,cs3c0,cs4c0,cs5c0,cs6c0,ce1c0,ce2c0,ce3c0,ce4c0, p)
[~, y] = ode23s(@(t,y) funcODE(t, y, p), t, [cs1c0,cs2c0,cs3c0,cs4c0,cs5c0,cs6c0,ce1c0,ce2c0,ce3c0,ce4c0]);
cs1 = y(:,1); % High Reactive Lignin
cs2 = y(:,2); % Low Reactive Lignin
cs3 = y(:,3); % Cellulose
cs4 = y(:,4); % Xylan
cs5 = y(:,5); % Glucomannan
cs6 = y(:,6); % Acetyls
ce1 = y(:,7); % Total dissolved lignin
ce2 = y(:,8); % Total dissolved carbohydrates
ce3 = y(:,9); % Total acetic acid
ce4 = y(:,10); % Total sulphite
pH = zeros(length(t), 1);
for i = 1 : length(t)
[~, pH(i)] = funcOtherVariables(cs1(i),cs2(i),cs3(i),cs4(i),cs5(i),cs6(i),ce1(i),ce2(i),ce3(i),ce4(i), p);
end
end
function dy_dt = funcODE(t, y, p)
cs1 = y(1); % High Reactive Lignin
cs2 = y(2); % Low Reactive Lignin
cs3 = y(3); % Cellulose
cs4 = y(4); % Xylan
cs5 = y(5); % Glucomannan
cs6 = y(6); % Acetyls
ce1 = y(7); % Total dissolved lignin
ce2 = y(8); % Total dissolved carbohydrates
ce3 = y(9); % Total acetic acid
ce4 = y(10); % Total sulphite
% Calculate the reaction rate.
r = funcOtherVariables(cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4, p);
% Mass balance (derivative terms)
dcs1_dt = r(1,1);
dcs2_dt = r(2,1);
dcs3_dt = r(3,1);
dcs4_dt = r(4,1);
dcs5_dt = r(5,1);
dcs6_dt = r(6,1);
dce1_dt = - p.ve1R1*r(1,1) - p.ve1R2*r(1,1) - p.ve1R3*r(2,1) - p.ve1R4*r(2,1);
dce2_dt = - p.ve2R5*r(3,1) - p.ve2R6*r(4,1) - p.ve2R7*r(5,1);
dce3_dt = - p.ve3R8*r(6,1);
dce4_dt = p.ve4R1*r(1,1) + p.ve4R2*r(1,1) + p.ve4R3*r(2,1) + p.ve4R4*r(2,1)...
+ p.ve4R5*r(3,1) + p.ve4R6*r(4,1) + p.ve4R7*r(5,1);
% Converting to vector form:
dy_dt = [dcs1_dt;dcs2_dt;dcs3_dt;dcs4_dt;dcs5_dt;dcs6_dt;dce1_dt;dce2_dt;dce3_dt;dce4_dt];
end
function [r, pH] = funcOtherVariables(cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4, p)
pH = fzero(@(pH) funcCharge(pH, cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4, p), 7);
cH = 10^-pH;
% Rate Equations
r1 = -p.k1*cs1^p.vs1R1 * (ce4*((ce4*p.KH2SO3*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R1 ...
-p.k2*cs1^p.vs1R2 * (ce4*((ce4*cH*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R2;
r2 = -p.k3*cs2^p.vs2R3 * (ce4*((ce4*p.KH2SO3*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R3 ...
-p.k4*cs2^p.vs2R4 * (ce4*((ce4*cH*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R4;
r3 = -p.k5*cs3^p.vs3R5 * (ce4*((ce4*p.KH2SO3*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R5;
r4 = -p.k6*cs4^p.vs4R6 * (ce4*((ce4*cH*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R6;
r5 = -p.k7*cs5^p.vs5R7 * (ce4*((ce4*p.KH2SO3*p.KHSO3)/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3)))^p.ve4R7;
r6 = -p.k8*cs6^p.vs6R8 * cH^p.vHR8 - p.k9*cs6^p.vs6R9 * (p.Kw/cH)^p.vOHR9;
r = [r1;r2;r3;r4;r5;r6];
end
function z = funcCharge(pH, p)
cH = 10^-pH;
% Left-hand side of the charge balance
LHS = p.ce6 + cH;
% Right-hand side of the charge balance
RHS = p.Kw/cH + 2*(s.ce4*p.KH2SO3*p.KHSO3/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*P.KHSO3)) +...
(s.ce4*cH*p.KH2SO3/(cH^2 + cH*p.KH2SO3 + p.KH2SO3*p.KHSO3))+ ...
2*(p.ce6*p.KH2CO3*p.KHCO3/(cH^2 + cH*p.KH2CO3 + p.KH2CO3*P.KHCO3)) +...
(s.ce4*cH*p.KH2CO3/(cH^2 + cH*p.KH2CO3 + p.KH2CO3*p.KHCO3))+ ...
(s.ce3*p.KAH)/(p.KAH + cH);
z = LHS - RHS; % Should be equal to zero
end
And I am running it using:
%% Parameters
% Rate constants (various units)
p.k1 = 2.95e-1;
p.k2 = 2.0e-2;
p.k3 = 1.0e-1;
p.k4 = 3.0e-1;
p.k5 = 3.0e-3;
p.k6 = 2.2;
p.k7 = 1.0e-2;
p.k8 = 5.0e-2;
p.k9 = 3.0e3;
% Equilibrium constants (various units)
p.KH2CO3 = 1e-3;
p.KHCO3 = 1e-7;
p.KH2SO3 = 1e-3;
p.KHSO3 = 1e-7;
p.KCH2OOH = 1e-4;
p.KCH3COOH = 1e-3;
p.Kw = 1e-8;
p.KAH = 1e-3;
% Constant concentration species (mol/L)
p.ce5 = 1; % Total carbonates
p.ce6 = 1; % Total sodium ions
% Stoichiometric coefficients
p.vs1R1 = 1;
p.vs1R2 = 1;
p.vs2R3 = 1;
p.vs2R4 = 1;
p.vs3R5 = 1;
p.vs4R6 = 1;
p.vs5R7 = 1;
p.vs6R8 = 1;
p.vs6R9 = 1;
p.ve1R1 = 1;
p.ve1R2 = 1;
p.ve1R3 = 1;
p.ve1R4 = 1;
p.ve2R5 = 1;
p.ve2R6 = 1;
p.ve2R7 = 1;
p.ve3R8 = 1;
p.ve4R1 = 1;
p.ve4R2 = 1;
p.ve4R3 = 1;
p.ve4R4 = 1;
p.ve4R5 = 1;
p.ve4R6 = 1;
p.ve4R7 = 1;
p.vHR8 = 1;
p.vOHR9 = 1;
%--------------------------------------------------------------------------
%% Initial conditions
cs1c0 = 1; % High Reactive Lignin
cs2c0 = 1; % Low Reactive Lignin
cs3c0 = 1; % Cellulose
cs4c0 = 1; % Xylan
cs5c0 = 1; % Glucomannan
cs6c0 = 1; % Acetyls
ce1c0 = 1; % Total dissolved lignin
ce2c0 = 1; % Total dissolved carbohydrates
ce3c0 = 1; % Total acetic acid
ce4c0 = 1; % Total sulphite
%------------------------------------------------------
%% Numerical integration
t = linspace(0, 100, 1);
[cs1,cs2,cs3,cs4,cs5,cs6,ce1,ce2,ce3,ce4,pH] = funcKinetics1(t, cs1c0,cs2c0,cs3c0,cs4c0,cs5c0,cs6c0,ce1c0,ce2c0,ce3c0,ce4c0, p);
I thank you.

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

采纳的回答

Bjorn Gustavsson
Bjorn Gustavsson 2021-3-26
编辑：Bjorn Gustavsson 2021-3-26
Your t variable will only have one component as you define it. My guess is that you want 100 elements between 0 and 1. If that's right you simply swapped the arguments to linspace:
t = linspace(0, 1, 100);
HTH

更多回答（0 个）

类别

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

标签

Community Treasure Hunt

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

Start Hunting!

Translated by