# ParametricNDSolve Behaviour

Let's consider a set of equations with initial conditions. Since equations of my interest are complicated, I want to solve those with ParametricNDSolve.
For simplicity let's discuss linear equations though.

testEquations = {x'[t] == a, y'[t] == b, x[0] == x0, y[0] == y0};
xy = ParametricNDSolveValue[testEquations, {x, y}, {t, -10, 10}, {a, b, x0, y0}];


The result I expect is a list of 2 interpolating functions (one for x, one for y), which should return me functions of t after substituting the parameters.
Now if we substitute something

subs = {a -> 1, b -> 1, x0 -> 1, y0 -> 2};
xy[a, b, x0, y0][[1]][t] /. subs


The output is

1[t]


The same output is for xy[a, b, x0, y0][[2]][t] and

xy[a, b, x0, y0][[1]] /. subs


Returns just 1.

xy[1, 1, 1, 1][[1]][0] /. subs


Gives output

ParametricNDSolveValue::ndnum: Encountered non-numerical value for a derivative at t == 0.. >>
1[0]


Does anyone know what did I miss?
I am using Mathematica 10.3.1

• xy[a, b, x0, y0][[1]][t] gives a[t], and this is the first that is evaluated. Then you apply the rules on this output - it's about the order of evaluation (like differentation - $f'(2)$ doesn't mean to differentiate $f(2)$). Why not just Set the values of parameters? Nov 14, 2016 at 10:04
• Oh, now I see it. I didn't use Set on the parameters because later I will need to use NonlinearModelFit, which also gives an error, so I tried to find a flaw in this piece of code first. The error message is looking similar to the output which I obtained by substitution and reads as '"The function value {1.[0.],-1.+1.[1.],-1.+1.[0.],-2.+1.[1.]} is not \ a list of real numbers with dimensions {4} at {b,c,n0,p0} = \ {1.,1.,1.,2.}."' May be I will have to open another question to be more specific. Thank you for your help anyway! Nov 14, 2016 at 12:08

As mentioned by corey979, this is a matter of evaluation order. The easiest workaround is probably:

Unevaluated@xy[a, b, x0, y0][[1]][t] /. subs


or

Hold@xy[a, b, x0, y0][[1]][t] /. subs // ReleaseHold


The expression xy[a, b, x0, y0][[1]][t] returns a[t] unless all the parameters are numbers.

The easiest way to avoid it, as it seems to me, is to use ParametricNDSolveValue 2 times to obtain each of the functions of interest, instead of trying to return a list of functions.

xSol = ParametricNDSolveValue[testEquations, x, {t, -10, 10}, {a, b, x0, y0}];
ySol = ParametricNDSolveValue[testEquations, y, {t, -10, 10}, {a, b, x0, y0}];
`