Plotting solution of ParametricNDSolve not possible for needed parameter range

I wish to numerically solve a parametric decay model. The system of equations is fairly simple. But plotting the solution for a wider parameter range is not possible (see code).

Mathematica throws 2 error messages:

1. ParametricNDSolve::nderr: Error test failure at t == 0.`; unable to continue.
2. InterpolatingFunction::dmval: Input value {2.04286*10^-11} lies outside the range of data in the interpolating function. Extrapolation will be used.

I tried to play around with the accuracy goal because of the small parameters (especialle k3), but somehow a smaller accuracy goal seems to work better (but not good enough) which to me is very strange...

Any ideas?

Here is my compilable code:

ClearAll[k1, k2, k3];

tmin = 0;
tmax = 1000*10^-9;

k1try = 1*10^5;
k2try = 1*10^-10;
k3try = 1*10^-28;

c = 6.2*10^16;

fitmodel = ParametricNDSolve[
{n'[t] == -k1*n[t] - k2*n[t]^2 - k3*n[t]^3,
PL[t] == (n[t]^2 )/(c^2),
n == c,
PL == 1},
PL, {t, tmin, tmax}, {k1, k2, k3}, AccuracyGoal -> 40]

Plot[Evaluate[
Table[PL[k1try, k2try, k3][t] /. fitmodel, {k3, 1*10^-28, 10*10^-26,
10*10^-28}]], {t, tmin, tmax}, PlotRange -> All]
• That seems to work but takes a long time: MaxSteps -> 10000000, MaxStepSize -> 10^-12 as options for ParametricNDSolve Jun 21 '18 at 10:06

tmin = 0;
tmax = 10^(-6);

k1try = 1*10^5;
k2try = 1*10^-10;
k3try = 1*10^-28;
c = 6.2*10^16;

fitmodel =
ParametricNDSolve[{n'[t] == -k1*n[t] - k2*c*n[t]^2 - k3*c*c*n[t]^3,
n == 1}, n, {t, tmin, tmax}, {k1, k2, k3}]

Plot[Evaluate[
Table[n[k1try, k2try, k3][t]^2 /. fitmodel, {k3, 1*10^-28,
10*10^-26, 10*10^-28}]], {t, tmin, tmax}, PlotRange -> All] • How did you arrive at this simplification? Jun 23 '18 at 15:49
• Made a replacement n[t]->c*n[t] Jun 25 '18 at 13:10
• I did not see the squared when you create the table. All clear now. Thanks! Jun 26 '18 at 8:27