1
$\begingroup$

I am trying to investigate some properties of DDEs with time-varying delays. As an approximation, I am trying to evaluate a basic DDE from time -10 to 10 with a specified constant delay, and then using the solution of that DDE from time 0 to 10 as the initial condition for a second DDE with a new constant delay over time 0 to 20, with the goal of repeating the process many times.

However, I cannot get NDSolve to accept the solution from the previous DDE as the initial condition of the new DDE. Here is what the code looks like.

Delay[x_] := -Exp[-x] + Pi/2; DDEstep[z_, y_] := NDSolve[{x'[t] == - x[t - Delay[10 z - 10]], x[t /; t <= 10 z - 10] == y}, x, {t, 10 z - 20, 10 z}]; step1 = DDEstep[1, 1]; step2 = DDEstep[2, x[t] /. step1[[1]]]

When I try to run it, I get the following error: "NDSolve: Initial history needs to be specified for all variables for delay-differential equations." NDSolve seems to find "x[t] /. step1[[1]]" to be an unacceptable form for the initial condition. For an ODE, it would be as easy as using Evaluate on the specific time of the initial condition (as has been described in other questions on the site) but since this is a DDE, I need my initial condition over a time interval. I know that user defined functions of t work when plugged in as initial conditions, so I would like someone to explain why "x[t] /. step1[[1]]" is not the same as a function of t.

Improve this question
$\endgroup$
4
  • 1
    $\begingroup$ If I change the initial time of integration to 10 z - 10 to match the time in the initial history, answers are computed. Do the solutions see right? $\endgroup$ Commented Jun 17, 2024 at 20:13
  • $\begingroup$ @MichaelE2 This has solved my problem. Thank you! One thing I am still confused about is that in the documentation for DDEs all of the examples have an initial time of integration that is less than or equal to the longest delay present in the equation. However, it would seem that is not necessary, and even caused errors for me. $\endgroup$ Commented Jun 17, 2024 at 22:28
  • $\begingroup$ You're welcome. I thought I had an answer, but after thinking about it, I'm not sure I understand either. Usually the domain $[a,b]$ as specified by {t, a, b} in NDSolve[] is independent of the DE. It does not even have to include initial or boundary conditions. I think there must be some internal test for the initial history that makes a mistake when it is given by an InterpolatingFunction. See the second method in my answer for a way to trick NDSolve into accepting it. If the solution seems wrong, let me know. $\endgroup$ Commented Jun 18, 2024 at 15:14
  • $\begingroup$ Thanks for the accept. I thought such a well-written question exposing a problem would have attracted more interest. In any case, I'm going to mark it a bug. If the community disagrees, we can open a meta Q&A. You can report the bug to WRI from within Mathematica, menu Help > Give Feedback.... $\endgroup$ Commented Jun 19, 2024 at 15:59

2 Answers 2

Reset to default
0
$\begingroup$

Here are two different ways of fixing it. The second way suggests that InterpolatingFunction is not handled or perhaps not allowed as a way to specify initial history in all cases. You might ask Wolfram about it, as it may be a bug. I'm not deeply familiar with how DDEs work in NDSolve[], so there may be some (un)documented restriction I don't know about.

I tend to follow the advice here, What are the most common pitfalls awaiting new users? (link to answer; see point #3), and avoid starting my variables and functions with capitals.

Method 1: If I change the initial time of integration to 10 z - 10 to match the time in the initial history, answers are computed. I don't think specifying the domain this way should be necessary.

delay[x_] := -Exp[-x] + Pi/2; ddeStep[z_, y_] := NDSolve[{x'[t] == -x[t - delay[10 z - 10]], x[t /; t <= 10 z - 10] == y}, x, {t, 10 z - 10, (* fix is here *) 10 z}]; step1 = ddeStep[1, 1] f[t_?NumericQ] = x[t] /. step1[[1]]; step2 = ddeStep[2, x[t] /. step1[[1]]]

Method 2: Hide InterpolatingFunction inside a _?NumericQ protected function; more advice here: What are the most common pitfalls awaiting new users? (link to another answer).

delay[x_] := -Exp[-x] + Pi/2; ddeStep[z_, y_] := NDSolve[{x'[t] == -x[t - delay[10 z - 10]], x[t /; t <= 10 z - 10] == y}, x, {t, 10 z - 20, 10 z}, "ExtrapolationHandler" -> (* for convenient plotting *) {Indeterminate &, "WarningMessage" -> False}]; step1 = ddeStep[1, 1] f[t_?NumericQ] = x[t] /. step1[[1]]; (* fix here and next line *) step2 = ddeStep[2, f[t]] Plot[x[t] /. Join[step1, step2] // Evaluate, {t, 0, 20}, PlotStyle -> {AbsoluteThickness[4], AbsoluteThickness[2]}, PlotRange -> 0.1, PlotLegends -> {HoldForm[step1], HoldForm[step2]}]

enter image description here

Apparently you can postprocess step2 to store the solution in step1 instead of the symbol f:

step2IF = step2 /. f -> (x /. step1[[1]])

Compare

step2[[1, 1, 2, 4, 1]] step2IF[[1, 1, 2, 4, 1]]

Caveat: I've played with InterpolatingFunction a lot, and it does not always respond so well to such open-heart surgery. However, I tested the above plot with step2IF and after Clear[f], and there were no problems. Clearing f did mean the original step2 no longer worked. (I did save all my work beforehand: By "not so well," I mean it has sometimes crashed both the kernel and the front end.)

Improve this answer
$\endgroup$
0
$\begingroup$

Here is how you can implement the approach:

ClearAll["Global`*"] (* Define the delay function *) Delay[x_] := -Exp[-x] + Pi/2; (* Function to solve a DDE from t0 to t1 and return the solution and initial condition *) DDEstep[z_, y_, t0_, t1_] := Module[{sol, init}, sol = First@NDSolve[{x'[t] == -x[t - Delay[10 z - 10]], x[t0] == y}, x, {t, t0, t1}]; init = x[t0] /. sol; {sol, init} ]; (* Initialize parameters *) zmax = 5; (* Maximum number of steps *) tmax = 20; (* Maximum final time for the entire process *) (* List to store solutions *) solutions = {}; (* Solve the sequence of DDEs *) Do[ If[i == 1, {sol, init} = DDEstep[i, 1, 10 i - 20, 10 i], {sol, init} = DDEstep[i, init, 10 i - 20, 10 i] ]; AppendTo[solutions, sol], {i, zmax} ]; (* Plot the solution of the last DDE *) Plot[Evaluate[x[t] /. Last[solutions]], {t, 0, tmax}, PlotRange -> All, PlotLabel -> "Solution of the last DDE with successive initial conditions"]

enter image description here

Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.