ode45 function matrix dimensions problem

2 views (last 30 days)
Carey n'eville
Carey n'eville on 8 Jan 2021
Edited: Star Strider on 8 Jan 2021
Hello friends I wrote a code. But this does not work. In this problem I have a parameter with range. I could not figure out how can I write. I tried something but it didn't work. I need to observe different Y_H values effects on system.
Parameter with range:
Y_H=[0.6 0.7]
CODE
clear all;
clc;
close all;
global mu_H Y_H Ks X_BH teta Ss_in
mu_Hmax=0.15; %days^-1
T=22.6; %centigrade
X_BH=2295; %mg/l
teta=1.25; %day
Y_H=[0.6 0.7];
Ks=250; %mg/l
%initial condition is Ss_in=1127 mg/l
Ss_in=1127; %mg/l
mu_H=mu_Hmax*1.072^(T-20); %day^-1
tspan=[30 32 34 36 50 52 57 59 61 63 65 72 75 76 79 83 85 87 89 96 99 106 109 116 118 121 124 125 126 127 131 133 134 135 136 137 140];
ode45(@substrate,tspan,Ss_in)
plot(t,Ss);
xlabel('Time (day)')
ylabel('Concentrations (mg/L)')
legend ('Substrate')
function Derivative = substrate (t,Ss)
global mu_H Y_H Ks X_BH teta Ss_in
Derivative = -((mu_H/Y_H)*(Ss/(Ks+Ss))*X_BH)+((1/teta)*(Ss_in-Ss)); %dSs/dt
Derivative=Derivative(:);
end
ERRORS I HAVE
Error using /
Matrix dimensions must agree.
Error in Solutionofmodeleqn31>substrate (line 26)
Derivative = -((mu_H/Y_H)*(Ss/(Ks+Ss))*X_BH)+((1/teta)*(Ss_in-Ss)); %dSs/dt
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in Solutionofmodeleqn31 (line 18)
ode45(@substrate,tspan,Ss_in)

Answers (2)

Alan Stevens
Alan Stevens on 8 Jan 2021
You need to loop through the Y_H values. See below
global mu_H Ks X_BH teta Ss_in
mu_Hmax=0.15; %days^-1
T=22.6; %centigrade
X_BH=2295; %mg/l
teta=1.25; %day
Y_H=[0.6 0.7];
Ks=250; %mg/l
%initial condition is Ss_in=1127 mg/l
Ss_in=1127; %mg/l
mu_H=mu_Hmax*1.072^(T-20); %day^-1
tspan=[30 32 34 36 50 52 57 59 61 63 65 72 75 76 79 83 85 87 89 96 99 ...
106 109 116 118 121 124 125 126 127 131 133 134 135 136 137 140];
Ss = zeros(numel(t),numel(Y_H));
% Loop through Y_H values
for i = 1:numel(Y_H)
[t,Ss(:,i)] =ode45(@(t,Ss) substrate(t,Ss,Y_H(i)),tspan,Ss_in);
end
plot(t,Ss);
xlabel('Time (day)')
ylabel('Concentrations (mg/L)')
legend ('Substrate')
function Derivative = substrate (~,Ss,Y_H)
global mu_H Ks X_BH teta Ss_in
Derivative = -((mu_H/Y_H)*(Ss/(Ks+Ss))*X_BH)+((1/teta)*(Ss_in-Ss)); %dSs/dt
Derivative=Derivative(:);
end

Star Strider
Star Strider on 8 Jan 2021
Edited: Star Strider on 8 Jan 2021
Do not use global variables!
Pass the additional parameters as arguments to the function, and use an appropriate anonymous function call in ode45.
Try this:
mu_Hmax=0.15; %days^-1
T=22.6; %centigrade
X_BH=2295; %mg/l
teta=1.25; %day
Y_H=[0.6 0.7];
Ks=250; %mg/l
%initial condition is Ss_in=1127 mg/l
Ss_in=1127; %mg/l
mu_H=mu_Hmax*1.072^(T-20); %day^-1
tspan=[30 32 34 36 50 52 57 59 61 63 65 72 75 76 79 83 85 87 89 96 99 106 109 116 118 121 124 125 126 127 131 133 134 135 136 137 140];
for k = 1:numel(Y_H)
[tv,Ssv(:,k)] = ode45(@(t,Ss)substrate(t, Ss, mu_H, Y_H(k), Ks, X_BH, teta, Ss_in),tspan,Ss_in);
end
figure
plot(tv,Ssv);
xlabel('Time (day)')
ylabel('Concentrations (mg/L)')
grid
legend (compose('Substrate: Y_H = %.1f',Y_H))
function Derivative = substrate (t,Ss,u_H, Y_H, Ks, X_BH, teta, Ss_in)
Derivative = -((mu_H/Y_H)*(Ss/(Ks+Ss))*X_BH)+((1/teta)*(Ss_in-Ss)); %dSs/dt
Derivative=Derivative(:);
end
EDIT — (8 Jan 2021 at 16:08)
Corrected typographical error.

Categories

Find more on Programming in Help Center and File Exchange

Products


Release

R2020b

Community Treasure Hunt

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

Start Hunting!