How to evaluate the "Domain" of a Parametric Function?

2

0

In a simple example I try to calculate the time the solution x[t] passes a second time x[t]==1.

X = ParametricNDSolveValue[{x''[t] + 2 .3 x'[t] + x[t] == 1, 
x[0] == x0, x'[0] == v0, flag[0] == 0
, WhenEvent[x[t] > 1 , {flag[t] -> 1}]
, WhenEvent [(x[t] < 1 ) && (flag[t] == 1) ,"StopIntegration"  ]}, 
x, {t, 0, 10}, {x0, v0}  ,DiscreteVariables -> {flag} ];
Plot[X[0, 1][t], {t, 0, X[0, 1]["Domain"][[1]][[-1]]},GridLines -> {{X[0, 1]["Domain"][[1]][[-1]]}, {1}}]

enter image description here

Depending on the initial conditions x0,v0 I now force the endtime to equal 4 which can be visualized in

ContourPlot[4 == X[x0, v0]["Domain"][[1]][[-1]], {x0, -.7, .7}, {v0, 0, 2},FrameLabel -> {x0, v0}]

enter image description here

Works! But if I try to calculate the contour for a special point

FindRoot[4 == X[0, v0]["Domain"][[1]][[-1]], {v0, 1}]
(*FindRoot[4 == X[0, v0]["Domain"][[1]][[-1]], {v0, 1}]*)

Mathematica can't evaluate!!!

What is wrong here? Thanks!

Ulrich Neumann

Posted 2019-02-08T13:48:50.850

Reputation: 25 490

I am a bit confused about ["Domain"] - it is given as a property for an InterpolatingFunction when I use PropertyList[ X[0,1]] for example. But Mathematica crashes (!) doing PropertyList[ {X[0,1], All } ] which seems to be how PropertyList should be used? – gwr – 2019-02-08T14:32:04.183

Instead of ["Domain"] you could also use [[1]]: X[0, 1][[1]] (*{{0., 4.27631}}*) – Ulrich Neumann – 2019-02-08T14:34:51.593

1I am getting Part::partd errors for the ContourPlot (which nontheless produces a plot) and also for the FindRoot which also gives FindRoot::nlnum error. – gwr – 2019-02-08T14:42:34.833

With Quiet@FindRoot[ With[{v = v0}, 4 == X[0, v]["Domain"][[1, -1]]], {v0, 1}] I am getting {v0 -> 1.49383}. – gwr – 2019-02-08T14:49:10.500

@gwr: Thanks, after testing your code FindRoot[4 == X[0, v0]["Domain"][[1, -1]], {v0, 1}] works, mysterious! – Ulrich Neumann – 2019-02-08T14:53:24.060

Answers

2

You need Evaluated -> False:

FindRoot[4 == X[0, v0]["Domain"][[1]][[-1]], {v0, 1}, Evaluated -> False]
(* {v0 -> 1.49383} *)

BTW, you can also add this option to ContourPlot to suppress the Part::partd warning.

xzczd

Posted 2019-02-08T13:48:50.850

Reputation: 44 878

1(+1) Do you have a different version of the documentation? :) In other words, where to find this? – gwr – 2019-02-09T08:53:27.450

1

@gwr The Evaluated option is indeed undocumented, but it's discussed quite a bit in this site, and I just happened to read several, for example this and this. You can search Evaluated in this site to find more. :)

– xzczd – 2019-02-09T08:59:23.850

Thanks. Seems to be undocumented since Version 6, when it was introduced according to MathGroup post from 2007. Mathematica is not for the common sense fraction of users I guess. ;-) – gwr – 2019-02-09T09:07:52.853

@xzczd Thank you for your very helpful answer, I was looking in the wrong direction for something like "SymbolicProcessing" -> 0 – Ulrich Neumann – 2019-02-09T09:52:41.650

1

This can all be done with DSolve

eqns = {x''[t] + 3/5 x'[t] + x[t] == 1, x[0] == x0, x'[0] == v0};

sol = DSolve[eqns, x, t][[1]];

Verifying the solution

eqns /. sol // Simplify

(* {True, True, True} *)

For x0 == 0

sol2 = NSolve[{x[4] == 1, x[t] == 1, 0 < t < 5} /. sol /. x0 -> 0, 
  {t, v0}]

(* {{v0 -> 1.49383, t -> 0.706716}, {v0 -> 1.49383, t -> 4.}} *)

ContourPlot[Evaluate[x[4] == 1 /. sol], {x0, -.7, .7}, {v0, 0, 2}, 
 FrameLabel -> (Style[#, 12, Bold] & /@ {x0, v0}),
 Epilog -> {Red, AbsolutePointSize[4], Point[{0, v0} /. sol2[[1, 1]]]}]

enter image description here

Plot[Evaluate@{1, x[t] /. sol /. x0 -> 0 /. sol2[[1, 1]]},
 {t, 0, 5},
 Epilog -> {Red, AbsolutePointSize[4], Point[{t, 1} /. sol2]}]

enter image description here

The interval for which x[t] >= 1 is

FunctionDomain[{Sqrt[x[t] - 1] /. sol /. x0 -> 0 /. sol2[[1, 1]], 
  0 < t < 5}, t]

(* 0.706716 <= t <= 4. *)

Bob Hanlon

Posted 2019-02-08T13:48:50.850

Reputation: 95 281

Thank you for your answer. My underlying problem is much more complex and I need to evaluate the parametric function as shown. – Ulrich Neumann – 2019-02-08T21:07:44.947