# System of differential equation with delay

I try to solve this differential system with delay:

DSolve[x'[t] == -0.12  x[t]  -
0.8 *  1.77   x[t] /(1 + x[t]^4)
E^(-0.1  *
2.8    ) (1.6 * 1.77  1/(1 + x[t - 2.8]^4) +  0.01) x[
t - 2.8] - (0.51/61)^-0.001  p[t]/(
p[t] + 1) ((2.5*10^4 - 2.5*10^4 * 0.4 * 0.01)/(2.5*10^4)^0.99) x[
t]^0.99, d'[t] == - 0.61  d[t] + 200 ,
p'[t] == -0.0412 p[t] + 0.61 d[t], {x[t], d[t], p[t]}, t]


In the render I wasn't able to understand the error message.

Is there any charitable sool to help me out.

• You need your three equations inside { and }. You have x'[t], y'[t] and p'[t] but are solving for x[t], d[t] and p[t]. There is still at least one more error I can't see at the moment. – Bill Feb 17 '15 at 7:22
• DSolve might not recognize that x[t - 2.8] is just the time delay of x[t]. By the way, can DSolve solve time delay differential equation? – Harry Feb 17 '15 at 7:32
• @Bill Excuse me. p[t] instead y[t], I give the correction. – Zbigniew Feb 17 '15 at 7:38
• @Harry In mathematica 10 it is possible to sove an EDO with delay. – Zbigniew Feb 17 '15 at 7:39
• @Zbigniew: First, you will have to enclose the equations in curly brackets, then you might need to specify the historic behavior, e.g. like x[t/;t<=0]==(*whatever*), and finally, you need to specify the range of t, e.g. {t,0,2}. – Jinxed Feb 17 '15 at 8:54

With the corrections suggested by Jinxed (with whatever set to 1) plus the additional correction of specifying the range of t (required for a delay differential equation), the code becomes

DSolve[{x'[t] == -0.12 x[t] - 0.8*1.77 x[t]/(1 + x[t]^4) E^(-0.1*2.8) (1.6*1.77 1/(1 +
x[t - 2.8]^4) + 0.01) x[t - 2.8] - (0.51/61)^-0.001 p[t]/(p[t] + 1)
((2.5*10^4 - 2.5*10^4*0.4*0.01)/(2.5*10^4)^0.99) x[t]^0.99,
d'[t] == -0.61 d[t] + 200,
p'[t] == -0.0412 p[t] + 0.61 d[t],
x[t /; t <= 0] == 1}, {x[t], d[t], p[t]}, ,{t, 0, 10}]


Since DSolve runs "forever", trying to solve this problem, I simplified the code by solving for d and p, which are independent of x.

ans1 = First@DSolve[{d'[t] == -0.61 d[t] + 200, p'[t] == -0.0412 p[t] + 0.61 d[t]}, {d[t],p[t]}, t]
(* {d[t] -> 327.869 + 1. E^(-0.61 t) C[1],
p[t] -> 1. E^(-0.0412 t) (5205.99 E^(0.0412 t) - 351.617 E^(0.61 t)) +
351.617 E^(-0.0412 t) (-1. E^(0.0412 t) + 1. E^(0.61 t)) +
1.07243 E^(-0.6512 t) (-1. E^(0.0412 t) + 1. E^(0.61 t)) C[1] + 1. E^(-0.0412 t) C[2]} *)


These then are substituted into the remaining differential equation

eq = Simplify[x'[t] == -0.12 x[t] - 0.8*1.77 x[t]/(1 + x[t]^4) E^(-0.1*2.8)
(1.6*1.77 1/(1 + x[t - 2.8]^4) + 0.01) x[t - 2.8] - (0.51/61)^-0.001 p[t]/(p[t] + 1)
((2.5*10^4 - 2.5*10^4*0.4*0.01)/(2.5*10^4)^0.99) x[t]^0.99 /. ans1]


The resulting expression, which is rather long to be displayed here, when inserted into DSolve

DSolve[{%, x[t /; t <= 0] == 1}, x[t], t]


still runs "forever", perhaps due to the x[t]^0.99 term.

This being the case, I used NDSolveValue. Doing so requires that the constants C[1] and C[2] be specified, and I chose 0. (Alternatively, one could use ParametricNDSolveValue with the two constants as parameters.)

Simplify[eq/.{C[1] -> 0, C[2] -> 0}]
sol = NDSolveValue[{%, x[t /; t <= 0] == 1}, x, {t, 0, 10}]
LogPlot[sol[t], {t, 0, 10}, AxesLabel -> {t, x}, PlotRange -> {10^-8, 1}]