Measurements.jl VS mpmath

Aside from the language features, some of the libraries in Julia make it really useful for statistical computing. One really cool library I am trying to use more and more in Julia is the Measurements library [1]. With the multiple dispatch system in Julia its super easy to integrate into most problems and can let you estimate error bounds on values programs produce. Super important for scientific applications.

I am hoping in the future that I can mix this in with some auto-diff problems to get uncertainty bounds on estimation problems with minimal fiddling with covariance matrices. Right now the performance is the only problem in integrating the library into pretty much any problem :(

[1] https://github.com/JuliaPhysics/Measurements.jl

What you've done here is tell SymPy to use extra precision for the intermediate (and final) output. This doesn't truly fix the problem of cancellation and loss of precision, but for many practical purposes it can postpone the problem long enough to give you a useful result.

Internally, SymPy uses mpmath (https://mpmath.org/) for representation of numbers to arbitrary precision. You could install and use the latter library directly, gaining extra precision without going through symbolic manipulation.

All that being said, it's still good practice to avoid loss of precision at the outset. Arbitrary-precision calculations are slow compared to hardware-native floating point operations. Using the example from mpmath's homepage in iPython:

    In [1]: import mpmath as mp; import scipy as sp; import numpy as np

Either use a fixed point system with enough precision (determined beforehand) or consider a library like https://mpmath.org.

Sure you can, check out projects for high precision numbers like https://mpmath.org/

Probably the most popular package for dealing with error propagation and arbitrary precision arithmetic in Python is mpmath, more specifically the mp.iv module. For more serious applications I'd take a look at MPFR and Arb, both in C. And there are tons of ball arithmetic and interval arithmetic libraries in Fortran.

Elaborating on /u/lanemik: if you're forced to do everything numerically and aren't able to use rationals, you can also use multiple-precision arithmetic. It's significantly slower, but it's as precise as you need it to be. Note that numpy will happily work with other objects that define arithmetic operations. I haven't messed with scipy enough to know how it does things.