A peculiar behavior of ParametricFunction

I am trying to use ParametricNDSolve to find a desired value for a parameter in an ODE, and the output seems to misbehave.

First, consider the following example given in the documentation:

sol1 = ParametricNDSolve[{x''[t] - x'[t] + x[t] == 0, x[0] == 1, x'[0] == s}, x, {t, 0, 30}, {s}]
(* ::Output:: *) (*
{x -> ParametricFunction[ <> ]}
*)
root1 = FindRoot[Evaluate[x[s][10] /. sol1], {s, 6}]
(* ::Output:: *) (*
{s -> 1.40296}
*)


Secondly, my equation with the Automatically chosen method:

sol2 = ParametricNDSolve[{X'[t] == t + (Sqrt[3] / Pi) Log[a/(1 - a)] Abs[t X[t]] - X[t], X[0] == 1}, {X}, {t, 0, 2}, a]
(* ::Output:: *) (*
{x -> ParametricFunction[ <> ]}
*)

root2 = FindRoot[Evaluate[X[a][2] /. sol2] == 3, {a, .5}]
(* ::Output:: *) (*
[!] InterpolatingFunction: Input {2} lies outside the range of data... (InterpolatingFunction::dmval)
{a -> 0.812787}
*)


In which I don't understand how the range is exceeded. Furthermore, specifying a method in the call yields:

sol3 = = ParametricNDSolve[{X'[t] == t + (Sqrt[3] / Pi) Log[a/(1 - a)] Abs[t X[t]] - X[t], X[0] == 1}, {X}, {t, 0, 2}, a, Method -> {"TimeIntegration" -> "ExplicitEuler"}]
(* ::Output:: *) (*
{x -> ParametricFunction[ <> ]}
*)

X[.7][2] /. sol3
(* ::Output:: *) (*
2.01701
*)

root3 = FindRoot[Evaluate[X[a][2] /. sol3] == 3, {a, .5}]
(* ::Output:: *) (*
[!] ParametricNDSolve: Encountered invalid NDSolveSensitivityMethod method data object at point t=0. (ParametricNDSolve::mdata)
[!] FindRoot: At {a}={0.5}, function value {-3+ParametricFunction[1,InternalBag[<1>],1,1,False,{{a$71966},<<5>>,{0}},{NDSolvebase$71973,NDSolveNDSolveParametricFunction[0,{ParametricNDSolve,InternalBag[<2>],None,ParametricNDSolve},{{{<<9>>},{<<9>>}},{0,{<<3>>},{<<2>>},{<<3>>}},None,{{<<7>>},{<<1>>},None,{}}},{X},<<4>>,{Cache->True,CacheTableLength->19,CacheTableWidth->7,CacheKeyMaxBytes->1000000,CacheResultMaxBytes->1000000,KeyComparison->None,ResultComparison->LessEqual},{},<<1>>]}][0.5][2]} is not a {1} dimensional list of numbers  (FindRoot::nlnum)
{a -> 0.5}
*)

X[.7][2] /. sol3
(* ::Output:: *) (*
ParametricFunction[ <> ][0.7][2]
*)


Note that t=2 is indeed in the interpolation range, and the call to the ParametricFunction instance fails after the unsuccessful FindRoot run.

How to understand this behavior, and is there a current workaround if it's a bug?

Related:

• In the first link, @Carl Woll's answer suggests that there is indeed something wrong with ParametricFunction, however I am not sure if the behavior in the current post is of the same origin. Root solving is a large numeric topic, and I haven't finished all the tech notes and SE QA's just yet; do be kind enough and point me to anything that I am supposed to read. Commented Dec 14, 2022 at 11:42
• I think "Findroot" makes a first guess using something like the tangent. This can lead to a bad first guess. Therefore, try a starting value closer to the searched root. E.g.: FindRoot[Evaluate[X[a][2] /. sol1] == 3, {a, .7}] Commented Dec 14, 2022 at 12:17
• @DanielHuber This indeed yields a good root and prevents the ParametricFunction from being spoiled Commented Dec 14, 2022 at 13:06
• To get a good starting value use Plot, i.e., Plot[(X[a][2] - 3) /. sol3, {a, 0, 1}] indicates the solution is in the vicinity of 0.8 Alternatively, avoid derivatives by using "FindRoot[lhs == rhs, {x, Subscript[x, 0], Subscript[x, 1]}] searches for a solution using Subscript[x, 0] and Subscript[x, 1] as the first two values of x, avoiding the use of derivatives." Commented Dec 14, 2022 at 15:27

Evaluated -> False

Clear[sol2, root2];
sol2 = ParametricNDSolve[{X'[t] ==
t + (Sqrt[3]/Pi) Log[a/(1 - a)] Abs[t X[t]] - X[t],
X[0] == 1}, {X}, {t, 0, 2}, a]
root2 = FindRoot[(X[a][2] /. sol2) == 3, {a, .5}, Evaluated -> False]


{a -> 0.812787}

• I find this indeed works. Would you explain what the option means and does? Thanks in advance Commented Dec 14, 2022 at 15:27
• Clear[sol]; sol = ParametricNDSolveValue[{X'[t] == t + (Sqrt[3]/Pi) Log[a/(1 - a)] Abs[t X[t]] - X[t], X[0] == 1}, {X}, {t, 0, 2}, {a}]; FindRoot[Indexed[sol[a], 1][2] == 3, {a, .5}] Commented Mar 12, 2023 at 6:23

Without Method -> {"TimeIntegration" -> "ExplicitEuler"} try more direct solution using ParametricNDSolveValue(only return X[2])

X3 = ParametricNDSolveValue[{X'[t] == t + (Sqrt[3]/Pi) Log[a/(1 - a)] Abs[t X[t]]- X[t], X[0] == 1},
X[2] , {t, 0, 2}, a  ]
root3 = FindRoot[ X3[a]   == 3, {a, 1/2 }]
(* {a -> 0.812787} *)


NMinimize evaluates without message " Input value {2} lies outside..."

NMinimize[ (X3[a]  - 3)^2, a]
(*{6.10426*10^-16, {a -> 0.812787}}  *)

• This seems to be my sol2? Sorry I pasted the wrong output though. Commented Dec 14, 2022 at 13:02
• No it isn't your sol2 Commented Dec 14, 2022 at 13:03
• Oh sorry I glanced over the second line... Commented Dec 14, 2022 at 13:08

and is there a current workaround

Some strange evaluation order? Here is a workaround

sol2 = ParametricNDSolve[{X'[t] == t + (Sqrt[3]/Pi) Log[a/(1 - a)] Abs[t X[t]] - X[t],
X[0] == 1}, {X}, {t, 0, 2}, a];
f[a_] := X[a] /. sol2
root2 = FindRoot[Evaluate[f[a][2] /. sol2] == 3, {a, .5}]


No beep and no warnings and no range exceeded messages.

V 13.1 on windows 10

Screen shot. No messages.

• v13.1 here too. Weirdly enough this doesn't seems to fix it for me. I'll clear the kernel and try again. Commented Dec 14, 2022 at 12:57
• @Gravifer no messages/warnings for me on V 13.1 on windows 10. Yes, make sure to Quit the kernel. Added screen shot Commented Dec 14, 2022 at 13:21
• I can now replicate your result. However, what dictates the result for me is /.sol2 in the FindRoot; with this replacement no exception is thrown. I don't fully understand why this is the case yet. Commented Dec 14, 2022 at 15:24