FiniteDifferences.jl
FiniteDifferences.jl copied to clipboard
Accuracy with Float32s is bad
@wesselb and I had a discussion about this a while ago, and I completely forgot to raise an issue about it. While FiniteDifferences accuracy for functions of Float64s is as you would expect, other types are a different matter. Of particular concern is Float32 as it's almost certainty the next most used type.
IIRC the issue is an issue of defaults. Specifically this one. float always yields a Float64, and eps(::Float64) is very different from eps(Float32). This changes the step-size calculation here.
If someone has time to look into this it would be greatly appreciated. I'm not entirely sure what the correct solution is, but a good approach to solving it might be
- recap how in principle to set the
epsparameter that is causing the problems (perhaps @wesselb could help here?) - figure out how to set
epsbased on this and make the appropriate changes - test that it works
If I'm not mistaken, the parameter eps is used to calculate the round-off error and should be set to the machine epsilon of the corresponding data type. The issue is that TINY is Float64, so bound always is Float64, even if the function outputs Float32s. In that case, eps is incorrectly set to the machine epsilon of Float64s, so the step size will be off.
For the particular case of Float64s and Float32s, most likely changing TINY to 1f-20 will fix most of the issues.
Might make sense to cast it explicitly wherever it's used to the element type of the thing that it's being added to. This way if (for some reason) we ever use Float16s or some extended-precision stuff we're protected.
Yes, that's a good suggestion. We should carefully check the code, to ensure that the assumption of Float64s is not also subtly baked in in other places.
Actually, it might make sense just to replace TINY with a function called add_tiny or something to prevent accidentally forgetting to cast stuff.
Yes, that's an even better solution.
Do we have a good example where the accuracy for Float32s is not what we would like it to be? I wonder if v0.12 resolves these issues.