21

Bug introduced in 8.0 and fixed in 9.0.0

I have an InterpolatingFunction based on irregularly-gridded data, like this:

int=Interpolation[{{0,0,1},{1,0,1},{0,1,1},{1,1,1},{0.5,0.5,0}},InterpolationOrder -> 1]

InterpolatingFunction[{{0.00075,0.00735},{0.00335,0.00555}}]

int[0.5,0.5]

0.

that I want to save in a file for later, using a recommendation from Szabolcs in a comment on this question:

Export["test.mmaz",Compress[int],"String"];

but importing it returns a broken function

intNew=Uncompress[Import["test.mmaz"]]

InterpolatingFunction[{{0.,1.},{0.,1.}},{4,4225,0,{5,0},{2,2},0,0,0,0,Automatic},{NDSolve`FEM`ElementMesh[{{0.,1.},{0.,1.}},{NDSolve`FEM`TriangleElement[<4>]}]},{1.,1.,1.,1.,0.},{Automatic}]

Not what I expected. Moreover, it does not evaluate with arguments:

intNew[0.5,0.5]

InterpolatingFunction[{{0.,1.},{0.,1.}},{4,4225,0,{5,0},{2,2},0,0,0,0,Automatic},{NDSolve`FEM`ElementMesh[{{0.,1.},{0.,1.}},{NDSolve`FEM`TriangleElement[<4>]}]},{1.,1.,1.,1.,0.},{Automatic}][0.5,0.5]

And yet

intNew == int

True

As far as I've seen this behaviour does not arise when using regularly-gridded data. I'm not aware of the full scope of the issue.

What breaks the InterpolatingFunction on import? Is there a fix for an existing saved function? What is the correct way to store the function safely?

Michael E2
  • 235,386
  • 17
  • 334
  • 747
JxB
  • 5,111
  • 1
  • 23
  • 40

1 Answers1

21

This is a bug in version 8 and has been fixed in the development version. For now, you have to export the data and reconstruct the interpolation once you have read in the data. What follows is way to recover your data. You should not use this on a day to day basis. The idea is to recover your data and store the data and then reconstruct the interpolation.

coords = intNew[[3, 1, 1]];
vals = Partition[intNew[[4]], 1];
data = Join[coords, vals, 2];
Interpolation[data, InterpolationOrder -> 1]

Update: Here is a slightly better fix. For large data the re-computation of the underlying mesh can be expensive. In this case, (and only in this case), you can use the following to avoid the expensive mesh creation.

mesh = intNew[[3, 1]];
vals = intNew[[4]];
iff = NDSolve`FEM`ElementInterpolation[{mesh}, vals]

Hope this helps.

  • +1, thanks. It helps to save the data. Remains a minor annoyance to have to reprocess the Interpolation between sessions when working with large irregularly gridded datasets, though. I'm still hoping someone finds a hack to store the interpolation, maybe by Blocking a key expression during the export/import... – JxB Mar 21 '12 at 17:58
  • @JxB, I doubt that is possible, since the bug was in the C code, but who knows perhaps someone has a clever idea. Why don't you send me an email of this site and we can have a chat. My email is ruebenko at the company name.com. –  Mar 21 '12 at 20:59
  • @ruebenko I just tried the workaround with NDSolve`FEM`ElementInterpolation[{mesh}, vals] and it returns Removed[$$Failure]. Are you sure it works (I use version 8.0.4)? – Alexey Popkov Aug 19 '12 at 11:35
  • @ruebenko Which versions of Mathematica are affected by this bug? – Alexey Popkov Aug 19 '12 at 14:13
  • @AlexeyPopkov, yes I am sure it works; V8. Can you post what you tried? –  Aug 20 '12 at 15:37
  • @ruebenko Try pts={{1,1,1},{2,1.2,.8},{3,1.3,.7},{3,2,1},{1.5,2,1},{1.5,1.5,1}};int=Interpolation[pts/.{a_,b_,c_}:>{{a,b},c},InterpolationOrder->1];DumpSave["int.mx",int];Get["int.mx"]; mesh=int[[3,1]]; vals=int[[4]]; iff=NDSolve`FEM`ElementInterpolation[{mesh},vals]; iff[1.1,1.5] – Alexey Popkov Aug 20 '12 at 18:15
  • 2
    @AlexeyPopkov, also fixed in the development version. For this slightly different case you can use mtemp = int[[3, 1]]; mesh = NDSolve\FEM`ElementMesh[ mtemp[[1]], {NDSolve`FEM`MeshElementCreate[ NDSolve`FEM`TriangleElement, mtemp[[2, 1, 1]]]}]` to recreate the underlying mesh. Hope this helps. –  Aug 20 '12 at 20:13