# How to Solve Functional Differential Equation

I tried to solve this equation numerically.

eqnEx = x''[t] + x[2 t] == 0;
NDSolve[{eqnEx, x[1] == 10, x'[1] == 0}, x[t], {t, 1, 10}]


Two important character of the equation are that function variable is 2t and initial condition is given at t=1 .

After Calculation, Mathematica gives this error message

NDSolve : The method currently implemented for delay differential equations does not support delays that depend directly on the time variable or dependent variables


After google, I find these kind of equation are called Functional Differential Equation and I could not obtain how to solve it. Is there a way to solve these kind of equation numerically?

• If you set $t=1$ the system remains ill defined since x[2] is not known? – chris May 2 at 6:37
• Maple 2021 Can solve this with numerics. – Mariusz Iwaniuk May 2 at 11:35
• @MariuszIwaniuk Please provide the answer from Maple? – bbgodfrey May 2 at 13:32
• @chris yes x[2] is not known. Edit : Is the problem solvable if x[2] provided? – OneMoreLovely May 3 at 8:53
• Welcome to Mathematica.SE! I hope you will become a regular contributor. To get started, 1) take the introductory tour now, 2) when you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge, 3) remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign, and 4) give help too, by answering questions in your areas of expertise. – bbgodfrey May 5 at 4:25

Based on the suggestion by Daniel Lichtblau, the functional differential equation in the question can be converted to a delay differential equation as follows. With

rult = t -> 2^-s;


2 t becomes 2^-(s-1), as desired. To proceed,

ruls = Simplify[Solve[Equal @@ rult, s] /. C[1] -> 0, t > 0][[1, 1]]
(* s -> -Log[t]/Log[2] *)


and the range of integration {t, 1, 10} becomes {s, smin, 0} with

smin = s /. ruls /. t -> 10
(* -Log[10]/Log[2] *)


The differential equation and boundary conditions themselves are transformed by

D[x[s /. ruls], t] /. Reverse[ruls] /. rult
(* -2^s x'[s]/Log[2] *)
Simplify[D[x[s /. ruls], t, t] /. Reverse[ruls] /. rult]
(* 4^s (Log[2] x'[s] + x''[s])/Log[2]^2 *)
eq = Simplify[(2^-s Log[2])^2 (%78 + x[s - 1])]
(* -4^-s Log[2]^2 x[-1 + s] + Log[2] x'[s] + x''[s] *)


and

x[0] == 10;
x'[0] == 0;


Unfortunately, the integration must begin at smin, as required by DSolve, and the boundary conditions are at 0. The problem can be solved, nonetheless.

sol = ParametricNDSolveValue[{eq == 0, x[s /; s < smin] == xp + (s - smin) xpp,
x'[s /; s < smin] == xpp}, {x[0], x'[0]}, {s, smin, 0}, {xp, xpp}];
Values[FindRoot[sol[xp0, xpp0] - {10, 0}, {xp0, 1}, {xpp0, 1}, Evaluated -> False]]
(* {-0.788537, -0.990015} *)


effectively transferring the boundary conditions to smin. Finally,

NDSolveValue[{eq == 0, x[s /; s < smin] == #1 + (s - smin) #2,
x'[s /; s < smin] == #2}, x[s], {s, smin, 0}] & @@ %;
ParametricPlot[{t /. rult, %}, {s, smin, 0}, AxesLabel -> {t, x}, PlotRange -> All,
LabelStyle -> {15, Bold, Black}, ImageSize -> Large, AspectRatio -> 1/GoldenRatio]


I close with the important caveat that the solution depends, as it must, on the assumed form of x[s] in the range {s, smin - 1, smin}, which is not specified in the question. Consequently, I made the simplest reasonable assumption, that x[s] varies linearly in that range.

A small amount of experimentation reveals that decreasing smin to -6 yields a solution insensitive to details of x[s] for s less than smin (i.e., for large t). Superimposing the two curves above shows very good, although not perfect, agreement over {t, 1, 10}.

• After considerate review of your solution, I think Equaiton (* -4^-s Log[2]^2 x[-1 + s] + Log[2] x'[s] + x''[s] *) should be (* 4^-s Log[2]^2 x[-1 + s] + Log[2] x'[s] + x''[s] *) not a minus sign. Then the solution will be same as mine. Still, your solution is so helpful and I appreciate it. – OneMoreLovely May 5 at 11:02
• @OneMoreLovely You are correct. I miscopied the original equation. I shall correct the answer now. Thanks. – bbgodfrey May 5 at 11:38

@bbgodfrey I tried to follow your method to solve the equation. The difference is that we choose different transformation of variables.

Step 1) make Delay Diff Eq using substitution by t = Exp[-s]

testfun = x''[t] + x[2 t]
testfun /. {x -> (x[-Log[#]] &)} /. t -> Exp[-s] // PowerExpand


Step 2) The substitution makes x[s=0]=10 and x'[s=0]=0are not the initial condition as you pointed out. The actual initial condition should be given at s=-Log[10]. I followed same logic as you. only difference is the result {-0.788537, -1.42829} not {28.1156, 49.9555} because of different parametrization.

Step 3) solve Delay Diff Eq with linear-assumed boundary condition.

sol2 = NDSolveValue[{eqn == 0, x[s /; s < smin] == -0.7885374724851758 + (s -
smin) (-1.428289649546688),
x'[s /; s < smin] == -1.428289649546688}, x, {s, smin, 0}]


Step 4) change variable s to t and plot. it looks like something different comparing your answer.

ParametricPlot[{t /. t -> Exp[-s], sol2[s]}, {s, smin, 0},
AxesLabel -> {t, x}, PlotRange -> All, LabelStyle -> {15, Bold, Black},
ImageSize -> Large, AspectRatio -> 1/GoldenRatio, Frame -> True]


Step 5) Check whether or not the obtained solution follows Diff equation in my question. We have checked with x[s] and Equation modified in Step 1. x[s - Log[2]] + E^(2 s) Derivative[1][x][s] + E^(2 s) (x^\[Prime]\[Prime])[s]==0

Derivative[2][sol2][-1] + Derivative[1][sol2][-1] + Exp[-2 (-1)] sol2[-1 - Log[2]] // N (* 4.48826*10^-8 *)
Derivative[2][sol2][-1.2] + Derivative[1][sol2][-1.2] + Exp[-2 (-1.2)] sol2[-1.2 - Log[2]] // N (* -2.452*10^-6 *)
Derivative[2][sol2][-1.6] + Derivative[1][sol2][-1.6] + Exp[-2 (-1.6)] sol2[-1.6 - Log[2]] // N (* -0.0000118291 *)


It seems our solution is consistent with the equation. I really concern why there is difference between our two solutions. Thanks a lot @bbgodfrey.

I have only numeric solution by Maple 2021.

Maple code:

sol := dsolve([diff(x(t), t \$ 2) + x(2*t) = 0, x(1) = 10, D(x)(1) = 0], numeric, delaymax = 1):
plots[odeplot](sol, [[t, x(t)], [t, D(x)(t)]], t = -1 .. 10, legend = [typeset(x(t)), typeset(diff(x(t), t))])


• Thank you for providing this solution from Maple. I had expected that values of x for t outside {1, 10} would be required as input, as they are for delay differential equations. – bbgodfrey May 2 at 14:13
• @bbgodfrey. I check this solution by Maple and is very wrong.Probably bug in dsolve. – Mariusz Iwaniuk May 2 at 17:15
• @Mariusz Iwaniuk. I checked your answer and comment carefully. Your comment is totally right. Since the second derivative in [7,10] is zero which can not be the right solution of this equation. Thank you again for giving me a clue to use Maple solving delay differential equation. Is there any method to solve the question with Mathematica or debug with Maple code? I guess the issue might come from the definition of delaymax. – OneMoreLovely May 3 at 10:43
• @OneMoreLovely Yours is not a delay differential equation, and NDSolve has no options for solving it. I do have some ideas but will not have time to pursue them until the end of the week. Am I correct in assuming that x[t] must be infinitely differentiable, for Infinity > t > 1? – bbgodfrey May 3 at 13:00
• @bgodfrey It can be transformed to a DDE via t=log(s). Problem is, one can only go downward in time. Something like this, if I did the algebra correctly (unlikely): NDSolveValue[{x[s + 1] == Log[2]^2*(x'[s]*s + x''[s]*s^2), x[2] == 10, x'[2] == 1/10}, {x[s]}, {s, 2, .1}]. – Daniel Lichtblau May 3 at 13:25