Math Library Quality Check

A little program I wrote can be used to analyze the quality of the math library implementation. It simply computes the real result with a high precision and determines from this the error of the value the libm implementation produced. Very easy.

To make the time waiting for the result a bit more interesting the program displays the error graphically in an updated window. Once the program is done it writes the plot in a GIF file which then can be used for later reference.

Below you can find a few examples of outputs for a few platforms I played with. Quite interesting results as I find.

Getting the Program

The program is available for your use. Of course under my (updated) extended GPL (which excludes use of the program on Mickeysoft and Be systems).

To install the program you have to have a few more packages available:

• the GNU mp library.
• libgd. This is a graphics library which can produce GIF files. You can get it from here.
• POSIX threads library.

The program is (so far) not written with portability in mind. No autoconf script, nothing. If you run into problems, fix them. If you make generally usable changes I'd like to see them. But: no whining!

Finally, here's the package.

Results

I've so far tested only a few implementations. They include:

1. GNU libc on ix86. This test uses (only?) assembler implementations of the functions. It therefore is rather a test of the FPU.
2. GNU libc on Alpha. This tests the C implementations of the functions in glibc.
3. libffm, a high-speed math library on Alpha, derived somehow from code used at Cray. We currently use version 0.28.
4. the original libfm from Cray, on Alpha.
5. Sun's libm in Solaris 2.7.
6. HP's libm in HPPA/HP-UX 11.0.
7. DEC's libm in DU 4.0 for Alpha.
8. Irix 6.5 on MIPS.

If you have more results, send them along. This is what I have so far.

Function	ix86/glibc	alpha/glibc	alpha/libffm	alpha/libfm	sparc/sol27	hppa/hpux11	alpha/du40	mips/irix65
exp
	0.50ulp	0.82ulp	1.50ulp	2.15ulp	0.80ulp	0.51ulp	0.50ulp	0.51ulp	max
	0.02%	5.76%	15.89%	20.31%	5.54%	0.25%	0.03%	0.12%	rate
cos
	0.88ulp	0.80ulp	2.73ulp	1.89ulp	0.76ulp	0.51ulp	0.53ulp	1.43ulp	max
	0.02%	1.35%	21.42%	12.63%	1.30%	0.18%	0.53%	11.50%	rate
sin
	0.50ulp	0.79ulp	2.76ulp	1.67ulp	0.82ulp	0.56ulp	0.53ulp	1.50ulp	max
	0.01%	1.46%	18.80%	11.70%	1.43%	0.25%	0.57%	10.21%	rate
tan
	1.11ulp	0.71ulp	3.38ulp	3.22ulp	0.53ulp	0.56ulp	0.54ulp	2.40ulp	max
	0.01%	1.53%	22.70%	22.27%	0.23%	0.38%	0.38%	19.60%	rate
cosh				N/A
	0.50ulp	1.19ulp	2.16ulp	N/A	1.21ulp	0.69ulp	0.50ulp	1.15ulp	max
	0.01%	8.89%	9.12%	N/A	8.83%	2.07%	0.02%	4.17%	rate
sinh				N/A
	0.50ulp	1.71ulp	1.85ulp	N/A	1.34ulp	0.92ulp	0.51ulp	1.18ulp	max
	0.04%	13.20%	17.02%	N/A	14.03%	3.45%	0.07%	5.80%	rate
tanh				N/A
	0.50ulp	1.81ulp	3.11ulp	N/A	1.96ulp	0.75ulp	0.57ulp	1.47ulp	max
	0.03%	13.43%	23.91%	N/A	13.04%	1.94%	0.55%	5.77%	rate
sqrt
	0.50ulp	0.50ulp	0.93ulp	2.81ulp	0.50ulp	0.50ulp	0.50ulp	0.50ulp	max
	0.01%	0.00%	7.25%	38.76%	0.00%	0.00%	0.00%	0.00%	rate

My interpretation of the result for Solaris and Irix is that Sun and SGI are selling a real C and math library separately. Can anybody confirm this? If this is true it's a really disgusting move.

To Do

There are two different directions in which the package will have to evolve:

• more functions must be added to libzem. Missing are among others:

  tanh

  Tangent Hyperbolicus

  acos

  Arc Cosine

  asin

  Arc Sine

  atan

  atan2

  Arc Tangent

  cbrt

  Cubic Root

• Port the program to more platforms.

Beside this somebody should perhaps run some checks whether the libzem functions actually are acurate. I'm pretty sure they are but who knows.