Differential equation in matlab

Question

A pathogen population ($P(t)$) under the effect of a drug($A(t)$) is modelled as

${dP\over dt}=rP-{E_{max}A^n\over C_{50}^n+A^n}P$ and $A(t)$ is modelled as

$A(t)={D\over V}(e^{-k(t-t_0)}-e^{-k_a(t-t_0)})$.

From literature I know the values of $D,V,t_0$(lag time), $k,k_a$. So, $A(t)$ can be obtained at the respective time points the ODE is run.
$r,E_{max},n,C_{50}$ are also known.

I want to solve the ODE ${dP\over dt}$ using ode45 of Matlab. (Without obtaining the analytical solution of $P(t)$ as this is only a simplified system and the original system cannot be solved analytically).

Can someone please let me know if the below code is right in solving $P(t)$. I am not sure about the way $A(t)$ is included into the ODE equation.

y0=10;
[t,y]=ode45(@Model,1:200,y0);

function s= Model(t,y)
r=0.4;
Emax=50;
n=1;
c_50=0.1;
k=0.25;
D=40;
V=6.6;
t0=1;
k_a=0.5;
s=zeros(1,1);
A=(D/V)*(exp(-k*(t-t0))-exp(-k_a*(t-t0)));
s(1)=r*y(1)-((Emax*(A^n))/((c_50^n)+(A^n)))*y(1);
end

Answer 1

Your time vector likely has too large steps for solver. Try replace

1:200

with

linspace(0,0.125,256)

Answer 2

I copied your code to Octave (it runs most matlab code) and the results seem wrong, as the solution "explodes". This could be a precision in ode45. I attach the Wolfram Mathematica solution to the same problem. The solution decays to zero almost immediately.