Detecting Floating-Point Errors via Atomic Conditions

Setup

Install the docker. Here are some guides for install docker for Ubuntu and docker for MacOS

Clone this repo to your local workspace:

$ git clone https://github.com/FP-Analysis/atomic-condition.git

The docker folder contains a Dockerfile that can be used to build our implementation. If an error about permission denied occurs, try to run docker in root sudo docker ....

$ cd atomic-condition/docker
$ docker build -t atomic .

It may takes a few minutes (5 minutes approx) for installing necessary packages and compiling the whole GSL Library.

Usage

Run this docker's container with interactive mode, the working directory is at /atom.

$ docker run -it atomic /bin/bash

Play with the Motivation Example

The small example foo in src/targetExample.c computes the motivation example in our paper: y = (1 - cos(x)) / (x * x)

Run on the small example foo, defined in src/targetExample.c:

( In docker's container: /atom )
$ make
$ bin/gslSolver.out example

The result shows:

...
Most suspicious input first:
-1.6006806510433069e-08         Output: 4.3331212264181757e-01
1.5707963267949008e+00          Output: 4.0528473456935066e-01
End of suspicious input list.
...

Then we use mpmath to calculate the accurate result of this two inputs:

(In python3)
>>> import mpmath
>>> mpmath.mp.prec=256
>>> def foo_h(x):
...     x = mpmath.mpf(x)
...     return (1-mpmath.cos(x))/(x*x)
...
>>> y1 = foo_h(-1.6006806510433069e-08)
>>> print(y1)
0.49999999999999998932425...
>>> y2 = foo_h(1.5707963267949008e+00)
>>> print(y2)
0.40528473456935062536197...

We can see that there is a significant error when x=-1.6006806510433069e-08, the first of our reported inputs.

If you want to play with your own function, just modify the foo function in src/targetExample.c, then

$ make
$ bin/gslSolver.out example

Run on GSL Functions

$ make
$ bin/gslSolver.out gsl <function_index>

Compute Relative Error (Only support GSL functions for now)

Using the oracle from mpmath to compute the relative error:

$ make
$ bin/gslSolver.out gsl <function_index> && python3 script/oracleMpmath.py

Example: GSL function gsl_sf_lngamma has the index 73,

$ bin/gslSolver.out gsl 73 && python3 script/oracleMpmath.py
...
Function Index: 73
Max Relative Error:
  Input:  -2.457024738220797e+00
  Output: 3.552713678800501e-15
  Oracle: 4.767633887335327e-15
        Relative Error: 2.54827e-01

Reproduce the Results in Paper

Getting the relative errors of GSL functions (Table 4 and Table 5)

You can check GSL functions one by one as mentioned in above section. You can also check all of the results in single run with:

$ make
$ bin/gslSolver.out gsl all
$ python3 script/oracleMpmath.py

Getting the statistics of the relative errors (Figure 3) and the ranking results (Figure 4)

(If you have already run these 3 commands, you don't need to run them again)
$ make
$ bin/gslSolver.out gsl all
$ python3 script/oracleMpmath.py

(Please do not insert other commands between this one and the above 3 ones, otherwise you need to run them once again to get the necessary intermediate results.)
$ python3 script/reproduce.py

Data sources for LSGA and DEMC

The data of LSGA's running time are collected from the LSGA paper, which claims (the running time)

... are approximately equal to 60 seconds.

The data of DEMC's running time are collected from the DEMC's data repository.

GSL Function List and mpmath Support

GSL Function Name	Index	mpmath support?
gsl_sf_airy_Ai	0	✔️
gsl_sf_airy_Bi	1	✔️
gsl_sf_airy_Ai_scaled	2	✔️
gsl_sf_airy_Bi_scaled	3	✔️
gsl_sf_airy_Ai_deriv	4	✔️
gsl_sf_airy_Bi_deriv	5	✔️
gsl_sf_airy_Ai_deriv_scaled	6	✔️
gsl_sf_airy_Bi_deriv_scaled	7	✔️
gsl_sf_bessel_J0	8	✔️
gsl_sf_bessel_J1	9	✔️
gsl_sf_bessel_Y0	10	✔️
gsl_sf_bessel_Y1	11	✔️
gsl_sf_bessel_I0	12	✔️
gsl_sf_bessel_I1	13	✔️
gsl_sf_bessel_I0_scaled	14	✔️
gsl_sf_bessel_I1_scaled	15	✔️
gsl_sf_bessel_K0	16	✔️
gsl_sf_bessel_K1	17	✔️
gsl_sf_bessel_K0_scaled	18	✔️
gsl_sf_bessel_K1_scaled	19	✔️
gsl_sf_bessel_j0	20	✔️
gsl_sf_bessel_j1	21	✔️
gsl_sf_bessel_j2	22	✔️
gsl_sf_bessel_y0	23	✔️
gsl_sf_bessel_y1	24	✔️
gsl_sf_bessel_y2	25	✔️
gsl_sf_bessel_i0_scaled	26	✔️
gsl_sf_bessel_i1_scaled	27	✔️
gsl_sf_bessel_i2_scaled	28	✔️
gsl_sf_bessel_k0_scaled	29	✔️
gsl_sf_bessel_k1_scaled	30	✔️
gsl_sf_bessel_k2_scaled	31	✔️
gsl_sf_clausen	32	✔️
gsl_sf_dawson	33
gsl_sf_debye_1	34
gsl_sf_debye_2	35
gsl_sf_debye_3	36
gsl_sf_debye_4	37
gsl_sf_debye_5	38
gsl_sf_debye_6	39
gsl_sf_dilog	40	✔️
gsl_sf_ellint_Kcomp	41	✔️
gsl_sf_ellint_Ecomp	42	✔️
gsl_sf_ellint_Dcomp	43
gsl_sf_erfc	44	✔️
gsl_sf_log_erfc	45	✔️
gsl_sf_erf	46	✔️
gsl_sf_erf_Z	47	✔️
gsl_sf_erf_Q	48	✔️
gsl_sf_hazard	49	✔️
gsl_sf_exp	50	✔️
gsl_sf_expm1	51	✔️
gsl_sf_exprel	52	✔️
gsl_sf_exprel_2	53	✔️
gsl_sf_expint_E1	54	✔️
gsl_sf_expint_E2	55	✔️
gsl_sf_expint_E1_scaled	56	✔️
gsl_sf_expint_E2_scaled	57	✔️
gsl_sf_expint_Ei	58	✔️
gsl_sf_expint_Ei_scaled	59	✔️
gsl_sf_Shi	60	✔️
gsl_sf_Chi	61	✔️
gsl_sf_expint_3	62
gsl_sf_Si	63	✔️
gsl_sf_Ci	64	✔️
gsl_sf_atanint	65
gsl_sf_fermi_dirac_m1	66	✔️
gsl_sf_fermi_dirac_0	67	✔️
gsl_sf_fermi_dirac_1	68	✔️
gsl_sf_fermi_dirac_2	69	✔️
gsl_sf_fermi_dirac_mhalf	70	✔️
gsl_sf_fermi_dirac_half	71	✔️
gsl_sf_fermi_dirac_3half	72	✔️
gsl_sf_lngamma	73	✔️
gsl_sf_gamma	74	✔️
gsl_sf_gammastar	75
gsl_sf_gammainv	76	✔️
gsl_sf_lambert_W0	77	✔️
gsl_sf_lambert_Wm1	78	✔️
gsl_sf_legendre_P1	79	✔️
gsl_sf_legendre_P2	80	✔️
gsl_sf_legendre_P3	81	✔️
gsl_sf_legendre_Q0	82	✔️
gsl_sf_legendre_Q1	83	✔️
gsl_sf_log	84	✔️
gsl_sf_log_abs	85	✔️
gsl_sf_log_1plusx	86	✔️
gsl_sf_log_1plusx_mx	87	✔️
gsl_sf_psi	88	✔️
gsl_sf_psi_1piy	89
gsl_sf_psi_1	90	✔️
gsl_sf_sin_pi	91
gsl_sf_cos_pi	92
gsl_sf_synchrotron_1	93
gsl_sf_synchrotron_2	94	✔️
gsl_sf_transport_2	95
gsl_sf_transport_3	96
gsl_sf_transport_4	97
gsl_sf_transport_5	98
gsl_sf_sin	99	✔️
gsl_sf_cos	100	✔️
gsl_sf_sinc	101	✔️
gsl_sf_lnsinh	102	✔️
gsl_sf_lncosh	103	✔️
gsl_sf_zeta	104	✔️
gsl_sf_zetam1	105	✔️
gsl_sf_eta	106	✔️