10

I've solved a differential equation using NDSolve with a large working precision prec. The output is a few interpolating functions which, when evaluated, return values with the same working precision specified before.

I'd like to export these functions to avoid rerunning NDSolve. I landed on using

Export["fun.m", f[x]]

for exporting the function and

f[x_] = Import["fun.m"]

to import the function. However, after importing, the values returned by f[x] are always given by the MachinePrecision. Taking Precision[f[x]] also returns MachinePrecision.

I've taken a look at the exported .m files and it looks like all of the data is there, so I am not sure what the issue is.

Here is a minimal example:

f = NDSolveValue[{y'[x] == y[x], y[0] == 1}, y, {x, 0, 2}, WorkingPrecision -> 20];

Export["fun.m", f];
fExported = Import["fun.m"];


f[1] // FullForm
(* 2.7182818122020232403695791388107538970319.68219897174271 *) fExported[1] // FullForm (* 2.718281812202023*)
Goofy
  • 2,767
  • 10
  • 12
MarcosMFlores
  • 315
  • 1
  • 6
  • 1
    fExported // InputForm and FilePrint["fun.m"] show the problem is probably with Import. The values at the end of "fun.m" are at the correct precision, but they are machine precision in fExported. It's either a bug or there's some option to Import you should use (don't know what that is, though). Using "fun.mx" seems to fix the problem, if that's any help. – Goofy Jan 02 '24 at 17:17
  • 2
    What seems to work is using DumpSave instead of Export. See also How do I save a variable or function definition to a file?. – Domen Jan 02 '24 at 17:23
  • 4
    Actually foo = First@Import[funFile, "HeldExpressions"] shows Import does its job right. But ReleaseHold[foo] shows that the evaluation of InterpolatingFunction converts the values to machine precision. I'd say it's a bug and should be reported to WRI. – Goofy Jan 02 '24 at 17:29
  • 2
    Please report this to WRI. – user21 Jan 02 '24 at 18:01

1 Answers1

6

This reproduces the behavior without using Export or Import:

f[1]//FullForm
(InterpolatingFunction@@f)[1]//FullForm
(*
2.7182818284590670945608961396241485466543375144201657067031`29.698970004336015
2.7182818284590673`
*)

The ".mx" file format, used by DumpSave, is a system-dependent file format that saves a binary version of the definitions of a symbol that does not have the bug above. The ".wxf" and ".wdx" formats seem to have the same evaluation of InterpolatingFunction problem as ".m".

Wolfram Language File Formats: https://reference.wolfram.com/language/guide/WolframLanguageFileFormats.html

Update: Potential workaround

Changing the interpolation method seems to avoid the bug, at a cost of some memory and a little time:

f = NDSolveValue[{y'[x] == y[x], y[0] == 1}, y, {x, 0, 2}, 
    WorkingPrecision -> 20]; // RepeatedTiming

(* {0.00133956, Null} *)


fAll = NDSolveValue[{y'[x] == y[x], y[0] == 1}, y, {x, 0, 2}, 
    WorkingPrecision -> 20, 
    InterpolationOrder -> All]; // RepeatedTiming


(* {0.0014855, Null} *)


ByteCount /@ {f, fAll}


(* {19656, 47656} *)


Export[FileNameJoin[{$TemporaryDirectory, "fun.m"}], fAll];
fAllExported = Import[%];
fAllExported[1]


(* 2.7182818308898678447 *)
Goofy
  • 2,767
  • 10
  • 12