summary refs log tree commit diff
path: root/sysdeps/libm-ieee754
diff options
context:
space:
mode:
authorUlrich Drepper <drepper@redhat.com>1997-06-19 19:11:22 +0000
committerUlrich Drepper <drepper@redhat.com>1997-06-19 19:11:22 +0000
commitd705269e75ff9846000abea981774b381b076482 (patch)
tree7bd717a8489ad113cb0d10b7628a88d6a55cb836 /sysdeps/libm-ieee754
parente7f1f046f51d95d5d05ada3abd3f2e88219118e9 (diff)
downloadglibc-d705269e75ff9846000abea981774b381b076482.tar.gz
glibc-d705269e75ff9846000abea981774b381b076482.tar.xz
glibc-d705269e75ff9846000abea981774b381b076482.zip
1997-06-19 19:38  Ulrich Drepper  <drepper@cygnus.com>

	* features.h: Define __STDC_IEC_559__ and _STDC_IEC_559_COMPLEX__.

	* elf/dl-minimal.c (__dcgettext): Remove assertion.

	* inet/rcmd.c: Correct a few typos.  Reported by Eric Troan.

	* manual/Makefile (distribute): Add dir.
	* manual/dir: New file.

	* math/Makefile (libm-support): Rename s_rinttol, s_rinttoll,
	s_roundtol, and s_roundtoll to s_lrint, s_llrint, s_lround,
	and s_llround respectively.
	(libm-calls): Add e_gamma_r.
	* math/libm-test.c (check_int_exc): New function.
	(signbit_test): Remove test for sign of NaN.
	(gamma_test): Clear exception after test of existence.
	Correct tests and and epsilons.
	(lgamma_test): Likewise.
	(ilogb_test): Correct all tests.
	(scalb_test): Rewrite.
	(rinttol_test): Rename to lrint_test and correct tests.
	(rinttoll_test): Likewise.
	(roundtol_test): Likewise.
	(roundtoll_test): Likewise.
	(main): Call lrint/lround functions instead of rinttol/roundtol.
	* math/math.h: Change prototypes for rinttol/roundtol.
	* math/mathcalls.h: Rearrange prototypes according to ISO C9X draft.
	* sysdeps/generic/mathbits.h: Define FP_ILOGB0 and FP_ILOGBNAN.
	* sysdeps/i386/mathbits.h: Likewise.
	* sysdeps/libm-i387/e_scalb.S: Handle special cases correctly.
	* sysdeps/libm-i387/e_scalbf.S: Likewise.
	* sysdeps/libm-i387/e_scalbl.S: Likewise.
	* sysdeps/libm-i387/s_asinh.S: Handle -inf correctly.
	* sysdeps/libm-i387/s_asinhf.S: Likewise.
	* sysdeps/libm-i387/s_asinhl.S: Likewise.
	* sysdeps/libm-i387/s_ilogb.S: Optimize.
	* sysdeps/libm-i387/s_ilogbf.S: Likewise.
	* sysdeps/libm-i387/s_ilogbl.S: Likewise.
	* sysdeps/libm-i387/s_rinttol.S: Rename to...
	* sysdeps/libm-i387/s_lrint.S: ...this.
	* sysdeps/libm-i387/s_rinttoll.S: Rename to...
	* sysdeps/libm-i387/s_llrint.S: ...this.
	* sysdeps/libm-i387/s_remquo.S: Correctly set sign of remainder.
	* sysdeps/libm-i387/s_remquof.S: Likewise.
	* sysdeps/libm-i387/s_remquol.S: Likewise.
	* sysdeps/libm-i387/e_gamma_r.c: New file.  Implementation of gamma
	function according to ISO C.
	* sysdeps/libm-i387/e_gammaf_r.c: New file.
	* sysdeps/libm-i387/e_gammal_r.c: New file.
	* sysdeps/libm-i387/e_lgamma_r.c: Don't let optimize compile the
	generation of exceptions away.
	* sysdeps/libm-i387/e_lgammaf_r.c: Likewise.
	* sysdeps/libm-i387/k_standard.c: Correct return value for infinity
	points of gamma function when not SVID mode.
	* sysdeps/libm-i387/s_rinttoll.c: Renamed to...
	* sysdeps/libm-i387/s_llrint.c: ...this.
	* sysdeps/libm-i387/s_rinttol.c: Renamed to...
	* sysdeps/libm-i387/s_lrint.c: ...this.
	* sysdeps/libm-i387/s_roundtoll.c: Renamed to...
	* sysdeps/libm-i387/s_llround.c: ...this.
	* sysdeps/libm-i387/s_roundtol.c: Renamed to...
	* sysdeps/libm-i387/s_lround.c: ..this.
	* sysdeps/libm-i387/s_scalbn.c: Change second parameter according to
	ISO C.
	* sysdeps/libm-i387/s_scalbnf.c: Likewise.
	* sysdeps/libm-i387/s_scalbnl.c: Likewise.
	* sysdeps/libm-i387/w_gamma.c: Call __ieee754_gamma_r if library
	mode is not _SVID_.
	* sysdeps/libm-i387/w_gammaf.c: Likewise.
	* sysdeps/libm-i387/w_gammal.c: Likewise.
	* sysdeps/m68k/fpu/__math.h: Rename __rinttol to __lrint and
	rinttol to lrint.
	* sysdeps/m68k/fpu/s_rinttol.c: Renamed to...
	* sysdeps/m68k/fpu/s_lrint.c: ...this.
	* sysdeps/m68k/fpu/s_rinttoll.c: Renamed to...
	* sysdeps/m68k/fpu/s_llrint.c: ...this.

	* md5-crypt/Makefile: Link md5test program with md5.o.

	* stdio-common/temptest.c: Don't use __stdio_gen_tempname which
	is not exported by the libc.so.

	* stdio-common/vfscanf.c: Correct scanning of strings after last
	change.

	* sysdeps/unix/sysv/linux/i386/sysdep.S: Use .comm to define errno.

1997-06-19 07:37  H.J. Lu  <hjl@gnu.ai.mit.edu>

	* time/tzfile.c (__tzfile_read): Store getc () return in int.

1997-06-13  Miles Bader  <miles@gnu.ai.mit.edu>

	* argp/argp-parse.c (argp_version_parser): Include `(PROGRAM ERROR)' in
	the no-version error text to indicate that something's fucked.
	[!_] (N_): New macro.
	(argp_default_options, argp_version_options): Wrap doc strings in N_().

1997-06-12  Miles Bader  <miles@gnu.ai.mit.edu>

	* argp/argp-parse.c (parser_parse_opt): Detect and report unhandled
	options here.
	(parser_parse_arg): Handle ARGP_KEY_ARGS here.
	Adjust NEXT pointer back if we fail to parse anything.
	(parser_parse_next): Simplify arg code.  Leave state NEXT frobbing
	to parser_parse_arg.

1997-06-11  Miles Bader  <miles@gnu.ai.mit.edu>

	* argp/argp.h (ARGP_KEY_ARGS, ARGP_KEY_FINI): New macros.
	* argp/argp-parse.c (parser_finalize): Do another pass over the
	parsers with ARGP_KEY_FINI.

1997-06-18  Miles Bader  <miles@gnu.ai.mit.edu>

	* string/Makefile (routines): Add argz-replace.

1997-06-16 00:16  Miles Bader  <miles@gnu.ai.mit.edu>

	* manual/string.texi (Argz Functions): Document argz_replace.

1997-06-12  Miles Bader  <miles@gnu.ai.mit.edu>

	* string/argz.h (__argz_replace, argz_replace): New declarations.
	* string/argz-replace.c: New file.

1997-06-11  Miles Bader  <miles@gnu.ai.mit.edu>

	* manual/argp.texi (Argp Special Keys): Document ARGP_KEY_ARGS.

1997-06-16 23:08  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* Makeconfig (libc-map): Remove definition.
	* Makerules (libc-map): Define it here, using the full name.
	(load-map-file): Remove case for empty $(..).
	($(common-objpfx)libc.so): Revert last change.

1997-06-17 22:18  Mark Kettenis  <kettenis@phys.uva.nl>

	* login/programs/utmpd.c (main): Improve signal handling.

	* login/programs/request.c (do_setutent, do_updwtmp):
	Allow arbitrary length filenames.
	* login/programs/utmpd.h (setutent_request, updwtmp_request):
	Get rid of fixed length file field.
	* login/utmp_daemon.c (do_setutent, do_updwtmp):
	Allow arbitrary	length filenames.

	* login/programs/request.c (do_pututline):
	Don't fail if connection->position is -1 on entry.

1997-06-15 16:32  Mark Kettenis  <kettenis@phys.uva.nl>

	* login/utmp_file.c (updwtmp_file): Use the same method for
	appending an entry as in pututline_file.

1997-06-11 18:59  Mark Kettenis  <kettenis@phys.uva.nl>

	* login/utmp_file.c (getutent_r_file):
	Use read lock instead of write lock.
	(getutline_r_file, internal_getut_r): Lock utmp file.
	(updwtmp_file): Use fcntl to lock file instead of flock.

1997-06-18 00:11  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* stdio-common/vfscanf.c (inchar, ungetc): Don't count EOF as
	character read in/put back.
	* stdio-common/tstscanf.c: Add test case for this.

1997-06-17 22:17  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* libc.map: Add more symbols.

1997-06-18 12:01  Ulrich Drepper  <drepper@cygnus.com>

	* manual/Makefile (install): Make sure `dir' file exists if we use
	install-info.

1997-06-17 19:32  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* manual/Makefile (info): Depend on dir-add.info.

1997-06-17 17:19  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* elf/ldd.bash.in: Fix spacing in message.

1997-06-17 14:28  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* elf/do-rel.h (elf_dynamic_do_rel): Always use version
	information if available.

1997-06-17 11:34  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sunrpc/Makefile ($(rpcsvc:%.x=$(objpfx)rpcsvc/%.h)): Make
	command non-empty to force make to recheck modification time.
	($(rpcsvc:%.x=$(objpfx)x%.c)): Likewise.

1997-06-17 00:26  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sysdeps/stub/e_acoshl.c: Set errno to ENOSYS.
	* sysdeps/stub/e_acosl.c: Likewise.
	* sysdeps/stub/e_asinl.c: Likewise.
	* sysdeps/stub/e_atan2l.c: Likewise.
	* sysdeps/stub/e_expl.c: Likewise.
	* sysdeps/stub/e_fmodl.c: Likewise.
	* sysdeps/stub/e_j0l.c: Likewise.
	* sysdeps/stub/e_j1l.c: Likewise.
	* sysdeps/stub/e_jnl.c: Likewise.
	* sysdeps/stub/e_lgammal_r.c: Likewise.
	* sysdeps/stub/e_log10l.c: Likewise.
	* sysdeps/stub/e_logl.c: Likewise.
	* sysdeps/stub/e_powl.c: Likewise.
	* sysdeps/stub/e_rem_pio2l.c: Likewise.
	* sysdeps/stub/e_sqrtl.c: Likewise.
	* sysdeps/stub/k_cosl.c: Likewise.
	* sysdeps/stub/k_rem_pio2l.c: Likewise.
	* sysdeps/stub/k_sinl.c: Likewise.
	* sysdeps/stub/k_tanl.c: Likewise.
	* sysdeps/stub/s_atanl.c: Likewise.
	* sysdeps/stub/s_erfl.c: Likewise.
	* sysdeps/stub/s_exp2.c: Likewise.
	* sysdeps/stub/s_exp2f.c: Likewise.
	* sysdeps/stub/s_exp2l.c: Likewise.
	* sysdeps/stub/s_expm1l.c: Likewise.
	* sysdeps/stub/s_log1pl.c: Likewise.
	* sysdeps/stub/s_log2l.c: Likewise.

1997-06-18 11:46  Ulrich Drepper  <drepper@cygnus.com>

	* sysdeps/libm-ieee754/w_gamma.c: If _LIB_VERSION is _SVID_ compute
	result as before last change.
	* sysdeps/libm-ieee754/w_gammaf.c: Likewise.
	* sysdeps/libm-ieee754/w_gammal.c: Likewise.

1997-06-16 23:37  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sysdeps/libm-ieee754/s_remquo.c: Fix off-by-one when computing
	quotient.
	* sysdeps/libm-ieee754/s_remquof.c: Likewise.
	* sysdeps/libm-ieee754/s_remquol.c: Likewise.

	* sysdeps/m68k/fpu/s_remquo.c: Remove FIXME and special case for
	quotient.

	* sysdeps/libm-ieee754/w_gamma.c: Add missing call to exp
	function.  Don't use global signgam.
	* sysdeps/libm-ieee754/w_gammaf.c: Likewise.
	* sysdeps/libm-ieee754/w_gammal.c: Likewise.

	* math/Makefile (libm-calls): Remove w_gamma_r.
	* sysdeps/libm-ieee754/w_gamma_r.c: Remove file.
	* sysdeps/libm-ieee754/w_gammaf_r.c: Remove file.
	* sysdeps/libm-ieee754/w_gammal_r.c: Remove file.

	* math/libm-test.c (atanh_test): Declare x only if needed.
	(signbit_test): Fix typo.
	(gamma_test): Check whether function is implemented.  Add
	epsilons.
	(lgamma_test): Likewise.
	(fmod_test): Add epsilons.
	(exp2_test): Use right function for existence test.

1997-06-07 09:20  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sysdeps/powerpc/Dist: Add fenv_const.c.
	* sysdeps/unix/sysv/linux/Dist: Add net/if_slip.h.
	* sysdeps/unix/sysv/linux/powerpc/Dist: Add init-first.h and
	syscall.h.
	* sysdeps/unix/sysv/linux/sparc/Dist: Add init-first.h.
	* string/Makefile (distribute): Add tst-svc.expect.
	* nis/Makefile (distribute): Add nis_intern.h and Banner.
	* elf/Makefile (distribute): Add dl-hash.h.
	* Rules (subdir_echo-distinfo): Add headers from $(distribute).

	* login/Makefile (others): Add utmpdump.
	* sysdeps/unix/sysv/linux/Makefile (sysdep_headers): Add
	net/if_slip.h.

	* manual/Makefile (dir-add.texi): Also look in indirectly included
	files.

1997-06-16 23:15  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sysdeps/unix/sysv/linux/configure.in: Remove commands to
	generate stdio_lim.h.
	* sysdeps/unix/sysv/linux/mk-stdiolim.c: Remove.
	* sysdeps/unix/sysv/linux/Makefile: Add rules to generate
	stdio_lim.h here.
	(common-generated): Add generated files.
	(inhibit-stdio_lim): Define.
	* sysdeps/posix/Makefile [$(inhibit-stdio_lim)=yes]: Disable rules
	to generate stdio_lim.h.

	* sysdeps/unix/sysv/linux/Makefile: Suppress inclusion of
	dependecy files if no_deps is set.
	($(objpfx)syscall-%.d): Add header file as target to dependency
	generation.

1997-06-14 19:19  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* sysdeps/unix/sysv/linux/arm/Dist: New file.

1997-06-14 17:59  Andreas Schwab  <schwab@issan.informatik.uni-dortmund.de>

	* extra-lib.mk (others): Depend on versioned shared library, not
	the unversioned one.
	* Makerules (build-shlib): Don't make the version link here.
	($(common-objpfx)libc.so$(libc.so-version)): New rule for libc
	version link.

1997-06-16 03:07  Ulrich Drepper  <drepper@cygnus.com>

	* sysdeps/libm-ieee754/k_standard.c: Undo change of Tue Aug 6
	01:13:56 1996.

	* argp/argp-help.c (argp_args_usage): Supply correct argp to
	filter_doc.
	* argp/argp-help.c (hol_add_cluster): Initialize CL->depth.
	* argp/argp-help.c (_help): Supply STATE to argp_args_usage.
	* argp/argp.h (ARGP_KEY_HELP_ARGS_DOC): New macro.
	* argp/argp-fmtstream.c: Add casts to prevent warnings.
	* argp/argp.h (OPTION_NO_USAGE): New macro.
	* argp/argp-help.c (usage_long_opt, usage_argful_short_opt,
	* argp/argp-fmtstream.c (__argp_fmtstream_update): Account for case
	* argp/argp-help.c <stddef.h>: New include.
	* argp/argp.h (argp_state_help, __argp_state_help, argp_usage,
	* argp/argp.h (argp_program_bug_address): Make const.
	* argp/argp-parse.c (argp_default_parser): Set STATE->name for OPT_PROGNAME.
	* argp/argp-help.c (__argp_error, __argp_failure, __argp_state_help):
	* argp/argp-parse.c (parser_init): Set PARSER->state.flags.
Diffstat (limited to 'sysdeps/libm-ieee754')
-rw-r--r--sysdeps/libm-ieee754/e_gamma_r.c47
-rw-r--r--sysdeps/libm-ieee754/e_gammaf_r.c46
-rw-r--r--sysdeps/libm-ieee754/e_gammal_r.c46
-rw-r--r--sysdeps/libm-ieee754/e_lgamma_r.c40
-rw-r--r--sysdeps/libm-ieee754/e_lgammaf_r.c22
-rw-r--r--sysdeps/libm-ieee754/s_llrint.c (renamed from sysdeps/libm-ieee754/s_rinttoll.c)6
-rw-r--r--sysdeps/libm-ieee754/s_llround.c (renamed from sysdeps/libm-ieee754/s_roundtoll.c)6
-rw-r--r--sysdeps/libm-ieee754/s_lrint.c (renamed from sysdeps/libm-ieee754/s_rinttol.c)6
-rw-r--r--sysdeps/libm-ieee754/s_lround.c (renamed from sysdeps/libm-ieee754/s_roundtol.c)6
-rw-r--r--sysdeps/libm-ieee754/s_remquo.c2
-rw-r--r--sysdeps/libm-ieee754/s_remquof.c2
-rw-r--r--sysdeps/libm-ieee754/s_remquol.c2
-rw-r--r--sysdeps/libm-ieee754/s_scalbn.c13
-rw-r--r--sysdeps/libm-ieee754/s_scalbnf.c20
-rw-r--r--sysdeps/libm-ieee754/s_scalbnl.c34
-rw-r--r--sysdeps/libm-ieee754/w_gamma.c22
-rw-r--r--sysdeps/libm-ieee754/w_gamma_r.c51
-rw-r--r--sysdeps/libm-ieee754/w_gammaf.c21
-rw-r--r--sysdeps/libm-ieee754/w_gammaf_r.c52
-rw-r--r--sysdeps/libm-ieee754/w_gammal.c21
-rw-r--r--sysdeps/libm-ieee754/w_gammal_r.c52
21 files changed, 261 insertions, 256 deletions
diff --git a/sysdeps/libm-ieee754/e_gamma_r.c b/sysdeps/libm-ieee754/e_gamma_r.c
new file mode 100644
index 0000000000..eb80809843
--- /dev/null
+++ b/sysdeps/libm-ieee754/e_gamma_r.c
@@ -0,0 +1,47 @@
+/* Implementation of gamma function according to ISO C.
+   Copyright (C) 1997 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Library General Public License as
+   published by the Free Software Foundation; either version 2 of the
+   License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Library General Public License for more details.
+
+   You should have received a copy of the GNU Library General Public
+   License along with the GNU C Library; see the file COPYING.LIB.  If not,
+   write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#include <math.h>
+#include <math_private.h>
+
+
+double
+__ieee754_gamma_r (double x, int *signgamp)
+{
+  /* We don't have a real gamma implementation now.  We'll use lgamma
+     and the exp function.  But due to the required boundary
+     conditions we must check some values separately.  */
+  int32_t hx;
+  u_int32_t lx;
+
+  EXTRACT_WORDS (hx, lx, x);
+
+  if ((hx & 0x7fffffff | lx) == 0)
+    /* Return value for x == 0 is NaN with invalid exception.  */
+    return x / x;
+  if (hx < 0 && (u_int32_t) hx < 0xfff00000 && __rint (x) == x)
+    {
+      /* Return value for integer x < 0 is NaN with invalid exception.  */
+      return (x - x) / (x - x);
+    }
+
+  /* XXX FIXME.  */
+  return __ieee754_exp (__ieee754_lgamma_r (x, signgamp));
+}
diff --git a/sysdeps/libm-ieee754/e_gammaf_r.c b/sysdeps/libm-ieee754/e_gammaf_r.c
new file mode 100644
index 0000000000..bf79ffeb57
--- /dev/null
+++ b/sysdeps/libm-ieee754/e_gammaf_r.c
@@ -0,0 +1,46 @@
+/* Implementation of gamma function according to ISO C.
+   Copyright (C) 1997 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Library General Public License as
+   published by the Free Software Foundation; either version 2 of the
+   License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Library General Public License for more details.
+
+   You should have received a copy of the GNU Library General Public
+   License along with the GNU C Library; see the file COPYING.LIB.  If not,
+   write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#include <math.h>
+#include <math_private.h>
+
+
+float
+__ieee754_gammaf_r (float x, int *signgamp)
+{
+  /* We don't have a real gamma implementation now.  We'll use lgamma
+     and the exp function.  But due to the required boundary
+     conditions we must check some values separately.  */
+  int32_t hx;
+
+  GET_FLOAT_WORD (hx, x);
+
+  if ((hx & 0x7fffffff) == 0)
+    /* Return value for x == 0 is NaN with invalid exception.  */
+    return x / x;
+  if (hx < 0 && (u_int32_t) hx < 0xff800000 && __rintf (x) == x)
+    {
+      /* Return value for integer x < 0 is NaN with invalid exception.  */
+      return (x - x) / (x - x);
+    }
+
+  /* XXX FIXME.  */
+  return __ieee754_expf (__ieee754_lgammaf_r (x, signgamp));
+}
diff --git a/sysdeps/libm-ieee754/e_gammal_r.c b/sysdeps/libm-ieee754/e_gammal_r.c
new file mode 100644
index 0000000000..4fd49886d1
--- /dev/null
+++ b/sysdeps/libm-ieee754/e_gammal_r.c
@@ -0,0 +1,46 @@
+/* Implementation of gamma function according to ISO C.
+   Copyright (C) 1997 Free Software Foundation, Inc.
+   This file is part of the GNU C Library.
+   Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+   The GNU C Library is free software; you can redistribute it and/or
+   modify it under the terms of the GNU Library General Public License as
+   published by the Free Software Foundation; either version 2 of the
+   License, or (at your option) any later version.
+
+   The GNU C Library is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   Library General Public License for more details.
+
+   You should have received a copy of the GNU Library General Public
+   License along with the GNU C Library; see the file COPYING.LIB.  If not,
+   write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
+   Boston, MA 02111-1307, USA.  */
+
+#include <math.h>
+#include <math_private.h>
+
+
+long double
+__ieee754_gammal_r (long double x, int *signgamp)
+{
+  /* We don't have a real gamma implementation now.  We'll use lgamma
+     and the exp function.  But due to the required boundary
+     conditions we must check some values separately.  */
+  u_int32_t es, hx, lx;
+
+  GET_LDOUBLE_WORDS (es, hx, lx, x);
+
+  if (((es & 0x7fff) | hx | lx) == 0)
+    /* Return value for x == 0 is NaN with invalid exception.  */
+    return x / x;
+  if ((hx & 0x8000) != 0 && (hx & 0x7fff) != 0x7fff && __rintl (x) == x)
+    {
+      /* Return value for integer x < 0 is NaN with invalid exception.  */
+      return (x - x) / (x - x);
+    }
+
+  /* XXX FIXME.  */
+  return __ieee754_expl (__ieee754_lgammal_r (x, signgamp));
+}
diff --git a/sysdeps/libm-ieee754/e_lgamma_r.c b/sysdeps/libm-ieee754/e_lgamma_r.c
index 1be4ddad8f..92e9556568 100644
--- a/sysdeps/libm-ieee754/e_lgamma_r.c
+++ b/sysdeps/libm-ieee754/e_lgamma_r.c
@@ -5,7 +5,7 @@
  *
  * Developed at SunPro, a Sun Microsystems, Inc. business.
  * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice 
+ * software is freely granted, provided that this notice
  * is preserved.
  * ====================================================
  */
@@ -15,12 +15,12 @@ static char rcsid[] = "$NetBSD: e_lgamma_r.c,v 1.7 1995/05/10 20:45:42 jtc Exp $
 #endif
 
 /* __ieee754_lgamma_r(x, signgamp)
- * Reentrant version of the logarithm of the Gamma function 
- * with user provide pointer for the sign of Gamma(x). 
+ * Reentrant version of the logarithm of the Gamma function
+ * with user provide pointer for the sign of Gamma(x).
  *
  * Method:
  *   1. Argument Reduction for 0 < x <= 8
- * 	Since gamma(1+s)=s*gamma(s), for x in [0,8], we may 
+ * 	Since gamma(1+s)=s*gamma(s), for x in [0,8], we may
  * 	reduce x to a number in [1.5,2.5] by
  * 		lgamma(1+s) = log(s) + lgamma(s)
  *	for example,
@@ -58,36 +58,36 @@ static char rcsid[] = "$NetBSD: e_lgamma_r.c,v 1.7 1995/05/10 20:45:42 jtc Exp $
  *	by
  *	  			    3       5             11
  *		w = w0 + w1*z + w2*z  + w3*z  + ... + w6*z
- *	where 
+ *	where
  *		|w - f(z)| < 2**-58.74
- *		
+ *
  *   4. For negative x, since (G is gamma function)
  *		-x*G(-x)*G(x) = pi/sin(pi*x),
  * 	we have
  * 		G(x) = pi/(sin(pi*x)*(-x)*G(-x))
  *	since G(-x) is positive, sign(G(x)) = sign(sin(pi*x)) for x<0
- *	Hence, for x<0, signgam = sign(sin(pi*x)) and 
+ *	Hence, for x<0, signgam = sign(sin(pi*x)) and
  *		lgamma(x) = log(|Gamma(x)|)
  *			  = log(pi/(|x*sin(pi*x)|)) - lgamma(-x);
- *	Note: one should avoid compute pi*(-x) directly in the 
+ *	Note: one should avoid compute pi*(-x) directly in the
  *	      computation of sin(pi*(-x)).
- *		
+ *
  *   5. Special Cases
  *		lgamma(2+s) ~ s*(1-Euler) for tiny s
  *		lgamma(1)=lgamma(2)=0
  *		lgamma(x) ~ -log(x) for tiny x
  *		lgamma(0) = lgamma(inf) = inf
  *	 	lgamma(-integer) = +-inf
- *	
+ *
  */
 
 #include "math.h"
 #include "math_private.h"
 
 #ifdef __STDC__
-static const double 
+static const double
 #else
-static double 
+static double
 #endif
 two52=  4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
 half=  5.00000000000000000000e-01, /* 0x3FE00000, 0x00000000 */
@@ -200,9 +200,9 @@ static double zero=  0.00000000000000000000e+00;
         }
 	switch (n) {
 	    case 0:   y =  __kernel_sin(pi*y,zero,0); break;
-	    case 1:   
+	    case 1:
 	    case 2:   y =  __kernel_cos(pi*(0.5-y),zero); break;
-	    case 3:  
+	    case 3:
 	    case 4:   y =  __kernel_sin(pi*(one-y),zero,0); break;
 	    case 5:
 	    case 6:   y = -__kernel_cos(pi*(y-1.5),zero); break;
@@ -226,9 +226,11 @@ static double zero=  0.00000000000000000000e+00;
 
     /* purge off +-inf, NaN, +-0, and negative arguments */
 	*signgamp = 1;
+	if ((unsigned int) hx==0xfff00000&&lx==0)
+	  return x-x;
 	ix = hx&0x7fffffff;
 	if(ix>=0x7ff00000) return x*x;
-	if((ix|lx)==0) return one/zero;
+	if((ix|lx)==0) return one/fabs(x);
 	if(ix<0x3b900000) {	/* |x|<2**-70, return -log(|x|) */
 	    if(hx<0) {
 	        *signgamp = -1;
@@ -237,9 +239,9 @@ static double zero=  0.00000000000000000000e+00;
 	}
 	if(hx<0) {
 	    if(ix>=0x43300000) 	/* |x|>=2**52, must be -integer */
-		return one/zero;
+		return x/zero;
 	    t = sin_pi(x);
-	    if(t==zero) return one/zero; /* -integer */
+	    if(t==zero) return one/fabsf(t); /* -integer */
 	    nadj = __ieee754_log(pi/fabs(t*x));
 	    if(t<zero) *signgamp = -1;
 	    x = -x;
@@ -275,7 +277,7 @@ static double zero=  0.00000000000000000000e+00;
 		p3 = t2+w*(t5+w*(t8+w*(t11+w*t14)));
 		p  = z*p1-(tt-w*(p2+y*p3));
 		r += (tf + p); break;
-	      case 2:	
+	      case 2:
 		p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*u5)))));
 		p2 = one+y*(v1+y*(v2+y*(v3+y*(v4+y*v5))));
 		r += (-0.5*y + p1/p2);
@@ -304,7 +306,7 @@ static double zero=  0.00000000000000000000e+00;
 	    y = z*z;
 	    w = w0+z*(w1+y*(w2+y*(w3+y*(w4+y*(w5+y*w6)))));
 	    r = (x-half)*(t-one)+w;
-	} else 
+	} else
     /* 2**58 <= x <= inf */
 	    r =  x*(__ieee754_log(x)-one);
 	if(hx<0) r = nadj - r;
diff --git a/sysdeps/libm-ieee754/e_lgammaf_r.c b/sysdeps/libm-ieee754/e_lgammaf_r.c
index 1d0122dbb2..f744d5320e 100644
--- a/sysdeps/libm-ieee754/e_lgammaf_r.c
+++ b/sysdeps/libm-ieee754/e_lgammaf_r.c
@@ -8,7 +8,7 @@
  *
  * Developed at SunPro, a Sun Microsystems, Inc. business.
  * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice 
+ * software is freely granted, provided that this notice
  * is preserved.
  * ====================================================
  */
@@ -21,9 +21,9 @@ static char rcsid[] = "$NetBSD: e_lgammaf_r.c,v 1.3 1995/05/10 20:45:47 jtc Exp
 #include "math_private.h"
 
 #ifdef __STDC__
-static const float 
+static const float
 #else
-static float 
+static float
 #endif
 two23=  8.3886080000e+06, /* 0x4b000000 */
 half=  5.0000000000e-01, /* 0x3f000000 */
@@ -136,9 +136,9 @@ static float zero=  0.0000000000e+00;
         }
 	switch (n) {
 	    case 0:   y =  __kernel_sinf(pi*y,zero,0); break;
-	    case 1:   
+	    case 1:
 	    case 2:   y =  __kernel_cosf(pi*((float)0.5-y),zero); break;
-	    case 3:  
+	    case 3:
 	    case 4:   y =  __kernel_sinf(pi*(one-y),zero,0); break;
 	    case 5:
 	    case 6:   y = -__kernel_cosf(pi*(y-(float)1.5),zero); break;
@@ -162,9 +162,11 @@ static float zero=  0.0000000000e+00;
 
     /* purge off +-inf, NaN, +-0, and negative arguments */
 	*signgamp = 1;
+	if ((unsigned int)hx==0xff800000)
+	  return x-x;
 	ix = hx&0x7fffffff;
 	if(ix>=0x7f800000) return x*x;
-	if(ix==0) return one/zero;
+	if(ix==0) return one/fabsf(x);
 	if(ix<0x1c800000) {	/* |x|<2**-70, return -log(|x|) */
 	    if(hx<0) {
 	        *signgamp = -1;
@@ -173,9 +175,9 @@ static float zero=  0.0000000000e+00;
 	}
 	if(hx<0) {
 	    if(ix>=0x4b000000) 	/* |x|>=2**23, must be -integer */
-		return one/zero;
+		return x/zero;
 	    t = sin_pif(x);
-	    if(t==zero) return one/zero; /* -integer */
+	    if(t==zero) return one/fabsf(t); /* -integer */
 	    nadj = __ieee754_logf(pi/fabsf(t*x));
 	    if(t<zero) *signgamp = -1;
 	    x = -x;
@@ -211,7 +213,7 @@ static float zero=  0.0000000000e+00;
 		p3 = t2+w*(t5+w*(t8+w*(t11+w*t14)));
 		p  = z*p1-(tt-w*(p2+y*p3));
 		r += (tf + p); break;
-	      case 2:	
+	      case 2:
 		p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*u5)))));
 		p2 = one+y*(v1+y*(v2+y*(v3+y*(v4+y*v5))));
 		r += (-(float)0.5*y + p1/p2);
@@ -240,7 +242,7 @@ static float zero=  0.0000000000e+00;
 	    y = z*z;
 	    w = w0+z*(w1+y*(w2+y*(w3+y*(w4+y*(w5+y*w6)))));
 	    r = (x-half)*(t-one)+w;
-	} else 
+	} else
     /* 2**58 <= x <= inf */
 	    r =  x*(__ieee754_logf(x)-one);
 	if(hx<0) r = nadj - r;
diff --git a/sysdeps/libm-ieee754/s_rinttoll.c b/sysdeps/libm-ieee754/s_llrint.c
index b2fccd17b7..faae106ece 100644
--- a/sysdeps/libm-ieee754/s_rinttoll.c
+++ b/sysdeps/libm-ieee754/s_llrint.c
@@ -33,7 +33,7 @@ static long double two52[2] =
 
 
 long long int
-__rinttoll (long double x)
+__llrint (long double x)
 {
   int32_t j0,sx;
   u_int32_t i0, i1, i;
@@ -141,7 +141,7 @@ static long double two63[2] =
 
 
 long long int
-__rinttoll (long double x)
+__llrint (long double x)
 {
   int32_t se,j0,sx;
   u_int32_t i0, i1, i;
@@ -239,4 +239,4 @@ __rinttoll (long double x)
 }
 #endif
 
-weak_alias (__rinttoll, rinttoll)
+weak_alias (__llrint, llrint)
diff --git a/sysdeps/libm-ieee754/s_roundtoll.c b/sysdeps/libm-ieee754/s_llround.c
index 49167236a6..aee0e31fc0 100644
--- a/sysdeps/libm-ieee754/s_roundtoll.c
+++ b/sysdeps/libm-ieee754/s_llround.c
@@ -27,7 +27,7 @@
 /* The `long double' is in fact the IEEE `double' type.  */
 
 long long int
-__roundtoll (long double x)
+__llround (long double x)
 {
   int32_t j0;
   u_int32_t i1, i0;
@@ -103,7 +103,7 @@ __roundtoll (long double x)
 }
 #else
 long long int
-__roundtoll (long double x)
+__llround (long double x)
 {
   int32_t j0;
   u_int32_t se, i1, i0;
@@ -176,4 +176,4 @@ __roundtoll (long double x)
   return se & 0x8000 ? -result : result;
 }
 #endif
-weak_alias (__roundtoll, roundtoll)
+weak_alias (__llround, llround)
diff --git a/sysdeps/libm-ieee754/s_rinttol.c b/sysdeps/libm-ieee754/s_lrint.c
index 257bdc744c..6779f974bd 100644
--- a/sysdeps/libm-ieee754/s_rinttol.c
+++ b/sysdeps/libm-ieee754/s_lrint.c
@@ -33,7 +33,7 @@ static long double two52[2] =
 
 
 long int
-__rinttol (long double x)
+__lrint (long double x)
 {
   int32_t j0,sx;
   u_int32_t i0,i1,i;
@@ -141,7 +141,7 @@ static long double two63[2] =
 
 
 long int
-__rinttol (long double x)
+__lrint (long double x)
 {
   int32_t se,j0,sx;
   u_int32_t i0,i1,i;
@@ -239,4 +239,4 @@ __rinttol (long double x)
 }
 #endif
 
-weak_alias (__rinttol, rinttol)
+weak_alias (__lrint, lrint)
diff --git a/sysdeps/libm-ieee754/s_roundtol.c b/sysdeps/libm-ieee754/s_lround.c
index bc0ceae0f8..0f91280a79 100644
--- a/sysdeps/libm-ieee754/s_roundtol.c
+++ b/sysdeps/libm-ieee754/s_lround.c
@@ -27,7 +27,7 @@
 /* The `long double' is in fact the IEEE `double' type.  */
 
 long int
-__roundtol (long double x)
+__lround (long double x)
 {
   int32_t j0;
   u_int32_t i1, i0;
@@ -101,7 +101,7 @@ __roundtol (long double x)
 }
 #else
 long int
-__roundtol (long double x)
+__lround (long double x)
 {
   int32_t j0;
   u_int32_t se, i1, i0;
@@ -174,4 +174,4 @@ __roundtol (long double x)
   return se & 0x8000 ? -result : result;
 }
 #endif
-weak_alias (__roundtol, roundtol)
+weak_alias (__lround, lround)
diff --git a/sysdeps/libm-ieee754/s_remquo.c b/sysdeps/libm-ieee754/s_remquo.c
index f0a32130a2..39f6ceab2e 100644
--- a/sysdeps/libm-ieee754/s_remquo.c
+++ b/sysdeps/libm-ieee754/s_remquo.c
@@ -82,6 +82,7 @@ __remquo (double x, double y, int *quo)
       if (x + x > y)
 	{
 	  x -= y;
+	  ++cquo;
 	  if (x + x >= y)
 	    {
 	      x -= y;
@@ -95,6 +96,7 @@ __remquo (double x, double y, int *quo)
       if (x > y_half)
 	{
 	  x -= y;
+	  ++cquo;
 	  if (x >= y_half)
 	    {
 	      x -= y;
diff --git a/sysdeps/libm-ieee754/s_remquof.c b/sysdeps/libm-ieee754/s_remquof.c
index 6fa02e47b3..b3870f4b02 100644
--- a/sysdeps/libm-ieee754/s_remquof.c
+++ b/sysdeps/libm-ieee754/s_remquof.c
@@ -81,6 +81,7 @@ __remquof (float x, float y, int *quo)
       if (x + x > y)
 	{
 	  x -= y;
+	  ++cquo;
 	  if (x + x >= y)
 	    {
 	      x -= y;
@@ -94,6 +95,7 @@ __remquof (float x, float y, int *quo)
       if (x > y_half)
 	{
 	  x -= y;
+	  ++cquo;
 	  if (x >= y_half)
 	    {
 	      x -= y;
diff --git a/sysdeps/libm-ieee754/s_remquol.c b/sysdeps/libm-ieee754/s_remquol.c
index 9ef424901b..b7835e6dc7 100644
--- a/sysdeps/libm-ieee754/s_remquol.c
+++ b/sysdeps/libm-ieee754/s_remquol.c
@@ -82,6 +82,7 @@ __remquol (long double x, long double p, int *quo)
       if (x + x > p)
 	{
 	  x -= p;
+	  ++cquo;
 	  if (x + x >= p)
 	    {
 	      x -= p;
@@ -95,6 +96,7 @@ __remquol (long double x, long double p, int *quo)
       if (x > p_half)
 	{
 	  x -= p;
+	  ++cquo;
 	  if (x >= p_half)
 	    {
 	      x -= p;
diff --git a/sysdeps/libm-ieee754/s_scalbn.c b/sysdeps/libm-ieee754/s_scalbn.c
index 439b966e69..0ce0ffdf2f 100644
--- a/sysdeps/libm-ieee754/s_scalbn.c
+++ b/sysdeps/libm-ieee754/s_scalbn.c
@@ -35,10 +35,10 @@ huge   = 1.0e+300,
 tiny   = 1.0e-300;
 
 #ifdef __STDC__
-	double __scalbn (double x, int n)
+	double __scalbn (double x, long int n)
 #else
 	double __scalbn (x,n)
-	double x; int n;
+	double x; long int n;
 #endif
 {
 	int32_t k,hx,lx;
@@ -49,17 +49,16 @@ tiny   = 1.0e-300;
 	    x *= two54;
 	    GET_HIGH_WORD(hx,x);
 	    k = ((hx&0x7ff00000)>>20) - 54;
-            if (n< -50000) return tiny*x; 	/*underflow*/
 	    }
         if (k==0x7ff) return x+x;		/* NaN or Inf */
         k = k+n;
-        if (k >  0x7fe) return huge*__copysign(huge,x); /* overflow  */
+        if (n> 50000 || k >  0x7fe)
+	  return huge*__copysign(huge,x); /* overflow  */
+	if (n< -50000) return tiny*__copysign(tiny,x); /*underflow*/
         if (k > 0) 				/* normal result */
 	    {SET_HIGH_WORD(x,(hx&0x800fffff)|(k<<20)); return x;}
         if (k <= -54)
-            if (n > 50000) 	/* in case integer overflow in n+k */
-		return huge*__copysign(huge,x);	/*overflow*/
-	    else return tiny*__copysign(tiny,x); 	/*underflow*/
+	  return tiny*__copysign(tiny,x); 	/*underflow*/
         k += 54;				/* subnormal result */
 	SET_HIGH_WORD(x,(hx&0x800fffff)|(k<<20));
         return x*twom54;
diff --git a/sysdeps/libm-ieee754/s_scalbnf.c b/sysdeps/libm-ieee754/s_scalbnf.c
index 3a83e54c56..4799c825c9 100644
--- a/sysdeps/libm-ieee754/s_scalbnf.c
+++ b/sysdeps/libm-ieee754/s_scalbnf.c
@@ -8,7 +8,7 @@
  *
  * Developed at SunPro, a Sun Microsystems, Inc. business.
  * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice 
+ * software is freely granted, provided that this notice
  * is preserved.
  * ====================================================
  */
@@ -31,10 +31,10 @@ huge   = 1.0e+30,
 tiny   = 1.0e-30;
 
 #ifdef __STDC__
-	float __scalbnf (float x, int n)
+	float __scalbnf (float x, long int n)
 #else
 	float __scalbnf (x,n)
-	float x; int n;
+	float x; long int n;
 #endif
 {
 	int32_t k,ix;
@@ -44,18 +44,18 @@ tiny   = 1.0e-30;
             if ((ix&0x7fffffff)==0) return x; /* +-0 */
 	    x *= two25;
 	    GET_FLOAT_WORD(ix,x);
-	    k = ((ix&0x7f800000)>>23) - 25; 
-            if (n< -50000) return tiny*x; 	/*underflow*/
+	    k = ((ix&0x7f800000)>>23) - 25;
 	    }
         if (k==0xff) return x+x;		/* NaN or Inf */
-        k = k+n; 
-        if (k >  0xfe) return huge*copysignf(huge,x); /* overflow  */
+        k = k+n;
+        if (n> 50000 || k >  0xfe)
+	  return huge*copysignf(huge,x); /* overflow  */
+	if (n< -50000)
+	  return tiny*copysignf(tiny,x);	/*underflow*/
         if (k > 0) 				/* normal result */
 	    {SET_FLOAT_WORD(x,(ix&0x807fffff)|(k<<23)); return x;}
         if (k <= -25)
-            if (n > 50000) 	/* in case integer overflow in n+k */
-		return huge*copysignf(huge,x);	/*overflow*/
-	    else return tiny*copysignf(tiny,x);	/*underflow*/
+	    return tiny*copysignf(tiny,x);	/*underflow*/
         k += 25;				/* subnormal result */
 	SET_FLOAT_WORD(x,(ix&0x807fffff)|(k<<23));
         return x*twom25;
diff --git a/sysdeps/libm-ieee754/s_scalbnl.c b/sysdeps/libm-ieee754/s_scalbnl.c
index 6825e9c7b3..3e80d85136 100644
--- a/sysdeps/libm-ieee754/s_scalbnl.c
+++ b/sysdeps/libm-ieee754/s_scalbnl.c
@@ -33,16 +33,16 @@ static const long double
 #else
 static long double
 #endif
-two54   =  1.80143985094819840000e+16, /* 0x4035, 0x00000000, 0x00000000 */
-twom54  =  5.55111512312578270212e-17, /* 0x3FC9, 0x00000000, 0x00000000 */
+two63   =  4.50359962737049600000e+15,
+twom63  =  1.08420217248550443400e-19;
 huge   = 1.0e+4900L,
 tiny   = 1.0e-4900L;
 
 #ifdef __STDC__
-	long double __scalbnl (long double x, int n)
+	long double __scalbnl (long double x, long int n)
 #else
 	long double __scalbnl (x,n)
-	long double x; int n;
+	long double x; long int n;
 #endif
 {
 	int32_t k,es,hx,lx;
@@ -50,22 +50,22 @@ tiny   = 1.0e-4900L;
         k = es&0x7fff;				/* extract exponent */
         if (k==0) {				/* 0 or subnormal x */
             if ((lx|(hx&0x7fffffff))==0) return x; /* +-0 */
-	    x *= two54;
-	    GET_HIGH_WORD(hx,x);
-	    k = ((hx&0x7ff00000)>>20) - 54;
-            if (n< -50000) return tiny*x; 	/*underflow*/
+	    x *= two63;
+	    GET_LDOUBLE_EXP(es,x);
+	    k = (hx&0x7fff) - 63;
 	    }
-        if (k==0x7ff) return x+x;		/* NaN or Inf */
+        if (k==0x7fff) return x+x;		/* NaN or Inf */
         k = k+n;
-        if (k >  0x7fe) return huge*__copysign(huge,x); /* overflow  */
+        if (n> 50000 || k > 0x7ffe)
+	  return huge*__copysignl(huge,x); /* overflow  */
+	if (n< -50000)
+	  return tiny*__copysignl(tiny,x);
         if (k > 0) 				/* normal result */
-	    {SET_HIGH_WORD(x,(hx&0x800fffff)|(k<<20)); return x;}
-        if (k <= -54)
-            if (n > 50000) 	/* in case integer overflow in n+k */
-		return huge*__copysign(huge,x);	/*overflow*/
-	    else return tiny*__copysign(tiny,x); 	/*underflow*/
+	    {SET_LDOUBLE_EXP(x,(es&0x8000)|k); return x;}
+        if (k <= -63)
+	    return tiny*__copysignl(tiny,x); 	/*underflow*/
         k += 54;				/* subnormal result */
-	SET_HIGH_WORD(x,(hx&0x800fffff)|(k<<20));
-        return x*twom54;
+	SET_LDOUBLE_EXP(x,(es&0x8000)|k);
+        return x*twom63;
 }
 weak_alias (__scalbnl, scalbnl)
diff --git a/sysdeps/libm-ieee754/w_gamma.c b/sysdeps/libm-ieee754/w_gamma.c
index 49e4bcf3c3..87a3408c3f 100644
--- a/sysdeps/libm-ieee754/w_gamma.c
+++ b/sysdeps/libm-ieee754/w_gamma.c
@@ -15,16 +15,13 @@ static char rcsid[] = "$NetBSD: w_gamma.c,v 1.7 1995/11/20 22:06:43 jtc Exp $";
 #endif
 
 /* double gamma(double x)
- * Return the logarithm of the Gamma function of x.
- *
- * Method: call gamma_r
+ * Return  the logarithm of the Gamma function of x or the Gamma function of x,
+ * depending on the library mode.
  */
 
 #include "math.h"
 #include "math_private.h"
 
-extern int signgam;
-
 #ifdef __STDC__
 	double __gamma(double x)
 #else
@@ -32,12 +29,19 @@ extern int signgam;
 	double x;
 #endif
 {
+	int signgam;
+        double y;
+	if (_LIB_VERSION == _SVID_)
+	  y = __ieee754_lgamma_r(x,&signgam);
+	else
+	  {
+	    y = __ieee754_gamma_r(x,&signgam);
+	    if (signgam < 0) y = -y;
 #ifdef _IEEE_LIBM
-	return __ieee754_lgamma_r(x,&signgam);
+	    return y;
 #else
-        double y;
-        y = __ieee754_lgamma_r(x,&signgam);
-        if(_LIB_VERSION == _IEEE_) return y;
+	    if(_LIB_VERSION == _IEEE_) return y;
+	  }
         if(!__finite(y)&&__finite(x)) {
             if(__floor(x)==x&&x<=0.0)
                 return __kernel_standard(x,x,41); /* gamma pole */
diff --git a/sysdeps/libm-ieee754/w_gamma_r.c b/sysdeps/libm-ieee754/w_gamma_r.c
deleted file mode 100644
index f9efc8cf2e..0000000000
--- a/sysdeps/libm-ieee754/w_gamma_r.c
+++ /dev/null
@@ -1,51 +0,0 @@
-/* @(#)wr_gamma.c 5.1 93/09/24 */
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * Developed at SunPro, a Sun Microsystems, Inc. business.
- * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice
- * is preserved.
- * ====================================================
- */
-
-#if defined(LIBM_SCCS) && !defined(lint)
-static char rcsid[] = "$NetBSD: w_gamma_r.c,v 1.7 1995/11/20 22:06:45 jtc Exp $";
-#endif
-
-/*
- * wrapper double gamma_r(double x, int *signgamp)
- */
-
-#include "math.h"
-#include "math_private.h"
-
-
-#ifdef __STDC__
-	double __gamma_r(double x, int *signgamp) /* wrapper lgamma_r */
-#else
-	double __gamma_r(x,signgamp)              /* wrapper lgamma_r */
-        double x; int *signgamp;
-#endif
-{
-#ifdef _IEEE_LIBM
-	return __ieee754_lgamma_r(x,signgamp);
-#else
-        double y;
-        y = __ieee754_lgamma_r(x,signgamp);
-        if(_LIB_VERSION == _IEEE_) return y;
-        if(!__finite(y)&&__finite(x)) {
-            if(__floor(x)==x&&x<=0.0)
-                return __kernel_standard(x,x,41); /* gamma pole */
-            else
-                return __kernel_standard(x,x,40); /* gamma overflow */
-        } else
-            return y;
-#endif
-}
-weak_alias (__gamma_r, gamma_r)
-#ifdef NO_LONG_DOUBLE
-strong_alias (__gamma_r, __gammal_r)
-weak_alias (__gamma_r, gammal_r)
-#endif
diff --git a/sysdeps/libm-ieee754/w_gammaf.c b/sysdeps/libm-ieee754/w_gammaf.c
index c2d21daf47..ecde183824 100644
--- a/sysdeps/libm-ieee754/w_gammaf.c
+++ b/sysdeps/libm-ieee754/w_gammaf.c
@@ -8,7 +8,7 @@
  *
  * Developed at SunPro, a Sun Microsystems, Inc. business.
  * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice 
+ * software is freely granted, provided that this notice
  * is preserved.
  * ====================================================
  */
@@ -20,8 +20,6 @@ static char rcsid[] = "$NetBSD: w_gammaf.c,v 1.4 1995/11/20 22:06:48 jtc Exp $";
 #include "math.h"
 #include "math_private.h"
 
-extern int signgam;
-
 #ifdef __STDC__
 	float __gammaf(float x)
 #else
@@ -29,12 +27,19 @@ extern int signgam;
 	float x;
 #endif
 {
+	int signgam;
+        float y;
+	if (_LIB_VERSION == _SVID_)
+	  y = __ieee754_lgammaf_r(x,&signgam);
+	else
+	  {
+	    y = __ieee754_gammaf_r(x,&signgam);
+	    if (signgam < 0) y = -y;
 #ifdef _IEEE_LIBM
-	return __ieee754_lgammaf_r(x,&signgam);
+	    return y;
 #else
-        float y;
-        y = __ieee754_lgammaf_r(x,&signgam);
-        if(_LIB_VERSION == _IEEE_) return y;
+	    if(_LIB_VERSION == _IEEE_) return y;
+	  }
         if(!__finitef(y)&&__finitef(x)) {
             if(__floorf(x)==x&&x<=(float)0.0)
 	        /* gammaf pole */
@@ -45,5 +50,5 @@ extern int signgam;
         } else
             return y;
 #endif
-}             
+}
 weak_alias (__gammaf, gammaf)
diff --git a/sysdeps/libm-ieee754/w_gammaf_r.c b/sysdeps/libm-ieee754/w_gammaf_r.c
deleted file mode 100644
index b15d5e84bb..0000000000
--- a/sysdeps/libm-ieee754/w_gammaf_r.c
+++ /dev/null
@@ -1,52 +0,0 @@
-/* w_gammaf_r.c -- float version of w_gamma_r.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * Developed at SunPro, a Sun Microsystems, Inc. business.
- * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice 
- * is preserved.
- * ====================================================
- */
-
-#if defined(LIBM_SCCS) && !defined(lint)
-static char rcsid[] = "$NetBSD: w_gammaf_r.c,v 1.4 1995/11/20 22:06:50 jtc Exp $";
-#endif
-
-/* 
- * wrapper float gammaf_r(float x, int *signgamp)
- */
-
-#include "math.h"
-#include "math_private.h"
-
-
-#ifdef __STDC__
-	float __gammaf_r(float x, int *signgamp) /* wrapper lgammaf_r */
-#else
-	float __gammaf_r(x,signgamp)              /* wrapper lgammaf_r */
-        float x; int *signgamp;
-#endif
-{
-#ifdef _IEEE_LIBM
-	return __ieee754_lgammaf_r(x,signgamp);
-#else
-        float y;
-        y = __ieee754_lgammaf_r(x,signgamp);
-        if(_LIB_VERSION == _IEEE_) return y;
-        if(!__finitef(y)&&__finitef(x)) {
-            if(__floorf(x)==x&&x<=(float)0.0)
-	        /* gammaf pole */
-                return (float)__kernel_standard((double)x,(double)x,141);
-            else
-	        /* gamma overflow */
-                return (float)__kernel_standard((double)x,(double)x,140);
-        } else
-            return y;
-#endif
-}             
-weak_alias (__gammaf_r, gammaf_r)
diff --git a/sysdeps/libm-ieee754/w_gammal.c b/sysdeps/libm-ieee754/w_gammal.c
index f72e2e7e44..605a4167c4 100644
--- a/sysdeps/libm-ieee754/w_gammal.c
+++ b/sysdeps/libm-ieee754/w_gammal.c
@@ -19,16 +19,12 @@ static char rcsid[] = "$NetBSD: $";
 #endif
 
 /* long double gammal(double x)
- * Return the logarithm of the Gamma function of x.
- *
- * Method: call gammal_r
+ * Return the Gamma function of x.
  */
 
 #include "math.h"
 #include "math_private.h"
 
-extern int signgam;
-
 #ifdef __STDC__
 	long double __gammal(long double x)
 #else
@@ -36,12 +32,19 @@ extern int signgam;
 	long double x;
 #endif
 {
+        long double y;
+	int signgam;
+	if (_LIB_VERSION == _SVID_)
+	  y = __ieee754_lgammal_r(x,&signgam);
+	else
+	  {
+	    y = __ieee754_gammal_r(x,&signgam);
+	    if (signgam < 0) y = -y;
 #ifdef _IEEE_LIBM
-	return __ieee754_lgammal_r(x,&signgam);
+	    return y;
 #else
-        long double y;
-        y = __ieee754_lgammal_r(x,&signgam);
-        if(_LIB_VERSION == _IEEE_) return y;
+	    if(_LIB_VERSION == _IEEE_) return y;
+	  }
         if(!__finitel(y)&&__finitel(x)) {
             if(__floorl(x)==x&&x<=0.0)
                 return __kernel_standard(x,x,241); /* gamma pole */
diff --git a/sysdeps/libm-ieee754/w_gammal_r.c b/sysdeps/libm-ieee754/w_gammal_r.c
deleted file mode 100644
index ad6a43ae9b..0000000000
--- a/sysdeps/libm-ieee754/w_gammal_r.c
+++ /dev/null
@@ -1,52 +0,0 @@
-/* w_gammal_r.c -- long double version of w_gamma_r.c.
- * Conversion to long double by Ulrich Drepper,
- * Cygnus Support, drepper@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * Developed at SunPro, a Sun Microsystems, Inc. business.
- * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice
- * is preserved.
- * ====================================================
- */
-
-#if defined(LIBM_SCCS) && !defined(lint)
-static char rcsid[] = "$NetBSD: $";
-#endif
-
-/*
- * wrapper long double gammal_r(long double x, int *signgamp)
- */
-
-#include "math.h"
-#include "math_private.h"
-
-
-#ifdef __STDC__
-	long double __gammal_r(long double x, int *signgamp)
-						/* wrapper lgammal_r */
-#else
-	long double __gammal_r(x,signgamp)	/* wrapper lgamma_r */
-        long double x; int *signgamp;
-#endif
-{
-#ifdef _IEEE_LIBM
-	return __ieee754_lgammal_r(x,signgamp);
-#else
-        long double y;
-        y = __ieee754_lgammal_r(x,signgamp);
-        if(_LIB_VERSION == _IEEE_) return y;
-        if(!__finitel(y)&&__finitel(x)) {
-            if(__floorl(x)==x&&x<=0.0)
-                return __kernel_standard(x,x,241); /* gamma pole */
-            else
-                return __kernel_standard(x,x,240); /* gamma overflow */
-        } else
-            return y;
-#endif
-}
-weak_alias (__gammal_r, gammal_r)