Problem in solving a set of 5 odes and getting Nan values

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

Answer 1

[Walter Roberson]

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.

[Ghs Jahanmir]

Thanks @Walter Roberson

But this is the case for eq(4).

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).

the same set of equations with same initial condition can be solved using symbolic version and converted to matlabfunction with same ode45.

I just wanted to make sure about my solution so I did the older approach of defining a separate funtion and call it.

C_vec is vector of 5 symolic y variables

eq is symbolic set of 5 equations

eqns = diff(C_vec,t)==eq;
[VF1,Sbs] = odeToVectorField(eqns);
VF=vpa(VF1);
odseSetfunctions = matlabFunction(VF,'Vars',{t,Y});

tpan=[0 10];
initial_conditions=zeros(W-1,1);
initial_conditions(5)=1000;  % ug/ml
[time,Ymatrix]=ode15s(@(T,Y)odseSetfunctions(T,Y),tpan,initial_conditions);

[Walter Roberson]

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

ans = NaN

5*(0/0)

ans = NaN

(0*5)/0

ans = NaN

0*(5/0)

ans = NaN

(0/0)*5

ans = NaN

0/0*0

ans = NaN

5*sin(0)/0

ans = NaN

syms x real
expr = 5*(sin(x)/x)

expr = 

limit(expr, x, 0, 'left')

ans = 

5

limit(expr, x, 0, 'right')

ans = 

5

5 * x / x

ans = 

5

That is, you can talk about limits, but the symbolic engine does not automatically calculate in terms of limits.

[Ghs Jahanmir]

@Walter Roberson

Yes, I got the point.

Thanks a lot for your time.