Variable inside an ode

15 次查看（过去 30 天）
Nora Rafael
Nora Rafael 2019-11-13
评论： Star Strider 2019-11-13
I am trying to solve this differential equation: dCb/dt = ((A*D)/(DV*h))*(Cs-Cb). I specify in Matlab the function:
f = @(t,Cb) ((A*D)/(h*DV))*(Cs-Cb)
I use ode45 to solve it:
[t,Cb] = ode45(f, tspan, Cb0);
Full code is below.
However, A is a variable and varies with time itself. It is actually A(t) =0.00155+(0.06043-0.0015)/(1+(t/4481.23119).^0.85949)
How do I add this A function to the below code?
Thanks
tspan = [0 9e4];
dv50=189/1000000000; %m,
dv90=795/1000000000;
r0=dv90/2; %m dv90
rho=1370; %kg/m3
W=18/1000000 %kg
N= W/((4/3)*pi*r0^3*rho);
A=(4*pi*r0^2*N); %m2
MW=402.88
D= 9.9e-9*MW^-0.453 %m2/s % from D(m2/s)=9.9e-9*MW^-0.453
h= r0 %m -
Cs= 0.032 %kg/m3
DV=0.0009 % m3, volume
Cb0=0;
f = @(t,Cb) A*D/(h*DV)*(Cs-Cb)
[t,Cb] = ode45(f, tspan, Cb0);
plot(t,Cb);
title('Bulk Concentration vs Time');
xlabel('Time (s)');
ylabel('Cb (kg/m3)');
legend('Nernst Brunner','Cb_exp');
  2 个评论
Adam
Adam 2019-11-13
Can you not just replace A in the definition of f, with its above representation in terms of t which is already an input to f?
Nora Rafael
Nora Rafael 2019-11-13
This seemed to work with no errors. But I am concerned as when I plot just this A function alone I get a strange curve.
This is the equation:
equation.JPG
and I implement it this way
SA=0.00155+(0.06043-0.0015)/(1+(t/4481.23119).^0.85949)
This is the curve that I get
Capture.JPG
This is a strange curve as this function should output an exponential decline.
I tested this function separately
function SA_=SurfaceArea(t)
SA_=0.00155+(0.06043-0.0015)/(1+(t/4481.23119).^0.85949);
end
>> t=1:9e4;
>> SA_=SurfaceArea(t)
and I get this error message:
Error using /
Matrix dimensions must agree.
Error in SurfaceArea (line 2)
SA_=0.00155+(0.06043-0.0015)/(1+(t/4481.23119).^0.85949);

请先登录，再进行评论。

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

采纳的回答

Star Strider
Star Strider 2019-11-13
Try this:
SA = @(t) 0.00155+(0.06043-0.0015)./(1+(t/4481.23119).^0.85949);
f = @(t,Cb) SA(t)*D/(h*DV)*(Cs-Cb)
The rest of your code is unchanged.
  4 个评论
显示 2更早的评论

请先登录，再进行评论。

更多回答（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