From beebe at math.utah.edu  Sun May  4 17:05:15 2003
From: beebe at math.utah.edu (beebe@math.utah.edu)
Date: Sun May  4 16:05:34 2003
Subject: [Rd] R-1.7.0 build feedback: NetBSD 1.6 (PR#2837): final report
Message-ID: <200305041405.h44E5F2b004719@pubhealth.ku.dk>

I've now done two rebuilds of R-1.7.0 on NetBSD 1.6, one with the
--without-zlib configure option, and one without.  Both builds use the
recently-installed gcc-3.2.3 compiler.

As before, the one built normally gets a segment violation, whereas
the one built with the --without-zlib option works.

The odd thing is that neither uses shared libraries for zlib:

	% ldd /usr/local/lib/R/bin/R.bin
	/usr/local/lib/R/bin/R.bin:
		 -lm.0 => /usr/lib/libm387.so.0
		 -lm.0 => /usr/lib/libm.so.0
		 -lg2c.0 => /usr/local/lib/libg2c.so.0
		 -lc.12 => /usr/lib/libc.so.12
		 -lgcc_s.1 => /usr/local/lib/libgcc_s.so.1

The other build shows the same library list.

The NetBSD 1.6 system ships with these libraries:

	/usr/lib/libz.a   /usr/lib/libz.so.0    /usr/lib/libz_p.a
	/usr/lib/libz.so  /usr/lib/libz.so.0.2  /usr/lib/libz_pic.a

but the pkg_info commands on NetBSD and OpenBSD produce no output.

There seems to be no simple way to confirm the zlib version on these
systems.  Running strings and nm on libz.a produces no useful
identification.  However, /usr/include/zlib.h reports "version 1.1.4,
March 11th, 2002", which is the same as the one that I installed in
/usr/local/, and is the latest stable release (I'm on the zlib beta
testers list).  That number is confirmed by this simple test program:

	% cat show-zlib-version.c
	#include <stdio.h>
	#include <stdlib.h>
	#include <zlib.h>

	int
	main()
	{
	    (void)printf("zlibVersion() reports [%s]\n", zlibVersion());
	    return (EXIT_SUCCESS);
	}

	% cc show-zlib-version.c -lz && ./a.out
	zlibVersion() reports [1.1.4]

The bad one only segment faults after installation.  If I run it in
the build tree as bin/R, it works, and passes all but the base-Ex
test.

Brian Ripley responds to the reg-tests-1.Rout.fail file found in the
working build:

>> ...
>> The error will be in the last command line, so
>> only the last few lines of the file are relevant: to wit
>>
>> > ## pweibull(large, log=T):
>> > stopifnot(pweibull(seq(1,50,len=1001), 2,3, log = TRUE) < 0)
>> Error: pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE) < 0 is not TRUE
>>
>> for which you would need to investigate the output of
>>
>> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
>>
>> and I guess that this is an accuracy problem in the runtime (libc).
>> ...

Here is what I get in that output:

	% R

	R : Copyright 2003, The R Development Core Team
	Version 1.7.0 Under development (unstable) (2003-04-11)

	R is free software and comes with ABSOLUTELY NO WARRANTY.
	You are welcome to redistribute it under certain conditions.
	Type `license()' or `licence()' for distribution details.

	R is a collaborative project with many contributors.
	Type `contributors()' for more information.

	Type `demo()' for some demos, `help()' for on-line help, or
	`help.start()' for a HTML browser interface to help.
	Type `q()' to quit R.

	> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
	  [1] -2.252266e+00 -2.162061e+00 -2.076474e+00 -1.995124e+00 -1.917675e+00
	   [6] -1.843830e+00 -1.773331e+00 -1.705943e+00 -1.641459e+00 -1.579693e+00
	  [11] -1.520477e+00 -1.463659e+00 -1.409102e+00 -1.356679e+00 -1.306276e+00
	  [16] -1.257788e+00 -1.211118e+00 -1.166177e+00 -1.122883e+00 -1.081160e+00
	  [21] -1.040937e+00 -1.002149e+00 -9.647332e-01 -9.286334e-01 -8.937957e-01
	  [26] -8.601698e-01 -8.277083e-01 -7.963667e-01 -7.661032e-01 -7.368779e-01
	  [31] -7.086534e-01 -6.813940e-01 -6.550661e-01 -6.296376e-01 -6.050778e-01
	  [36] -5.813577e-01 -5.584494e-01 -5.363265e-01 -5.149634e-01 -4.943358e-01
	  [41] -4.744204e-01 -4.551946e-01 -4.366369e-01 -4.187266e-01 -4.014437e-01
	  [46] -3.847688e-01 -3.686833e-01 -3.531692e-01 -3.382092e-01 -3.237864e-01
	  [51] -3.098845e-01 -2.964877e-01 -2.835807e-01 -2.711486e-01 -2.591770e-01
	  [56] -2.476519e-01 -2.365597e-01 -2.258870e-01 -2.156210e-01 -2.057491e-01
	  [61] -1.962592e-01 -1.871392e-01 -1.783776e-01 -1.699631e-01 -1.618847e-01
	  [66] -1.541315e-01 -1.466932e-01 -1.395596e-01 -1.327205e-01 -1.261664e-01
	  [71] -1.198878e-01 -1.138753e-01 -1.081199e-01 -1.026130e-01 -9.734578e-02
	  [76] -9.231000e-02 -8.749747e-02 -8.290025e-02 -7.851057e-02 -7.432090e-02
	  [81] -7.032386e-02 -6.651231e-02 -6.287926e-02 -5.941794e-02 -5.612173e-02
	  [86] -5.298423e-02 -4.999919e-02 -4.716054e-02 -4.446239e-02 -4.189903e-02
	  [91] -3.946490e-02 -3.715463e-02 -3.496298e-02 -3.288491e-02 -3.091551e-02
	  [96] -2.905005e-02 -2.728393e-02 -2.561271e-02 -2.403212e-02 -2.253801e-02
	 [101] -2.112637e-02 -1.979336e-02 -1.853526e-02 -1.734848e-02 -1.622957e-02
	 [106] -1.517522e-02 -1.418223e-02 -1.324752e-02 -1.236817e-02 -1.154132e-02
	 [111] -1.076427e-02 -1.003442e-02 -9.349278e-03 -8.706452e-03 -8.103661e-03
	 [116] -7.538723e-03 -7.009554e-03 -6.514161e-03 -6.050648e-03 -5.617201e-03
	 [121] -5.212098e-03 -4.833694e-03 -4.480428e-03 -4.150814e-03 -3.843441e-03
	 [126] -3.556967e-03 -3.290122e-03 -3.041701e-03 -2.810561e-03 -2.595622e-03
	 [131] -2.395859e-03 -2.210307e-03 -2.038050e-03 -1.878228e-03 -1.730026e-03
	 [136] -1.592676e-03 -1.465456e-03 -1.347685e-03 -1.238723e-03 -1.137968e-03
	 [141] -1.044855e-03 -9.588518e-04 -8.794615e-04 -8.062168e-04 -7.386800e-04
	 [146] -6.764415e-04 -6.191181e-04 -5.663515e-04 -5.178069e-04 -4.731716e-04
	 [151] -4.321541e-04 -3.944823e-04 -3.599030e-04 -3.281801e-04 -2.990941e-04
	 [156] -2.724409e-04 -2.480308e-04 -2.256875e-04 -2.052476e-04 -1.865595e-04
	 [161] -1.694827e-04 -1.538870e-04 -1.396520e-04 -1.266662e-04 -1.148267e-04
	 [166] -1.040384e-04 -9.421341e-05 -8.527081e-05 -7.713589e-05 -6.973986e-05
	 [171] -6.301937e-05 -5.691614e-05 -5.137658e-05 -4.635146e-05 -4.179554e-05
	 [176] -3.766734e-05 -3.392878e-05 -3.054498e-05 -2.748400e-05 -2.471657e-05
	 [181] -2.221595e-05 -1.995767e-05 -1.791939e-05 -1.608070e-05 -1.442298e-05
	 [186] -1.292925e-05 -1.158404e-05 -1.037325e-05 -9.284062e-06 -8.304808e-06
	 [191] -7.424880e-06 -6.634643e-06 -5.925350e-06 -5.289063e-06 -4.718585e-06
	 [196] -4.207393e-06 -3.749581e-06 -3.339801e-06 -2.973218e-06 -2.645460e-06
	 [201] -2.352578e-06 -2.091005e-06 -1.857524e-06 -1.649233e-06 -1.463518e-06
	 [206] -1.298022e-06 -1.150627e-06 -1.019425e-06 -9.027020e-07 -7.989171e-07
	 [211] -7.066874e-07 -6.247715e-07 -5.520563e-07 -4.875440e-07 -4.303409e-07
	 [216] -3.796467e-07 -3.347456e-07 -2.949976e-07 -2.598306e-07 -2.287339e-07
	 [221] -2.012514e-07 -1.769765e-07 -1.555466e-07 -1.366387e-07 -1.199652e-07
	 [226] -1.052701e-07 -9.232584e-08 -8.093002e-08 -7.090295e-08 -6.208508e-08
	 [231] -5.433485e-08 -4.752674e-08 -4.154950e-08 -3.630461e-08 -3.170488e-08
	 [236] -2.767316e-08 -2.414125e-08 -2.104887e-08 -1.834283e-08 -1.597615e-08
	 [241] -1.390740e-08 -1.210008e-08 -1.052202e-08 -9.144876e-09 -7.943739e-09
	 [246] -6.896685e-09 -5.984448e-09 -5.190104e-09 -4.498796e-09 -3.897489e-09
	 [251] -3.374750e-09 -2.920564e-09 -2.526156e-09 -2.183845e-09 -1.886912e-09
	 [256] -1.629483e-09 -1.406425e-09 -1.213253e-09 -1.046054e-09 -9.014167e-10
	 [261] -7.763638e-10 -6.683025e-10 -5.749754e-10 -4.944174e-10 -4.249193e-10
	 [266] -3.649955e-10 -3.133551e-10 -2.688775e-10 -2.305899e-10 -1.976489e-10
	 [271] -1.693233e-10 -1.449798e-10 -1.240699e-10 -1.061191e-10 -9.071710e-11
	 [276] -7.750911e-11 -6.618894e-11 -5.649181e-11 -4.818967e-11 -4.108569e-11
	 [281] -3.501033e-11 -2.981737e-11 -2.538114e-11 -2.159339e-11 -1.836120e-11
	 [286] -1.560441e-11 -1.325440e-11 -1.125233e-11 -9.547585e-12 -8.096857e-12
	 [291] -6.862844e-12 -5.813794e-12 -4.922507e-12 -4.165557e-12 -3.523182e-12
	 [296] -2.978284e-12 -2.516320e-12 -2.124856e-12 -1.793343e-12 -1.512790e-12
	 [301] -1.275424e-12 -1.074696e-12 -9.050538e-13 -7.618350e-13 -6.409318e-13
	 [306] -5.389023e-13 -4.528600e-13 -3.803624e-13 -3.193001e-13 -2.678968e-13
	 [311] -2.247091e-13 -1.882938e-13 -1.577627e-13 -1.321165e-13 -1.105782e-13
	 [316] -9.248158e-14 -7.727152e-14 -6.461498e-14 -5.395684e-14 -4.496403e-14
	 [321] -3.752554e-14 -3.130829e-14 -2.609024e-14 -2.176037e-14 -1.809664e-14
	 [326] -1.498801e-14 -1.254552e-14 -1.043610e-14 -8.659740e-15 -7.216450e-15
	 [331] -5.995204e-15 -4.884981e-15 -4.107825e-15 -3.330669e-15 -2.775558e-15
	 [336] -2.331468e-15 -1.887379e-15 -1.554312e-15 -1.332268e-15 -1.110223e-15
	 [341] -8.881784e-16 -7.771561e-16 -5.551115e-16 -5.551115e-16 -4.440892e-16
	 [346] -3.330669e-16 -3.330669e-16 -2.220446e-16 -2.220446e-16 -1.110223e-16
	 [351] -1.110223e-16 -1.110223e-16 -1.110223e-16 -1.110223e-16 -1.110223e-16
	 [356]  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00
	...all lines are zero from here on...
	 [996]  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00
	[1001]  0.000000e+00

When I compare it with output from FreeBSD 5.0 on the same platform, I
get 4 columns instead of 5, and the numbers are nonzero right up to
the last:

	> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
	   [1]  -2.252266e+00  -2.162061e+00  -2.076474e+00  -1.995124e+00
	   [5]  -1.917675e+00  -1.843830e+00  -1.773331e+00  -1.705943e+00
	...
	 [993] -1.765317e-119 -1.028280e-119 -5.986435e-120 -3.483321e-120
	 [997] -2.025755e-120 -1.177467e-120 -6.840359e-121 -3.971708e-121
	[1001] -2.304857e-121

On OpenBSD 3.2 on the same platform (all of these systems run on top
of VMware on top of GNU/Linux 7.2 or 8.0 on Pentium III 600MHz dual
CPU servers), I get 5-column output, and trailing zeros:

	> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
	   [1] -2.252266e+00 -2.162061e+00 -2.076474e+00 -1.995124e+00 -1.917675e+00
	   [6] -1.843830e+00 -1.773331e+00 -1.705943e+00 -1.641459e+00 -1.579693e+00
	  [11] -1.520477e+00 -1.463659e+00 -1.409102e+00 -1.356679e+00 -1.306276e+00
	  [16] -1.257788e+00 -1.211118e+00 -1.166177e+00 -1.122883e+00 -1.081160e+00
	  [21] -1.040937e+00 -1.002149e+00 -9.647332e-01 -9.286334e-01 -8.937957e-01
	...

	 [351] -1.110223e-16 -1.110223e-16 -1.110223e-16 -1.110223e-16 -1.110223e-16
	 [356]  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00
	...
	 [996]  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00  0.000000e+00
	[1001]  0.000000e+00

On Solaris 9 x86 on the same platform, I get 4-column output:

	> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
	   [1]  -2.252266e+00  -2.162061e+00  -2.076474e+00  -1.995124e+00
	   [5]  -1.917675e+00  -1.843830e+00  -1.773331e+00  -1.705943e+00
	   [9]  -1.641459e+00  -1.579693e+00  -1.520477e+00  -1.463659e+00
	  [13]  -1.409102e+00  -1.356679e+00  -1.306276e+00  -1.257788e+00
	  [17]  -1.211118e+00  -1.166177e+00  -1.122883e+00  -1.081160e+00
	...
	 [353]  -8.543005e-17  -7.001648e-17  -5.735327e-17  -4.695528e-17
	 [357]  -3.842190e-17  -3.142256e-17  -2.568459e-17  -2.098321e-17
	 [361]  -1.713324e-17  -1.398220e-17  -1.140459e-17  -9.297196e-18
	...
	 [993] -1.765317e-119 -1.028280e-119 -5.986435e-120 -3.483321e-120
	 [997] -2.025755e-120 -1.177467e-120 -6.840359e-121 -3.971708e-121
	[1001] -2.304857e-121

On the underlying GNU/Linux systems, I get 4-column output:

	> pweibull(seq(1, 50, len = 1001), 2, 3, log = TRUE)
	   [1]  -2.252266e+00  -2.162061e+00  -2.076474e+00  -1.995124e+00
	   [5]  -1.917675e+00  -1.843830e+00  -1.773331e+00  -1.705943e+00
	   [9]  -1.641459e+00  -1.579693e+00  -1.520477e+00  -1.463659e+00
	  [13]  -1.409102e+00  -1.356679e+00  -1.306276e+00  -1.257788e+00
	  [17]  -1.211118e+00  -1.166177e+00  -1.122883e+00  -1.081160e+00
	...
	 [353]  -8.543005e-17  -7.001648e-17  -5.735327e-17  -4.695528e-17
	 [357]  -3.842190e-17  -3.142256e-17  -2.568459e-17  -2.098321e-17
	 [361]  -1.713324e-17  -1.398220e-17  -1.140459e-17  -9.297196e-18
	...
	 [997] -2.025755e-120 -1.177467e-120 -6.840359e-121 -3.971708e-121
	[1001] -2.304857e-121

I'd like to get to the bottom of these differences.  My numerical
calculator project, extended hoc

	http://www.math.utah.edu/hoc

has an extensive battery of validation tests, which have turned up
bugs in the math library on a half-dozen systems, including a couple
in glibc.  However, although the tests show some differences between
the four guest O/Ses, the chi-square and incomplete gamma function
tests pass.

The tests reveal that OpenBSD has these anomalies: exp(-Inf) -> NaN,
exp(Inf) -> Inf, but exp(-MAXNORMAL) -> 0 and exp(MAXNORMAL) -> +Inf,
as expected.

Examination of the R-1.7.0/src/nmath/pweibull.c file shows that
pweibull() calls pow(), log(), log1p(), and R_D_exp(), and the latter
is defined in ./src/nmath/dpq.h as

#define R_D_exp(x)       (log_p  ?  (x)   : exp(x))      /* exp(x) */

I therefore tried some experiments with hoc, using log1p() as the most
likely candidate for errors, since it is quite new, and much less used
than exp() and log().  And voil, log1p() is at least part of the
problem:

FreeBSD, Solaris 9 x86, GNU/Linux, Solaris SPARC:

	hoc64> for (k = 0; k <= 100; ++k) println k, log1p(2^-k)
	0 0.69314718055994529
	1 0.40546510810816438
	...
	50 8.8817841970012484e-16
	51 4.4408920985006252e-16
	52 2.2204460492503128e-16
	53 1.1102230246251565e-16
	...
	99 1.5777218104420236e-30
	100 7.8886090522101181e-31

NetBSD, OpenBSD:

	hoc64> for (k = 0; k <= 100; ++k) println k, log1p(2^-k)
	0 0.69314718055994529
	1 0.40546510810816438
	...
	50 8.8817841970012484e-16
	51 4.4408920985006252e-16
	52 2.2204460492503128e-16
	53 0
	54 0
	...
	99 0
	100 0

OpenBSD 3.3 was announced earlier this week, so if my colleague who
installs these systems can find the time, we may be able to test it as
well.

log1p will now get some additional tests in the hoc validation suite.

-------------------------------------------------------------------------------
- Nelson H. F. Beebe                    Tel: +1 801 581 5254                  -
- Center for Scientific Computing       FAX: +1 801 581 4148                  -
- University of Utah                    Internet e-mail: beebe@math.utah.edu  -
- Department of Mathematics, 110 LCB        beebe@acm.org  beebe@computer.org -
- 155 S 1400 E RM 233                       beebe@ieee.org                    -
- Salt Lake City, UT 84112-0090, USA    URL: http://www.math.utah.edu/~beebe  -


[Rd] R-1.7.0 build feedback: NetBSD 1.6 (PR#2837): final report
beebe@math.utah.edu