Problem in solving a set of 5 odes and getting Nan values
3 views (last 30 days)
Show older comments
Ghs Jahanmir
on 15 Sep 2021
Commented: Ghs Jahanmir
on 15 Sep 2021
Dear All,
I am trying to solve a set of 5 ode equaions using ode45 solver by using a separate independent function (named "Multi_func" contains my equations. However, they are solved with NaN values in final "eq" matrix.
any help is appreciated.
%%%%% Parameters
Q=[1100 350 1240 930 3620];
V_i=[3180 1696 1065 27000 2670];
V_t=[2290 1610 300 26649 0 ];
k_el=[0 0 90 0 0];
g=0;
tpan=[0 1];
Vm=[0 1.19e-4 2e-5 0 0];
K_M=[0 27 32 0 0];
W=6;
initial_conditions=zeros(1,W-1);
initial_conditions(5)=10; % ug/ml
[t,eqs]=ode45(@(t,eqs) Multi_func(t,eqs,Q,V_i,k_el,V_t,Vm,K_M,g,W),tpan,initial_conditions);
Time2=t./60; %convert to hr
cc=jet(W);
for i=1:W-1
plot(Time2,eqs(:,i),'color', cc(i,:));
hold on
end
my Multi_func function is:
function eq = Multi_func(t,eqs,Q,V_i,k_el,V_t,Vm,K_M,g,W)
eq=zeros(W-1,1);
eq(1)=((Q(1)*eqs(5)-Q(1)*eqs(1))/V_i(1))-(k_el(1)*eqs(5)/V_i(1))-V_t(1)*(((Vm(1)*eqs(1)))/(K_M(1)+eqs(1)))/V_i(1);
eq(2)=((Q(2)*eqs(5)-Q(2)*eqs(2))/V_i(2))+(Q(1)*eqs(1)/V_i(2))-(k_el(2)*eqs(5)/V_i(2))-V_t(2)*(((Vm(2)*eqs(2)))/(K_M(2)+eqs(2)))/V_i(2);
eq(3)=((Q(3)*eqs(5)-Q(3)*eqs(3))/V_i(3))-(k_el(3)*eqs(5)/V_i(3))-V_t(3)*(((Vm(3)*eqs(3)))/(K_M(3)+eqs(3)))/V_i(3);
eq(4)=((Q(4)*eqs(5)-Q(4)*eqs(4))/V_i(4))-(k_el(4)*eqs(5)/V_i(4))-V_t(4)*(((Vm(4)*eqs(4)))/(K_M(4)+eqs(4)))/V_i(4);
eq(5)=(Q(1)*eqs(1)+Q(2)*eqs(2)+Q(3)*eqs(3)+Q(4)*eqs(4)-Q(5)*eqs(5))/V_i(5)+V_t(5)*(((Vm(5)*eqs(5)))/(K_M(5)+eqs(5)))/V_i(5)+g/V_i(5);
end
0 Comments
Accepted Answer
Walter Roberson
on 15 Sep 2021
K_M=[0 27 32 0 0];
First K_M is 0.
W=6;
initial_conditions=zeros(1,W-1);
initial_conditions(5)=10; % ug/ml
First boundary condition is 0
eq(1)=((Q(1)*eqs(5)-Q(1)*eqs(1))/V_i(1))-(k_el(1)*eqs(5)/V_i(1))-V_t(1)*(((Vm(1)*eqs(1)))/(K_M(1)+eqs(1)))/V_i(1); ^^^^^^^^^^^^^^^^
Notice the /(K_M(1)+eqs(1)) . As pointed out above, K_M(1) is 0 an eqs(1) starts out 0, and sum of those two is 0, so you have a division by 0 there.
3 Comments
Walter Roberson
on 15 Sep 2021
Vm(1) is also zero. It shuld be able to cross off the (((Vm(1)*eqs(1)))/(K_M(1)+eqs(1)))/V_i(1).
Try it.
(5*0)/0
5*(0/0)
(0*5)/0
0*(5/0)
(0/0)*5
0/0*0
5*sin(0)/0
syms x real
expr = 5*(sin(x)/x)
limit(expr, x, 0, 'left')
limit(expr, x, 0, 'right')
5 * x / x
That is, you can talk about limits, but the symbolic engine does not automatically calculate in terms of limits.
More Answers (0)
See Also
Categories
Find more on Symbolic Math Toolbox in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!