3

I want to evaluate the solution of a system of non-linear ODEs using ParametricNDSolve. The output of ParametricNDSolve is a ParametricFunction object. Let's call it $f_\theta(t)$, where $\theta$ is a list of all parameters used in ParametricNDSolve.

For a given choice of the parameters $\hat\theta$, the ParametricFunction becomes a InterpolatingFunction object. I want to evaluate this function at multiple points $\{t_1,t_2,...,t_n\}$, which is, obtaining the value $f_{\hat\theta}(t_i)$.

Does Mathematica fully solve the system of equations each time I call $f_{\hat\theta}(t_i)$ for a different $i$? Or does it cache the InterpolatingFunction after the first call (let's say, $f_{\hat\theta}(t_1)$) and uses it to obtain the value of $f_{\hat\theta}$ at $t_2,...,t_n$?

anonymous
  • 421
  • 2
  • 9

1 Answers1

4

To extend this Q&A beyond one that is "easily found in the documentation," where some basic facts about "ParametricCaching" are given, I'll add some undocumented suboptions to Method -> {"ParametricCaching" -> {Automatic, ...opts...}}:

sol = ParametricNDSolveValue[{x''[t] + a x[t] == 0, x[0] == b, 
    x'[0] == 0}, x, {t, 0, 1000}, {a, b}, MaxSteps -> ∞];
sol[[-1, -1, -3]]
(*
{"Cache" -> True, "CacheTableLength" -> 19, "CacheTableWidth" -> 7, 
 "CacheKeyMaxBytes" -> 1000000, "CacheResultMaxBytes" -> 1000000, 
 "KeyComparison" -> None, "ResultComparison" -> LessEqual}
*)

For instance, we can use the timings to verify when a cached solution is reused. This confirms that the "TableCache*" suboptions seem to set the dimensions of the cache table:

sol1 = ParametricNDSolveValue[{x''[t] + a x[t] == 0, x[0] == b, 
    x'[0] == 0}, x, {t, 0, 1000}, {a, b}, MaxSteps -> ∞, 
   Method -> {"ParametricCaching" -> {Automatic, 
       "CacheTableLength" -> 1, "CacheTableWidth" -> 1}}]; 
sol1[1, 1][1] // AbsoluteTiming
sol1[1, 1][2] // AbsoluteTiming
sol1[1, 2][2] // AbsoluteTiming
sol1[1, 1][2] // AbsoluteTiming
(*
  {0.017357, 0.540302}
  {0.000014, -0.416147}
  {0.020204, -0.832294}
  {0.01881, -0.416147}   <-- N.B. Recomputed
*)

sol2 = ParametricNDSolveValue[{x''[t] + a x[t] == 0, x[0] == b, 
   x'[0] == 0}, x, {t, 0, 1000}, {a, b}, MaxSteps -> ∞, 
  Method -> {"ParametricCaching" -> {Automatic, 
      "CacheTableLength" -> 1, "CacheTableWidth" -> 2}}]; 
sol2[1, 1][1] // AbsoluteTiming
sol2[1, 1][2] // AbsoluteTiming
sol2[1, 2][2] // AbsoluteTiming
sol2[1, 1][2] // AbsoluteTiming
(*
  {0.019945, 0.540302}
  {0.000017, -0.416147}
  {0.01882, -0.832294}
  {0.000022, -0.416147}  <-- N.B. Reused
*)

The memory-related options should be self-explanatory and useful. We can test "CacheResultMaxBytes" but I'm not sure how to test "CacheKeyMaxBytes".

(* bytes of a result *)
sol1[1, 1] // ByteCount
(*  305704  *)

sol3 = ParametricNDSolveValue[{x''[t] + a x[t] == 0, x[0] == b, 
   x'[0] == 0}, x, {t, 0, 1000}, {a, b}, MaxSteps -> ∞, 
  Method -> {"ParametricCaching" -> {Automatic, 
      "CacheResultMaxBytes" -> 400000}}]; 
sol3[1, 1][1] // AbsoluteTiming
sol3[1, 1][2] // AbsoluteTiming
(*
  {0.015709, 0.540302}
  {0.000013, -0.416147}  <-- N.B. Reused
*)


sol4 = ParametricNDSolveValue[{x''[t] + a x[t] == 0, x[0] == b, 
   x'[0] == 0}, x, {t, 0, 1000}, {a, b}, MaxSteps -> ∞, 
  Method -> {"ParametricCaching" -> {Automatic, 
      "CacheResultMaxBytes" -> 1000}}];
sol4[1, 1][1] // AbsoluteTiming
sol4[1, 1][2] // AbsoluteTiming
(*
  {0.018927, 0.540302}
  {0.022862, -0.416147}  <-- N.B. Recomputed
*)

I don't about keys and their function in ParametricNDSolve caching, nor how "ResultComparison" is used.

(Note it's important to rerun ParametricNDSolveValue each time you want to test to reset the cache table.)

Michael E2
  • 235,386
  • 17
  • 334
  • 747