# Solving coupled differential equations with DSolve/NDSolve

I am trying to solve system of coupled differential equations as follows, but somehow its not returning me any solution, with or without parametric values using DSolve/NDSolve. It doesn't return me anything.

The same model is working in OCTAVE solver.

Any help? Following is the code:

     kapc = 0.04;
drugdosage = 0;
katc = 2.079;
kpmc = 0.56545;
maxdrug = 250;
kiap = 5.4;

kiev = 1;
kt = 1;
kgir = 1;
katg = 1;

kdtc = 0.3466; ktapc = 0.2605; g1 = 2.7840; ktatc = 0.2605; g2 = 2.7840; kx = 0.1; tdtc = 500; Mt = 1; g3 = 3.9913; kink = 3.9913; knkc = 3.1623; HLA = 1; nk = 1; vac = 1; dose = 0; dd = 0.002; gamma = 0.2; a = 0.5; k = 0.005;

eq1 := tcell'[t] ==
kapc*(1 + (kgir*(drugdosage))) *(((katg*m[t] + (drugdosage))^
g1)/(ktapc^g1 + (katg* m[t] + (drugdosage)^g1))) + (katc*
tcell[t]*(((katg*m[t])^
g2)/(ktatc^g2 + (katg*(m[t]))^g2))) -  (kdtc * ((tcell[
t]^3) + kx* (tcell[t]^0.1) * ((tcell[t])^0.9)));
eq2 := m'[
t] == (kpmc *m[t] * (Mt - m[t])) - (
kiap* tcell[
t] * ( (katg*m[t]) / ((1 + (kiev) )*m[t] ) )) - (kink*
m[t]*((katg^(-g3))/((knkc^g3) + (katg^(-g3)))));
eq3 := d'[t] ==  -(Gamma*d[t]) + maxdrug;
sol = NDSolve[{eq1, eq2, eq3, tcell[0] == 0, m[0] == 0.00081,
d[0] == 0, d[1.3] == maxdrug, d[21.3] == maxdrug,
d[42.3] == maxdrug, d[63.3] == maxdrug, d[84.3] == maxdrug,
d[105.3] == maxdrug, d[126.3] == maxdrug, d[147.3] == maxdrug,
d[168.3] == maxdrug, d[189.3] == maxdrug,
d[210.3] == maxdrug}, {tcell[t], m[t], d[t]}, {t, 0, 250}]
{sol1, sol2} = sol[[1, All, 2]];
Plot[{sol1, sol2}, {t, 0, 250}, PlotRange -> All]

• Analytical solution not possible with DSolve. Try NDSolve.
– zhk
Mar 17, 2017 at 4:06
• What is m[2]?
– zhk
Mar 17, 2017 at 4:09
• its actually m[t], typo while editing code here. NDSolve gives me following result does it makes sense? {tcell, m} /. DSolve[{Derivative[1][tcell][t] == (0.04 m[t]^2.784)/( 0.0236379 + m[t]) + (2.079 m[t]^2.784 tcell[t])/( 0.0236379 + m[t]^2.784) - 0.3466 (0.1 tcell[t]^1. + tcell[t]^3), Derivative[1][m][t] == -0.0399105 m[t] + 0.56545 (1 - m[t]) m[t] - (5.4 m[t] tcell[t])/(1 + m[t]), True, True}, {tcell, m}, t] Mar 17, 2017 at 4:18

After correcting the typo m[2], DSolve is still unable to find an analytical solution. So, the obvious choice is to find numerical one using NDSolve.

eq1 = tcell'[t] == kapc*(1 + (kgir*(drugdosage)))*(((katg*m[t] + (drugdosage))^
g1)/(ktapc^g1 + (katg*m[t] + (drugdosage)^g1))) + (katc*
tcell[t]*(((katg*m[t])^g2)/(ktatc^g2 + (katg*m[t])^g2))) - (kdtc*((tcell[t]^3) +
kx*(tcell[t]^0.1)*((tcell[t])^0.9)));

eq2 = m'[t] == (kpmc*m[t]*(Mt - m[t])) - (kiap*tcell[t]*((katg*m[t])/(1 + (kiev*m[t]))))
- (kink*m[t]*((katg^(-g3))/((knkc^g3) + (katg^(-g3)))));

eq3 = d'[t] == -(g1*d[t]) + maxdrug;

sol = NDSolve[{eq1, eq2, eq3, tcell[0] == 0, m[0] == 0.00081}, {tcell[t], m[t], d[t]},
{t, 0, 250}];

{tcsol, msol,dsol} = sol[[1, All, 2]];

Plot[{tcsol, msol, dsol}, {t, 0, 250}, PlotRange -> All]


• Thanks, but it does not work when I update the model and incorporate third equation.I am interested in plotting the values of tcell[t] and m[t] over the scale of time from t : 0 - 250 but it somehow doesn' t work.I have updated the model in the main question ... .. Mar 17, 2017 at 5:41
• @Osman Don't use := and Gamma. You just need one initial condition on d.
– zhk
Mar 17, 2017 at 6:57
• @ Maple SE - Area 51 Proposal The problem is that I want to introduce the values of d at certain intervals which will decay with time. Is there any other way I can do that? or using multiple conditions is wrong? Mar 18, 2017 at 17:18