I notice some inconsistencies in the fact that NonlinearModelFit works with the output of ParametricNDSolve. Here's an example that works (starting with a new kernel):

```
eqs = {a'(t) == -k1 a(t) - k2 a(t)^2,
b'(t) == k1 a(t) + k2 a(t)^2,
a(0) == a0, b(0) == b0};
fixedparams = {k1 -> 1.2, b0 -> 0};
fns = {a, b};
params = {k2, a0};
solution = ParametricNDSolve(eqs /. fixedparams, fns, {t, 0, 5}, params)
fitfn = a /. solution;
paramsForDataSet = {k2 -> 1.263, a0 -> 0.0321};
dataset = {#, ((fitfn @@ params) /. paramsForDataSet)(#) +
RandomVariate(NormalDistribution(0, 0.0002))} & /@ Range(0, 5, 0.01);
ListPlot(dataset, PlotRange -> Full)
```

```
initialGuess = {k2 -> 2.0, a0 -> 0.3};
tmp = Values@initialGuess;
Dynamic@Column({Show(ListPlot(dataset, PlotRange -> Full),
Plot((fitfn @@ tmp)(t), {t, 0, 5},
PlotRange -> Full, PlotStyle -> Red),
PlotRange -> Full, ImageSize -> Large),
ListPlot({#1, #2 - (fitfn @@ tmp)(#1)} & @@@ dataset,
PlotRange -> Full, AspectRatio -> 0.2,
ImageSize -> Large)})
```

This last bit gives me a dynamic update plot of my fit and residuals when it converges. Here is the editing procedure:

```
result = NonlinearModelFit(dataset, (fitfn @@ params)(t),
Evaluate(List @@@ initialGuess), t,
StepMonitor :> (tmp = params))
tmp = Values@result("BestFitParameters")
```

It looks great! But when I slightly complicate the model, the kernel crashes. Again from a new kernel:

```
eqs = {a'(t) == -k1 a(t) - k2 a(t)^2, b'(t) == k1 a(t) + k2 a(t)^2,
c(t) == q a(t) + r b(t), c(0) == q a0 + r b0, a(0) == a0,
b(0) == b0};
fixedparams = {k1 -> 1.2, b0 -> 0};
fns = {a, b, c};
params = {k2, a0, q, r};
solution = ParametricNDSolve(eqs /. fixedparams, fns, {t, 0, 5}, params)
fitfn = c /. solution;
paramsForDataSet = {k2 -> 1.263, a0 -> 0.0321, q -> 0.341,
r -> 0.8431};
dataset = {#, ((fitfn @@ params) /. paramsForDataSet)(#) +
RandomVariate(NormalDistribution(0, 0.0002))} & /@ Range(0, 5, 0.01);
ListPlot(dataset, PlotRange -> Full)
```

```
initialGuess = {k2 -> 2.0, a0 -> 0.3, q -> 0.32, r -> 0.88};
tmp = Values@initialGuess;
Dynamic@Column({Show(ListPlot(dataset, PlotRange -> Full),
Plot((fitfn @@ tmp)(t), {t, 0, 5}, PlotRange -> Full,
PlotStyle -> Red),
PlotRange -> Full, ImageSize -> Large),
ListPlot({#1, #2 - (fitfn @@ tmp)(#1)} & @@@ dataset,
PlotRange -> Full, AspectRatio -> 0.2,
ImageSize -> Large)})
result = NonlinearModelFit(dataset, (fitfn @@ params)(t),
Evaluate(List @@@ initialGuess), t,
StepMonitor :> (tmp = params))
tmp = Values@result("BestFitParameters")
```

The only differences are:

- add c (t) and c (0) to equations
- adding c to fns
- add q and r to params
- adding values for q and r to paramsForDataSet and initialGuess
- change fitfn in c instead of

Everything else is the same, but this time the kernel hangs. All suggestions will be welcome.

(In case there is a bug in Mathematica, I have submitted a bug report to Wolfram, but I do not want to rule out that I can do something wrong, that is That's why I'm asking the question here too.)