Skip to content
New issue

Have a question about this project? # for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “#”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? # to your account

Rounding floating point errors #1381

Open
matchatealeaf opened this issue Dec 11, 2024 · 2 comments
Open

Rounding floating point errors #1381

matchatealeaf opened this issue Dec 11, 2024 · 2 comments

Comments

@matchatealeaf
Copy link

I'm doing some calculations which result in an expression that looks like this:

1.0000000000000004b - 3.0000000000000013N*b + 3.000000000000001(N^2)*b + 1.1102230246251565e-16(N^3)*b + (6.123233995736769e-17 - 6.123233995736769e-17N + 6.123233995736769e-17(N^2) - 6.123233995736769e-17N*b + 1.8369701987210307e-16(N^2)*b - 2.4492935982947074e-16(N^3)*b + 1.2246467991473537e-16(N^4)*b)*im

As one can see, there are lot of terms that should really be zero as they are in the orders of 1e-16 despite using simplify.
Similarly, there are terms that should really have a factor of whole numbers rather than 1.000... and 3.000...

How does one simplify the expression to get rid of these float point errors?

@ChrisRackauckas
Copy link
Member

Did you use something like 1/5? Instead of using floating point numbers, just use rational numbers 1//5.

@matchatealeaf
Copy link
Author

Thanks for the reply.

I don't think using rational numbers will help me as the expression is a result of calculations involving square roots that promotes them to floats.

In this specific case, a lot of terms just happen to cancel out which gives a nice form involving integer factors, but I do not know that before hand and in different cases it might result in an irrational factor.

Currently I am using the following rewrite rule which I adapted from this discussion in Discourse.

using Symbolics.SymbolicUtils.Rewriters

rule_too_small =
    @rule ~x::(xs -> xs isa AbstractFloat && abs(round(xs) - xs) <= 1e-14) => round(~x)
rule = Prewalk(PassThrough(rule_too_small))

This seems to work but I don't know if there are any edge cases where it will fail.

It would be good if there is something like sympy's nsimplify instead.

# for free to join this conversation on GitHub. Already have an account? # to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants