JehanneOS/newlib
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

2009-03-25 Craig Howland <howland@LGSInnovations.com>

Browse Source 
* libc/include/math.h:  (llround, llroundf): Declare.
	* libm/common/s_llround.c: New file, implementing llround().
	* libm/common/sf_llround.c: New file, implementing llroundf().
	* libm/common/sf_lround.c: Remove spurious cast in _DOUBLE_IS_32BITS
	version of function.
	* libm/common/sf_lrint.c: Ditto.
	* libm/common/sf_logb.c:  Corrected return for subnormal argument
	by replacing existing function with a version created from sf_ilogb.c.
	* libm/common/s_logb.c: Ditto, except starting point s_ilogb.c.  Also
	added documentation for logb() and logbf().
	* libm/common/s_signbit.c:  Add signbit() documentation.
	* libm/common/s_log2.c: Update return values to match what w_log2.c has,
	since log2 uses log(); add note about being derived instead of direct.
	* libm/common/sf_fma.c: Add casts to attempt to get correct results,
	as well as comments pointing out problems with the implementation.
	* libm/common/s_fma.c: Add fma() and fmaf() documentation.
	* libm/common/sf_remquo.c: Incorrect quotient returns for large values
	corrected by discarding existing function and replacing with Sun
	verion, with some enhancements.
	* libm/common/s_remquo.c: Ditto.  Add remquo() and remquof()
	documentation.
	* libm/common/s_fmax.c: Add fmax() and fmaxf() documentation.
	* libm/common/s_fmin.c: Add fmin() and fminf() documentation.
	* libm/common/s_fdim.c: Return NAN for NAN arg, add fdim() and fdimf()
	documentation.
	* libm/common/sf_fdim.c: Return NAN for NAN arg, HUGE_VALF for inf arg.
	* libm/common/s_trunc.c: Add trunc() and truncf() documentation.
	* libm/common/s_rint.c: Add rint() and rintf() documentation.
	* libm/common/s_round.c: Add round() and roundf() documentation.
	* libm/common/s_scalbn.c: Add scalbln() and scalblnf() documentation.
	* libm/common/s_infinity.c: Add infinity() and infinityf()
	documentation.
	* libm/common/s_lround.c: Add lround(), lroundf(), llround(), and
	llroundf() documentation.
	* libm/common/s_lrint.c: Add lrint(), lrintf(), llrint(), and llrintf()
	documentation.
	* libm/common/isgreater.c: New file for documenting math.h function-like
	macros isgreater(), isgreaterequal(), isless(), islessequal(),
	islessgreater(), and isunordered().
	* libm/common/s_isnan.c: Add documentation for function-like macros
	fpclassify(), isfinite(), isinf(), isnan(), and isnormal().
	* libm/common/s_nearbyint.c: Add nearbyint() and nearbyintf()
	documentation.
	* libm/common/Makefile.am: Add s_llround.c (src); sf_llround.c (fsrc);
	s_fdim.def, s_fma.def, s_fmax.def, s_fmin.def,
	s_logb.def, s_lrint.def, s_lround.def, s_nearbyint.def, s_remquo.def,
        s_rint.def, s_round.def, s_signbit.def, s_trunc.def, and
        isgreater.def (chobj);
	re-name all existing chew files (chobj) to match source file base
	names (put in underscores), delete all special targets for chew files
	(leaving all to be generated by rule).
	* libm/common/Makefile.in: regenerate.
	* libm/math/w_exp2.c: Add "base 2" to documentation description (and
	delete TRAD_SYNOPSIS).
	* libm/math/w_gamma.c: Add tgamma() and tgammaf() documentation, along
	with some history behind the function names.
	* libm/math/math.tex: Add includes for newly-added documentation (see
	.def additions to common/Makefile.am and math/Makefile.am in this
	ChangeLog list), adjusted existing .def file names to match source file
	base names (added underscores); add mention of HUGE_VALF; rename
	"Version of library" section to "Error Handling" and add some text
	about floating-point exception; added section "Standards Compliance And
	Portability".
	* libm/math/Makefile.am: Add w_exp2.def (chobj);
	re-name all existing chew files (chobj) to match source file base
	names, delete all special targets for chew files (leaving all to be
	generated by rule).
	* libm/math/Makefile.in: regenerated
	* doc/makedoc.c: Change silent ignoring of commands < 5 characters
	to a failure when reading macro file for commands < 4 characters;
	add -v (verbose) option for printing some debugging information;
	get rid of spurious translation of "@*" to "*" (no source files used @*,
	so no existing doc pages were affected); clean up some compiler
	warnings.
	* doc/doc.str: add BUGS and SEEALSO sections (to match texi2pod.pl
	which has them); Remove ITEM command (redundant with makedoc built-in
	"o", not used in any present source file so nothing is lost, anyway).
	* HOWTO: New file to hold information for maintainers regarding how
	to do things.  Initial sections on documentation and ELIX levels.
This commit is contained in:
Jeff Johnston
2009-03-25 19:13:24 +00:00
parent baea11cd6a
commit 139f923bb4
35 changed files with 1451 additions and 652 deletions
Download Patch File Download Diff File Expand all files Collapse all files

82
newlib/ChangeLog
View File

@@ -1,3 +1,85 @@
2009-03-25 Craig Howland <howland@LGSInnovations.com>
* libc/include/math.h: (llround, llroundf): Declare.
* libm/common/s_llround.c: New file, implementing llround().
* libm/common/sf_llround.c: New file, implementing llroundf().
* libm/common/sf_lround.c: Remove spurious cast in _DOUBLE_IS_32BITS
version of function.
* libm/common/sf_lrint.c: Ditto.
* libm/common/sf_logb.c: Corrected return for subnormal argument
by replacing existing function with a version created from sf_ilogb.c.
* libm/common/s_logb.c: Ditto, except starting point s_ilogb.c. Also
added documentation for logb() and logbf().
* libm/common/s_signbit.c: Add signbit() documentation.
* libm/common/s_log2.c: Update return values to match what w_log2.c has,
since log2 uses log(); add note about being derived instead of direct.
* libm/common/sf_fma.c: Add casts to attempt to get correct results,
as well as comments pointing out problems with the implementation.
* libm/common/s_fma.c: Add fma() and fmaf() documentation.
* libm/common/sf_remquo.c: Incorrect quotient returns for large values
corrected by discarding existing function and replacing with Sun
verion, with some enhancements.
* libm/common/s_remquo.c: Ditto. Add remquo() and remquof()
documentation.
* libm/common/s_fmax.c: Add fmax() and fmaxf() documentation.
* libm/common/s_fmin.c: Add fmin() and fminf() documentation.
* libm/common/s_fdim.c: Return NAN for NAN arg, add fdim() and fdimf()
documentation.
* libm/common/sf_fdim.c: Return NAN for NAN arg, HUGE_VALF for inf arg.
* libm/common/s_trunc.c: Add trunc() and truncf() documentation.
* libm/common/s_rint.c: Add rint() and rintf() documentation.
* libm/common/s_round.c: Add round() and roundf() documentation.
* libm/common/s_scalbn.c: Add scalbln() and scalblnf() documentation.
* libm/common/s_infinity.c: Add infinity() and infinityf()
documentation.
* libm/common/s_lround.c: Add lround(), lroundf(), llround(), and
llroundf() documentation.
* libm/common/s_lrint.c: Add lrint(), lrintf(), llrint(), and llrintf()
documentation.
* libm/common/isgreater.c: New file for documenting math.h function-like
macros isgreater(), isgreaterequal(), isless(), islessequal(),
islessgreater(), and isunordered().
* libm/common/s_isnan.c: Add documentation for function-like macros
fpclassify(), isfinite(), isinf(), isnan(), and isnormal().
* libm/common/s_nearbyint.c: Add nearbyint() and nearbyintf()
documentation.
* libm/common/Makefile.am: Add s_llround.c (src); sf_llround.c (fsrc);
s_fdim.def, s_fma.def, s_fmax.def, s_fmin.def,
s_logb.def, s_lrint.def, s_lround.def, s_nearbyint.def, s_remquo.def,
s_rint.def, s_round.def, s_signbit.def, s_trunc.def, and
isgreater.def (chobj);
re-name all existing chew files (chobj) to match source file base
names (put in underscores), delete all special targets for chew files
(leaving all to be generated by rule).
* libm/common/Makefile.in: regenerate.
* libm/math/w_exp2.c: Add "base 2" to documentation description (and
delete TRAD_SYNOPSIS).
* libm/math/w_gamma.c: Add tgamma() and tgammaf() documentation, along
with some history behind the function names.
* libm/math/math.tex: Add includes for newly-added documentation (see
.def additions to common/Makefile.am and math/Makefile.am in this
ChangeLog list), adjusted existing .def file names to match source file
base names (added underscores); add mention of HUGE_VALF; rename
"Version of library" section to "Error Handling" and add some text
about floating-point exception; added section "Standards Compliance And
Portability".
* libm/math/Makefile.am: Add w_exp2.def (chobj);
re-name all existing chew files (chobj) to match source file base
names, delete all special targets for chew files (leaving all to be
generated by rule).
* libm/math/Makefile.in: regenerated
* doc/makedoc.c: Change silent ignoring of commands < 5 characters
to a failure when reading macro file for commands < 4 characters;
add -v (verbose) option for printing some debugging information;
get rid of spurious translation of "@*" to "*" (no source files used @*,
so no existing doc pages were affected); clean up some compiler
warnings.
* doc/doc.str: add BUGS and SEEALSO sections (to match texi2pod.pl
which has them); Remove ITEM command (redundant with makedoc built-in
"o", not used in any present source file so nothing is lost, anyway).
* HOWTO: New file to hold information for maintainers regarding how
to do things. Initial sections on documentation and ELIX levels.
2009-03-25 Richard Earnshaw <rearnsha@arm.com>
* libc/machine/arm/strcmp.c (strcmp_unaligned): Correctly

2
newlib/libc/include/math.h
View File

@@ -229,6 +229,7 @@ extern long int lrint _PARAMS((double));
extern _LONG_LONG_TYPE int llrint _PARAMS((double));
extern double round _PARAMS((double));
extern long int lround _PARAMS((double));
extern long long int llround _PARAMS((double));
extern double trunc _PARAMS((double));
extern double remquo _PARAMS((double, double, int *));
extern double copysign _PARAMS((double, double));
@@ -297,6 +298,7 @@ extern long int lrintf _PARAMS((float));
extern _LONG_LONG_TYPE llrintf _PARAMS((float));
extern float roundf _PARAMS((float));
extern long int lroundf _PARAMS((float));
extern long long int llroundf _PARAMS((float));
extern float truncf _PARAMS((float));
extern float remquof _PARAMS((float, float, int *));
extern float copysignf _PARAMS((float, float));

81
newlib/libm/common/Makefile.am
View File

@@ -10,7 +10,7 @@ src = s_finite.c s_copysign.c s_modf.c s_scalbn.c \
s_log1p.c s_nan.c s_nextafter.c s_pow10.c \
s_rint.c s_logb.c s_log2.c s_matherr.c s_lib_ver.c \
s_fdim.c s_fma.c s_fmax.c s_fmin.c s_fpclassify.c s_lrint.c \
s_lround.c s_nearbyint.c s_remquo.c s_round.c s_scalbln.c \
s_lround.c s_llround.c s_nearbyint.c s_remquo.c s_round.c s_scalbln.c \
s_signbit.c s_trunc.c
fsrc = sf_finite.c sf_copysign.c sf_modf.c sf_scalbn.c \
@@ -19,7 +19,7 @@ fsrc = sf_finite.c sf_copysign.c sf_modf.c sf_scalbn.c \
sf_log1p.c sf_nan.c sf_nextafter.c sf_pow10.c \
sf_rint.c sf_logb.c sf_log2.c \
sf_fdim.c sf_fma.c sf_fmax.c sf_fmin.c sf_fpclassify.c sf_lrint.c \
sf_lround.c sf_nearbyint.c sf_remquo.c sf_round.c \
sf_lround.c sf_llround.c sf_nearbyint.c sf_remquo.c sf_round.c \
sf_scalbln.c sf_trunc.c
libcommon_la_LDFLAGS = -Xcompiler -nostdlib
@@ -37,9 +37,13 @@ endif # USE_LIBTOOL
include $(srcdir)/../../Makefile.shared
chobj = scbrt.def scopysign.def sexp10.def sexpm1.def silogb.def \
sinfinity.def sisnan.def slog1p.def smatherr.def smodf.def \
snan.def snextafter.def spow10.def sscalbn.def
chobj = s_cbrt.def s_copysign.def s_exp10.def s_expm1.def s_ilogb.def \
s_infinity.def s_isnan.def s_log1p.def s_matherr.def s_modf.def \
s_nan.def s_nextafter.def s_pow10.def s_scalbn.def \
s_fdim.def s_fma.def s_fmax.def s_fmin.def \
s_logb.def s_log2.def s_lrint.def s_lround.def s_nearbyint.def \
s_remquo.def s_rint.def s_round.def s_signbit.def s_trunc.def \
isgreater.def
SUFFIXES = .def
@@ -55,73 +59,6 @@ doc: $(chobj)
CLEANFILES = $(chobj) *.ref
# Texinfo does not appear to support underscores in file names, so we
# name the .def files without underscores.
scopysign.def: s_copysign.c
$(CHEW) < $(srcdir)/s_copysign.c >$@ 2>/dev/null
touch stmp-def
scbrt.def: s_cbrt.c
$(CHEW) < $(srcdir)/s_cbrt.c >$@ 2>/dev/null
touch stmp-def
serf.def: s_erf.c
$(CHEW) < $(srcdir)/s_serf.c >$@ 2>/dev/null
touch stmp-def
sexp10.def: s_exp10.c
$(CHEW) < $(srcdir)/s_exp10.c >$@ 2>/dev/null
touch stmp-def
sexpm1.def: s_expm1.c
$(CHEW) < $(srcdir)/s_expm1.c >$@ 2>/dev/null
touch stmp-def
silogb.def: s_ilogb.c
$(CHEW) < $(srcdir)/s_ilogb.c >$@ 2>/dev/null
touch stmp-def
sinfinity.def: s_infinity.c
$(CHEW) < $(srcdir)/s_infinity.c >$@ 2>/dev/null
touch stmp-def
sisnan.def: s_isnan.c
$(CHEW) < $(srcdir)/s_isnan.c >$@ 2>/dev/null
touch stmp-def
slog1p.def: s_log1p.c
$(CHEW) < $(srcdir)/s_log1p.c >$@ 2>/dev/null
touch stmp-def
slog2.def: s_log2.c
$(CHEW) < $(srcdir)/s_log2.c >$@ 2>/dev/null
touch stmp-def
smodf.def: s_modf.c
$(CHEW) < $(srcdir)/s_modf.c >$@ 2>/dev/null
touch stmp-def
smatherr.def: s_matherr.c
$(CHEW) < $(srcdir)/s_matherr.c >$@ 2>/dev/null
touch stmp-def
snan.def: s_nan.c
$(CHEW) < $(srcdir)/s_nan.c >$@ 2>/dev/null
touch stmp-def
snextafter.def: s_nextafter.c
$(CHEW) < $(srcdir)/s_nextafter.c >$@ 2>/dev/null
touch stmp-def
spow10.def: s_pow10.c
$(CHEW) < $(srcdir)/s_pow10.c >$@ 2>/dev/null
touch stmp-def
sscalbn.def: s_scalbn.c
$(CHEW) < $(srcdir)/s_scalbn.c >$@ 2>/dev/null
touch stmp-def
# A partial dependency list.
$(lib_a_OBJECTS): $(srcdir)/../../libc/include/math.h fdlibm.h

116
newlib/libm/common/Makefile.in
View File

@@ -72,9 +72,10 @@ am__objects_1 = lib_a-s_finite.$(OBJEXT) lib_a-s_copysign.$(OBJEXT) \
lib_a-s_fma.$(OBJEXT) lib_a-s_fmax.$(OBJEXT) \
lib_a-s_fmin.$(OBJEXT) lib_a-s_fpclassify.$(OBJEXT) \
lib_a-s_lrint.$(OBJEXT) lib_a-s_lround.$(OBJEXT) \
lib_a-s_nearbyint.$(OBJEXT) lib_a-s_remquo.$(OBJEXT) \
lib_a-s_round.$(OBJEXT) lib_a-s_scalbln.$(OBJEXT) \
lib_a-s_signbit.$(OBJEXT) lib_a-s_trunc.$(OBJEXT)
lib_a-s_llround.$(OBJEXT) lib_a-s_nearbyint.$(OBJEXT) \
lib_a-s_remquo.$(OBJEXT) lib_a-s_round.$(OBJEXT) \
lib_a-s_scalbln.$(OBJEXT) lib_a-s_signbit.$(OBJEXT) \
lib_a-s_trunc.$(OBJEXT)
am__objects_2 = lib_a-sf_finite.$(OBJEXT) lib_a-sf_copysign.$(OBJEXT) \
lib_a-sf_modf.$(OBJEXT) lib_a-sf_scalbn.$(OBJEXT) \
lib_a-sf_cbrt.$(OBJEXT) lib_a-sf_exp10.$(OBJEXT) \
@@ -88,9 +89,10 @@ am__objects_2 = lib_a-sf_finite.$(OBJEXT) lib_a-sf_copysign.$(OBJEXT) \
lib_a-sf_fdim.$(OBJEXT) lib_a-sf_fma.$(OBJEXT) \
lib_a-sf_fmax.$(OBJEXT) lib_a-sf_fmin.$(OBJEXT) \
lib_a-sf_fpclassify.$(OBJEXT) lib_a-sf_lrint.$(OBJEXT) \
lib_a-sf_lround.$(OBJEXT) lib_a-sf_nearbyint.$(OBJEXT) \
lib_a-sf_remquo.$(OBJEXT) lib_a-sf_round.$(OBJEXT) \
lib_a-sf_scalbln.$(OBJEXT) lib_a-sf_trunc.$(OBJEXT)
lib_a-sf_lround.$(OBJEXT) lib_a-sf_llround.$(OBJEXT) \
lib_a-sf_nearbyint.$(OBJEXT) lib_a-sf_remquo.$(OBJEXT) \
lib_a-sf_round.$(OBJEXT) lib_a-sf_scalbln.$(OBJEXT) \
lib_a-sf_trunc.$(OBJEXT)
@USE_LIBTOOL_FALSE@am_lib_a_OBJECTS = $(am__objects_1) \
@USE_LIBTOOL_FALSE@ $(am__objects_2)
lib_a_OBJECTS = $(am_lib_a_OBJECTS)
@@ -102,15 +104,16 @@ am__objects_3 = s_finite.lo s_copysign.lo s_modf.lo s_scalbn.lo \
s_log1p.lo s_nan.lo s_nextafter.lo s_pow10.lo s_rint.lo \
s_logb.lo s_log2.lo s_matherr.lo s_lib_ver.lo s_fdim.lo \
s_fma.lo s_fmax.lo s_fmin.lo s_fpclassify.lo s_lrint.lo \
s_lround.lo s_nearbyint.lo s_remquo.lo s_round.lo s_scalbln.lo \
s_signbit.lo s_trunc.lo
s_lround.lo s_llround.lo s_nearbyint.lo s_remquo.lo s_round.lo \
s_scalbln.lo s_signbit.lo s_trunc.lo
am__objects_4 = sf_finite.lo sf_copysign.lo sf_modf.lo sf_scalbn.lo \
sf_cbrt.lo sf_exp10.lo sf_expm1.lo sf_ilogb.lo sf_infinity.lo \
sf_isinf.lo sf_isinff.lo sf_isnan.lo sf_isnanf.lo sf_log1p.lo \
sf_nan.lo sf_nextafter.lo sf_pow10.lo sf_rint.lo sf_logb.lo \
sf_log2.lo sf_fdim.lo sf_fma.lo sf_fmax.lo sf_fmin.lo \
sf_fpclassify.lo sf_lrint.lo sf_lround.lo sf_nearbyint.lo \
sf_remquo.lo sf_round.lo sf_scalbln.lo sf_trunc.lo
sf_fpclassify.lo sf_lrint.lo sf_lround.lo sf_llround.lo \
sf_nearbyint.lo sf_remquo.lo sf_round.lo sf_scalbln.lo \
sf_trunc.lo
@USE_LIBTOOL_TRUE@am_libcommon_la_OBJECTS = $(am__objects_3) \
@USE_LIBTOOL_TRUE@ $(am__objects_4)
libcommon_la_OBJECTS = $(am_libcommon_la_OBJECTS)
@@ -280,7 +283,7 @@ src = s_finite.c s_copysign.c s_modf.c s_scalbn.c \
s_log1p.c s_nan.c s_nextafter.c s_pow10.c \
s_rint.c s_logb.c s_log2.c s_matherr.c s_lib_ver.c \
s_fdim.c s_fma.c s_fmax.c s_fmin.c s_fpclassify.c s_lrint.c \
s_lround.c s_nearbyint.c s_remquo.c s_round.c s_scalbln.c \
s_lround.c s_llround.c s_nearbyint.c s_remquo.c s_round.c s_scalbln.c \
s_signbit.c s_trunc.c
fsrc = sf_finite.c sf_copysign.c sf_modf.c sf_scalbn.c \
@@ -289,7 +292,7 @@ fsrc = sf_finite.c sf_copysign.c sf_modf.c sf_scalbn.c \
sf_log1p.c sf_nan.c sf_nextafter.c sf_pow10.c \
sf_rint.c sf_logb.c sf_log2.c \
sf_fdim.c sf_fma.c sf_fmax.c sf_fmin.c sf_fpclassify.c sf_lrint.c \
sf_lround.c sf_nearbyint.c sf_remquo.c sf_round.c \
sf_lround.c sf_llround.c sf_nearbyint.c sf_remquo.c sf_round.c \
sf_scalbln.c sf_trunc.c
libcommon_la_LDFLAGS = -Xcompiler -nostdlib
@@ -300,9 +303,13 @@ libcommon_la_LDFLAGS = -Xcompiler -nostdlib
@USE_LIBTOOL_FALSE@noinst_LIBRARIES = lib.a
@USE_LIBTOOL_FALSE@lib_a_SOURCES = $(src) $(fsrc)
@USE_LIBTOOL_FALSE@lib_a_CFLAGS = $(AM_CFLAGS)
chobj = scbrt.def scopysign.def sexp10.def sexpm1.def silogb.def \
sinfinity.def sisnan.def slog1p.def smatherr.def smodf.def \
snan.def snextafter.def spow10.def sscalbn.def
chobj = s_cbrt.def s_copysign.def s_exp10.def s_expm1.def s_ilogb.def \
s_infinity.def s_isnan.def s_log1p.def s_matherr.def s_modf.def \
s_nan.def s_nextafter.def s_pow10.def s_scalbn.def \
s_fdim.def s_fma.def s_fmax.def s_fmin.def \
s_logb.def s_log2.def s_lrint.def s_lround.def s_nearbyint.def \
s_remquo.def s_rint.def s_round.def s_signbit.def s_trunc.def \
isgreater.def
SUFFIXES = .def
CHEW = ../../doc/makedoc -f $(srcdir)/../../doc/doc.str
@@ -555,6 +562,12 @@ lib_a-s_lround.o: s_lround.c
lib_a-s_lround.obj: s_lround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-s_lround.obj `if test -f 's_lround.c'; then $(CYGPATH_W) 's_lround.c'; else $(CYGPATH_W) '$(srcdir)/s_lround.c'; fi`
lib_a-s_llround.o: s_llround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-s_llround.o `test -f 's_llround.c' || echo '$(srcdir)/'`s_llround.c
lib_a-s_llround.obj: s_llround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-s_llround.obj `if test -f 's_llround.c'; then $(CYGPATH_W) 's_llround.c'; else $(CYGPATH_W) '$(srcdir)/s_llround.c'; fi`
lib_a-s_nearbyint.o: s_nearbyint.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-s_nearbyint.o `test -f 's_nearbyint.c' || echo '$(srcdir)/'`s_nearbyint.c
@@ -753,6 +766,12 @@ lib_a-sf_lround.o: sf_lround.c
lib_a-sf_lround.obj: sf_lround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-sf_lround.obj `if test -f 'sf_lround.c'; then $(CYGPATH_W) 'sf_lround.c'; else $(CYGPATH_W) '$(srcdir)/sf_lround.c'; fi`
lib_a-sf_llround.o: sf_llround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-sf_llround.o `test -f 'sf_llround.c' || echo '$(srcdir)/'`sf_llround.c
lib_a-sf_llround.obj: sf_llround.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-sf_llround.obj `if test -f 'sf_llround.c'; then $(CYGPATH_W) 'sf_llround.c'; else $(CYGPATH_W) '$(srcdir)/sf_llround.c'; fi`
lib_a-sf_nearbyint.o: sf_nearbyint.c
$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(lib_a_CFLAGS) $(CFLAGS) -c -o lib_a-sf_nearbyint.o `test -f 'sf_nearbyint.c' || echo '$(srcdir)/'`sf_nearbyint.c
@@ -943,73 +962,6 @@ objectlist.awk.in: $(noinst_LTLIBRARIES)
doc: $(chobj)
# Texinfo does not appear to support underscores in file names, so we
# name the .def files without underscores.
scopysign.def: s_copysign.c
$(CHEW) < $(srcdir)/s_copysign.c >$@ 2>/dev/null
touch stmp-def
scbrt.def: s_cbrt.c
$(CHEW) < $(srcdir)/s_cbrt.c >$@ 2>/dev/null
touch stmp-def
serf.def: s_erf.c
$(CHEW) < $(srcdir)/s_serf.c >$@ 2>/dev/null
touch stmp-def
sexp10.def: s_exp10.c
$(CHEW) < $(srcdir)/s_exp10.c >$@ 2>/dev/null
touch stmp-def
sexpm1.def: s_expm1.c
$(CHEW) < $(srcdir)/s_expm1.c >$@ 2>/dev/null
touch stmp-def
silogb.def: s_ilogb.c
$(CHEW) < $(srcdir)/s_ilogb.c >$@ 2>/dev/null
touch stmp-def
sinfinity.def: s_infinity.c
$(CHEW) < $(srcdir)/s_infinity.c >$@ 2>/dev/null
touch stmp-def
sisnan.def: s_isnan.c
$(CHEW) < $(srcdir)/s_isnan.c >$@ 2>/dev/null
touch stmp-def
slog1p.def: s_log1p.c
$(CHEW) < $(srcdir)/s_log1p.c >$@ 2>/dev/null
touch stmp-def
slog2.def: s_log2.c
$(CHEW) < $(srcdir)/s_log2.c >$@ 2>/dev/null
touch stmp-def
smodf.def: s_modf.c
$(CHEW) < $(srcdir)/s_modf.c >$@ 2>/dev/null
touch stmp-def
smatherr.def: s_matherr.c
$(CHEW) < $(srcdir)/s_matherr.c >$@ 2>/dev/null
touch stmp-def
snan.def: s_nan.c
$(CHEW) < $(srcdir)/s_nan.c >$@ 2>/dev/null
touch stmp-def
snextafter.def: s_nextafter.c
$(CHEW) < $(srcdir)/s_nextafter.c >$@ 2>/dev/null
touch stmp-def
spow10.def: s_pow10.c
$(CHEW) < $(srcdir)/s_pow10.c >$@ 2>/dev/null
touch stmp-def
sscalbn.def: s_scalbn.c
$(CHEW) < $(srcdir)/s_scalbn.c >$@ 2>/dev/null
touch stmp-def
# A partial dependency list.
$(lib_a_OBJECTS): $(srcdir)/../../libc/include/math.h fdlibm.h

75
newlib/libm/common/isgreater.c Normal file
View File

@@ -0,0 +1,75 @@
/* isgreater.c: This file contains no source code, but rather only the
* man-page comments. All of the documented "functions" are actually macros
* defined in math.h (q.v.). */
/*
FUNCTION
<<isgreater>>, <<isgreaterequal>>, <<isless>>, <<islessequal>>, <<islessgreater>>, and <<isunordered>>--comparison macros
INDEX
isgreater
INDEX
isgreaterequal
INDEX
isless
INDEX
islessequal
INDEX
islessgreater
INDEX
isunordered
ANSI_SYNOPSIS
#include <math.h>
int isgreater(real-floating <[x]>, real-floating <[y]>);
int isgreaterequal(real-floating <[x]>, real-floating <[y]>);
int isless(real-floating <[x]>, real-floating <[y]>);
int islessequal(real-floating <[x]>, real-floating <[y]>);
int islessgreater(real-floating <[x]>, real-floating <[y]>);
int isunordered(real-floating <[x]>, real-floating <[y]>);
DESCRIPTION
<<isgreater>>, <<isgreaterequal>>, <<isless>>, <<islessequal>>,
<<islessgreater>>, and <<isunordered>> are macros defined for use in
comparing floating-point numbers without raising any floating-point
exceptions.
The relational operators (i.e. <, >, <=, and >=) support the usual mathematical
relationships between numeric values. For any ordered pair of numeric
values exactly one of the relationships--less, greater, and equal--is
true. Relational operators may raise the "invalid" floating-point
exception when argument values are NaNs. For a NaN and a numeric value, or
for two NaNs, just the unordered relationship is true (i.e., if one or both
of the arguments a NaN, the relationship is called unordered). The specified
macros are quiet (non floating-point exception raising) versions of the
relational operators, and other comparison macros that facilitate writing
efficient code that accounts for NaNs without suffering the "invalid"
floating-point exception. In the synopses shown, "real-floating" indicates
that the argument is an expression of real floating type.
Please note that saying that the macros do not raise floating-point
exceptions, it is referring to the function that they are performing. It
is certainly possible to give them an expression which causes an exception.
For example:
o+
o NaN < 1.0
causes an "invalid" exception,
o isless(NaN, 1.0)
does not, and
o isless(NaN*0., 1.0)
causes an exception due to the "NaN*0.", but not from the
resultant reduced comparison of isless(NaN, 1.0).
o-
RETURNS
@comment Formatting note: "$@" forces a new line
No floating-point exceptions are raised for any of the macros.@*
The <<isgreater>> macro returns the value of (x) > (y).@*
The <<isgreaterequal>> macro returns the value of (x) >= (y).@*
The <<isless>> macro returns the value of (x) < (y).@*
The <<islessequal>> macro returns the value of (x) <= (y).@*
The <<islessgreater>> macro returns the value of (x) < (y) || (x) > (y).@*
The <<isunordered>> macro returns 1 if either of its arguments is NaN and 0 otherwise.
PORTABILITY
C99, POSIX.
*/

37
newlib/libm/common/s_fdim.c
View File

@@ -3,6 +3,39 @@
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
*/
/*
FUNCTION
<<fdim>>, <<fdimf>>--positive difference
INDEX
fdim
INDEX
fdimf
ANSI_SYNOPSIS
#include <math.h>
double fdim(double <[x]>, double <[y]>);
float fdimf(float <[x]>, float <[y]>);
DESCRIPTION
The <<fdim>> functions determine the positive difference between their
arguments, returning:
. <[x]> - <[y]> if <[x]> > <[y]>, or
@ifnottex
. +0 if <[x]> <= <[y]>, or
@end ifnottex
@tex
. +0 if <[x]> $\leq$ <[y]>, or
@end tex
. NAN if either argument is NAN.
A range error may occur.
RETURNS
The <<fdim>> functions return the positive difference value.
PORTABILITY
ANSI C, POSIX.
*/
#include "fdlibm.h"
@@ -17,7 +50,9 @@
#endif
{
int c = __fpclassifyd(x);
if (c == FP_NAN || c == FP_INFINITE)
if (c == FP_NAN) return(x);
if (__fpclassifyd(y) == FP_NAN) return(y);
if (c == FP_INFINITE)
return HUGE_VAL;
return x > y ? x - y : 0.0;

38
newlib/libm/common/s_fma.c
View File

@@ -1,3 +1,41 @@
/*
FUNCTION
<<fma>>, <<fmaf>>--floating multiply add
INDEX
fma
INDEX
fmaf
ANSI_SYNOPSIS
#include <math.h>
double fma(double <[x]>, double <[y]>, double <[z]>);
float fmaf(float <[x]>, float <[y]>, float <[z]>);
DESCRIPTION
The <<fma>> functions compute (<[x]> * <[y]>) + <[z]>, rounded as one ternary
operation: they compute the value (as if) to infinite precision and round once
to the result format, according to the rounding mode characterized by the value
of FLT_ROUNDS. That is, they are supposed to do this: see below.
RETURNS
The <<fma>> functions return (<[x]> * <[y]>) + <[z]>, rounded as one ternary
operation.
BUGS
This implementation does not provide the function that it should, purely
returning "(<[x]> * <[y]>) + <[z]>;" with no attempt at all to provide the
simulated infinite precision intermediates which are required. DO NOT USE THEM.
If double has enough more precision than float, then <<fmaf>> should provide
the expected numeric results, as it does use double for the calculation. But
since this is not the case for all platforms, this manual cannot determine
if it is so for your case.
PORTABILITY
ANSI C, POSIX.
*/
#include "fdlibm.h"
#ifndef _DOUBLE_IS_32BITS

25
newlib/libm/common/s_fmax.c
View File

@@ -3,6 +3,31 @@
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
*/
/*
FUNCTION
<<fmax>>, <<fmaxf>>--maximum
INDEX
fmax
INDEX
fmaxf
ANSI_SYNOPSIS
#include <math.h>
double fmax(double <[x]>, double <[y]>);
float fmaxf(float <[x]>, float <[y]>);
DESCRIPTION
The <<fmax>> functions determine the maximum numeric value of their arguments.
NaN arguments are treated as missing data: if one argument is a NaN and the
other numeric, then the <<fmax>> functions choose the numeric value.
RETURNS
The <<fmax>> functions return the maximum numeric value of their arguments.
PORTABILITY
ANSI C, POSIX.
*/
#include "fdlibm.h"

25
newlib/libm/common/s_fmin.c
View File

@@ -3,6 +3,31 @@
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
*/
/*
FUNCTION
<<fmin>>, <<fminf>>--minimum
INDEX
fmin
INDEX
fminf
ANSI_SYNOPSIS
#include <math.h>
double fmin(double <[x]>, double <[y]>);
float fminf(float <[x]>, float <[y]>);
DESCRIPTION
The <<fmin>> functions determine the minimum numeric value of their arguments.
NaN arguments are treated as missing data: if one argument is a NaN and the
other numeric, then the <<fmin>> functions choose the numeric value.
RETURNS
The <<fmin>> functions return the minimum numeric value of their arguments.
PORTABILITY
ANSI C, POSIX.
*/
#include "fdlibm.h"

15
newlib/libm/common/s_infinity.c
View File

@@ -5,7 +5,7 @@
/*
FUNCTION
<<infinity>>, <<infinityf>>---representation of infinity
<<infinity>>, <<infinityf>>--representation of infinity
INDEX
infinity
@@ -17,17 +17,18 @@ ANSI_SYNOPSIS
double infinity(void);
float infinityf(void);
TRAD_SYNOPSIS
#include <math.h>
double infinity();
float infinityf();
DESCRIPTION
<<infinity>> and <<infinityf>> return the special number IEEE
infinity in double- and single-precision arithmetic
respectively.
PORTABILITY
<<infinity>> and <<infinityf>> are neither standard C nor POSIX. C and
POSIX require macros HUGE_VAL and HUGE_VALF to be defined in math.h, which
Newlib defines to be infinities corresponding to these archaic infinity()
and infinityf() functions in floating-point implementations which do have
infinities.
QUICKREF
infinity - pure

113
newlib/libm/common/s_isnan.c
View File

@@ -13,8 +13,21 @@
/*
FUNCTION
<<isnan>>, <<isnanf>>, <<isinf>>, <<isinff>>, <<finite>>, <<finitef>>---test for exceptional numbers
<<fpclassify>>, <<isfinite>>, <<isinf>>, <<isnan>>, and <<isnormal>>--floating-point classification macros; <<finite>>, <<finitef>>, <<isinf>>, <<isinff>>, <<isnan>>, <<isnanf>>--test for exceptional numbers
@c C99 (start
INDEX
fpclassify
INDEX
isfinite
INDEX
isinf
INDEX
isnan
INDEX
isnormal
@c C99 end)
@c SUSv2 (start
INDEX
isnan
INDEX
@@ -28,8 +41,18 @@ INDEX
isinff
INDEX
finitef
@c SUSv2 end)
ANSI_SYNOPSIS
[C99 standard macros:]
#include <math.h>
int fpclassify(real-floating <[x]>);
int isfinite(real-floating <[x]>);
int isinf(real-floating <[x]>);
int isnan(real-floating <[x]>);
int isnormal(real-floating <[x]>);
[Archaic SUSv2 functions:]
#include <ieeefp.h>
int isnan(double <[arg]>);
int isinf(double <[arg]>);
@@ -38,30 +61,64 @@ ANSI_SYNOPSIS
int isinff(float <[arg]>);
int finitef(float <[arg]>);
TRAD_SYNOPSIS
#include <ieeefp.h>
int isnan(<[arg]>)
double <[arg]>;
int isinf(<[arg]>)
double <[arg]>;
int finite(<[arg]>);
double <[arg]>;
int isnanf(<[arg]>);
float <[arg]>;
int isinff(<[arg]>);
float <[arg]>;
int finitef(<[arg]>);
float <[arg]>;
DESCRIPTION
These functions provide information on the floating-point
<<fpclassify>>, <<isfinite>>, <<isinf>>, <<isnan>>, and <<isnormal>> are macros
defined for use in classifying floating-point numbers. This is a help because
of special "values" like NaN and infinities. In the synopses shown,
"real-floating" indicates that the argument is an expression of real floating
type. These function-like macros are C99 and POSIX-compliant, and should be
used instead of the now-archaic SUSv2 functions.
The <<fpclassify>> macro classifies its argument value as NaN, infinite, normal,
subnormal, zero, or into another implementation-defined category. First, an
argument represented in a format wider than its semantic type is converted to
its semantic type. Then classification is based on the type of the argument.
The <<fpclassify>> macro returns the value of the number classification macro
appropriate to the value of its argument:
o+
o FP_INFINITE
<[x]> is either plus or minus infinity;
o FP_NAN
<[x]> is "Not A Number" (plus or minus);
o FP_NORMAL
<[x]> is a "normal" number (i.e. is none of the other special forms);
o FP_SUBNORMAL
<[x]> is too small be stored as a regular normalized number (i.e. loss of precision is likely); or
o FP_ZERO
<[x]> is 0 (either plus or minus).
o-
The "<<is>>" set of macros provide a useful set of shorthand ways for
classifying floating-point numbers, providing the following equivalent
relations:
o+
o <<isfinite>>(<[x]>)
returns non-zero if <[x]> is finite. (It is equivalent to
(<<fpclassify>>(<[x]>) != FP_INFINITE && <<fpclassify>>(<[x]>) != FP_NAN).)
o <<isinf>>(<[x]>)
returns non-zero if <[x]> is infinite. (It is equivalent to
(<<fpclassify>>(<[x]>) == FP_INFINITE).)
o <<isnan>>(<[x]>)
returns non-zero if <[x]> is NaN. (It is equivalent to
(<<fpclassify>>(<[x]>) == FP_NAN).)
o <<isnormal>>(<[x]>)
returns non-zero if <[x]> is normal. (It is equivalent to
(<<fpclassify>>(<[x]>) == FP_NORMAL).)
o-
The archaic SUSv2 functions provide information on the floating-point
argument supplied.
There are five major number formats:
There are five major number formats ("exponent" referring to the
biased exponent in the binary-encoded number):
o+
o zero
A number which contains all zero bits.
A number which contains all zero bits, excluding the sign bit.
o subnormal
A number with a zero exponent but a nonzero fraction.
o normal
@@ -85,7 +142,21 @@ DESCRIPTION
and <<isinf>> are macros that operate on multiple types of
floating-point. The SUSv2 standard declares <<isnan>> as
a function taking double. Newlib has decided to declare
them both as macros in math.h and as functions in ieeefp.h.
them both as macros in math.h and as functions in ieeefp.h to
maintain backward compatibility.
RETURNS
@comment Formatting note: "$@" forces a new line
The fpclassify macro returns the value corresponding to the appropriate FP_ macro.@*
The isfinite macro returns nonzero if <[x]> is finite, else 0.@*
The isinf macro returns nonzero if <[x]> is infinite, else 0.@*
The isnan macro returns nonzero if <[x]> is an NaN, else 0.@*
The isnormal macro returns nonzero if <[x]> has a normal value, else 0.
PORTABILITY
math.h macros are C99, POSIX.
ieeefp.h funtions are outdated and should be avoided.
QUICKREF
isnan - pure

68
newlib/libm/common/s_llround.c Normal file
View File

@@ -0,0 +1,68 @@
/* lround adapted to be llround for Newlib, 2009 by Craig Howland. */
/*
* ====================================================
* 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.
* ====================================================
*/
#include "fdlibm.h"
#ifndef _DOUBLE_IS_32BITS
long long int
llround(double x)
{
__int32_t sign, exponent_less_1023;
/* Most significant word, least significant word. */
__uint32_t msw, lsw;
long long int result;
EXTRACT_WORDS(msw, lsw, x);
/* Extract sign. */
sign = ((msw & 0x80000000) ? -1 : 1);
/* Extract exponent field. */
exponent_less_1023 = ((msw & 0x7ff00000) >> 20) - 1023;
msw &= 0x000fffff;
msw |= 0x00100000;
if (exponent_less_1023 < 20)
{
if (exponent_less_1023 < 0)
{
if (exponent_less_1023 < -1)
return 0;
else
return sign;
}
else
{
msw += 0x80000 >> exponent_less_1023;
result = msw >> (20 - exponent_less_1023);
}
}
else if (exponent_less_1023 < (8 * sizeof (long long int)) - 1)
{
if (exponent_less_1023 >= 52)
result = ((long long int) msw << (exponent_less_1023 - 20)) | (lsw << (exponent_less_1023 - 52));
else
{
unsigned int tmp = lsw + (0x80000000 >> (exponent_less_1023 - 20));
if (tmp < lsw)
++msw;
result = ((long long int) msw << (exponent_less_1023 - 20)) | (tmp >> (52 - exponent_less_1023));
}
}
else
/* Result is too large to be represented by a long long int. */
return (long long int)x;
return sign * result;
}
#endif /* _DOUBLE_IS_32BITS */

46
newlib/libm/common/s_log2.c
View File

@@ -14,7 +14,7 @@
/*
FUNCTION
<<log2>>, <<log2f>>---base 2 logarithm
<<log2>>, <<log2f>>--base 2 logarithm
INDEX
log2
INDEX
@@ -25,28 +25,40 @@ ANSI_SYNOPSIS
double log2(double <[x]>);
float log2f(float <[x]>);
TRAD_SYNOPSIS
#include <math.h>
double log2(<[x]>);
double <[x]>;
float log2f(<[x]>);
float <[x]>;
DESCRIPTION
<<log2>> returns the base 2 logarithm of <[x]>.
The <<log2>> functions compute the base-2 logarithm of <[x]>. A domain error
occurs if the argument is less than zero. A range error occurs if the
argument is zero.
<<log2f>> is identical, save that it takes and returns <<float>> values.
The Newlib implementations are not full, intrinisic calculations, but
rather are derivatives based on <<log>>. (Accuracy might be slightly off from
a direct calculation.) In addition to functions, they are also implemented as
macros defined in math.h:
. #define log2(x) (log (x) / _M_LOG2_E)
. #define log2f(x) (logf (x) / (float) _M_LOG2_E)
To use the functions instead, just undefine the macros first.
You can use the (non-ANSI) function <<matherr>> to specify error
handling for these functions, indirectly through the respective <<log>>
function.
RETURNS
On success, <<log2>> and <<log2f>> return the calculated value.
If the result underflows, the returned value is <<0>>. If the
result overflows, the returned value is <<HUGE_VAL>>. In
either case, <<errno>> is set to <<ERANGE>>.
The <<log2>> functions return
@ifnottex
<<log base-2(<[x]>)>>
@end ifnottex
@tex
$log_2(x)$
@end tex
on success.
When <[x]> is zero, the
returned value is <<-HUGE_VAL>> and <<errno>> is set to <<ERANGE>>.
When <[x]> is negative, the returned value is NaN (not a number) and
<<errno>> is set to <<EDOM>>. You can control the error behavior via
<<matherr>>.
PORTABILITY
<<log2>> and <<log2f>> are required by ISO/IEC 9899:1999 and the
System V Interface Definition (Issue 6).
C99, POSIX, System V Interface Definition (Issue 6).
*/
/*

104
newlib/libm/common/s_logb.c
View File

@@ -1,5 +1,5 @@
/* @(#)s_logb.c 5.1 93/09/24 */
/* 2009 for Newlib: Sun's s_ilogb.c converted to be s_logb.c. */
/* @(#)s_ilogb.c 5.1 93/09/24 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
@@ -10,33 +10,101 @@
* is preserved.
* ====================================================
*/
/*
* double logb(x)
* IEEE 754 logb. Included to pass IEEE test suite. Not recommend.
* Use ilogb instead.
FUNCTION
<<logb>>, <<logbf>>--get exponent of floating-point number
INDEX
logb
INDEX
logbf
ANSI_SYNOPSIS
#include <math.h>
double logb(double <[x]>);
float logbf(float <[x]>);
DESCRIPTION
The <<logb>> functions extract the exponent of <[x]>, as a signed integer value
in floating-point format. If <[x]> is subnormal it is treated as though it were
normalized; thus, for positive finite <[x]>,
@ifnottex
1 <= (<[x]> * FLT_RADIX to the power (-logb(<[x]>))) < FLT_RADIX.
@end ifnottex
@tex
$1 \leq ( x \cdot FLT\_RADIX ^ {-logb(x)} ) < FLT\_RADIX$.
@end tex
A domain error may occur if the argument is zero.
In this floating-point implementation, FLT_RADIX is 2. Which also means
that for finite <[x]>, <<logb>>(<[x]>) = <<floor>>(<<log2>>(<<fabs>>(<[x]>))).
All nonzero, normal numbers can be described as
@ifnottex
<[m]> * 2**<[p]>, where 1.0 <= <[m]> < 2.0.
@end ifnottex
@tex
$m \cdot 2^p$, where $1.0 \leq m < 2.0$.
@end tex
The <<logb>> functions examine the argument <[x]>, and return <[p]>.
The <<frexp>> functions are similar to the <<logb>> functions, but
returning <[m]> adjusted to the interval [.5, 1) or 0, and <[p]>+1.
RETURNS
@comment Formatting note: "$@" forces a new line
When <[x]> is:@*
+inf or -inf, +inf is returned;@*
NaN, NaN is returned;@*
0, -inf is returned, and the divide-by-zero exception is raised;@*
otherwise, the <<logb>> functions return the signed exponent of <[x]>.
PORTABILITY
ANSI C, POSIX
SEEALSO
frexp, ilogb
*/
/* double logb(double x)
* return the binary exponent of non-zero x
* logb(0) = -inf, raise divide-by-zero floating point exception
* logb(+inf|-inf) = +inf (no signal is raised)
* logb(NaN) = NaN (no signal is raised)
* Per C99 recommendation, a NaN argument is returned unchanged.
*/
#include "fdlibm.h"
#ifndef _DOUBLE_IS_32BITS
double
#ifdef __STDC__
double logb(double x)
logb(double x)
#else
double logb(x)
double x;
logb(x)
double x;
#endif
{
__int32_t lx,ix;
EXTRACT_WORDS(ix,lx,x);
ix &= 0x7fffffff; /* high |x| */
if((ix|lx)==0) return -1.0/fabs(x);
if(ix>=0x7ff00000) return x*x;
if((ix>>=20)==0) /* IEEE 754 logb */
return -1022.0;
else
return (double) (ix-1023);
__int32_t hx,lx,ix;
EXTRACT_WORDS(hx,lx,x);
hx &= 0x7fffffff; /* high |x| */
if(hx<0x00100000) { /* 0 or subnormal */
if((hx|lx)==0) {
double xx;
/* arg==0: return -inf and raise divide-by-zero exception */
INSERT_WORDS(xx,hx,lx); /* +0.0 */
return -1./xx; /* logb(0) = -inf */
}
else /* subnormal x */
if(hx==0) {
for (ix = -1043; lx>0; lx<<=1) ix -=1;
} else {
for (ix = -1022,hx<<=11; hx>0; hx<<=1) ix -=1;
}
return (double) ix;
}
else if (hx<0x7ff00000) return (hx>>20)-1023; /* normal # */
else if (hx>0x7ff00000 || lx) return x; /* x==NaN */
else return HUGE_VAL; /* x==inf (+ or -) */
}
#endif /* _DOUBLE_IS_32BITS */

38
newlib/libm/common/s_lrint.c
View File

@@ -10,6 +10,44 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<lrint>>, <<lrintf>>, <<llrint>>, <<llrintf>>--round to integer
INDEX
lrint
INDEX
lrintf
INDEX
llrint
INDEX
llrintf
ANSI_SYNOPSIS
#include <math.h>
long int lrint(double <[x]>);
long int lrintf(float <[x]>);
long long int llrint(double <[x]>);
long long int llrintf(float <[x]>);
DESCRIPTION
The <<lrint>> and <<llrint>> functions round their argument to the nearest
integer value, using the current rounding direction. If the rounded value is
outside the range of the return type, the numeric result is unspecified. A
range error may occur if the magnitude of <[x]> is too large.
The "inexact" floating-point exception is raised in implementations that
support it when the result differs in value from the argument (i.e., when
a fraction actually has been truncated).
RETURNS
<[x]> rounded to an integral value, using the current rounding direction.
SEEALSO
<<lround>>
PORTABILITY
ANSI C, POSIX
*/
/*
* lrint(x)

38
newlib/libm/common/s_lround.c
View File

@@ -8,6 +8,44 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<lround>>, <<lroundf>>, <<llround>>, <<llroundf>>--round to integer, to nearest
INDEX
lround
INDEX
lroundf
INDEX
llround
INDEX
llroundf
ANSI_SYNOPSIS
#include <math.h>
long int lround(double <[x]>);
long int lroundf(float <[x]>);
long long int llround(double <[x]>);
long long int llroundf(float <[x]>);
DESCRIPTION
The <<lround>> and <<llround>> functions round their argument to the
nearest integer value, rounding halfway cases away from zero, regardless
of the current rounding direction. If the rounded value is outside the
range of the return type, the numeric result is unspecified (depending
upon the floating-point implementation, not the library). A range
error may occur if the magnitude of x is too large.
RETURNS
<[x]> rounded to an integral value as an integer.
SEEALSO
See the <<round>> functions for the return being the same floating-point type
as the argument. <<lrint>>, <<llrint>>.
PORTABILITY
ANSI C, POSIX
*/
#include "fdlibm.h"

37
newlib/libm/common/s_nearbyint.c
View File

@@ -8,6 +8,43 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<nearbyint>>, <<nearbyintf>>--round to integer
INDEX
nearbyint
INDEX
nearbyintf
ANSI_SYNOPSIS
#include <math.h>
double nearbyint(double <[x]>);
float nearbyintf(float <[x]>);
DESCRIPTION
The <<nearbyint>> functions round their argument to an integer value in
floating-point format, using the current rounding direction and
(supposedly) without raising the "inexact" floating-point exception.
See the <<rint>> functions for the same function with the "inexact"
floating-point exception being raised when appropriate.
BUGS
Newlib does not support the floating-point exception model, so that
the floating-point exception control is not present and thereby what may
be seen will be compiler and hardware dependent in this regard.
The Newlib <<nearbyint>> functions are identical to the <<rint>>
functions with respect to the floating-point exception behavior, and
will cause the "inexact" exception to be raised for most targets.
RETURNS
<[x]> rounded to an integral value, using the current rounding direction.
PORTABILITY
ANSI C, POSIX
SEEALSO
<<rint>>, <<round>>
*/
#include <math.h>
#include "fdlibm.h"

219
newlib/libm/common/s_remquo.c
View File

@@ -1,39 +1,208 @@
/* Copyright (C) 2002 by Red Hat, Incorporated. All rights reserved.
/* Adapted for Newlib, 2009. (Allow for int < 32 bits; return *quo=0 during
* errors to make test scripts easier.) */
/* @(#)e_fmod.c 1.3 95/01/18 */
/*-
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<remquo>>, <<remquof>>--remainder and part of quotient
INDEX
remquo
INDEX
remquof
ANSI_SYNOPSIS
#include <math.h>
double remquo(double <[x]>, double <[y]>, int *<[quo]>);
float remquof(float <[x]>, float <[y]>, int *<[quo]>);
DESCRIPTION
The <<remquo>> functions compute the same remainder as the <<remainder>>
functions; this value is in the range -<[y]>/2 ... +<[y]>/2. In the object
pointed to by <<quo>> they store a value whose sign is the sign of <<x>>/<<y>>
and whose magnitude is congruent modulo 2**n to the magnitude of the integral
quotient of <<x>>/<<y>>. (That is, <<quo>> is given the n lsbs of the
quotient, not counting the sign.) This implementation uses n=31 if int is 32
bits or more, otherwise, n is 1 less than the width of int.
For example:
. remquo(-29.0, 3.0, &<[quo]>)
returns -1.0 and sets <[quo]>=10, and
. remquo(-98307.0, 3.0, &<[quo]>)
returns -0.0 and sets <[quo]>=-32769, although for 16-bit int, <[quo]>=-1. In
the latter case, the actual quotient of -(32769=0x8001) is reduced to -1
because of the 15-bit limitation for the quotient.
RETURNS
When either argument is NaN, NaN is returned. If <[y]> is 0 or <[x]> is
infinite (and neither is NaN), a domain error occurs (i.e. the "invalid"
floating point exception is raised or errno is set to EDOM), and NaN is
returned.
Otherwise, the <<remquo>> functions return <[x]> REM <[y]>.
BUGS
IEEE754-2008 calls for <<remquo>>(subnormal, inf) to cause the "underflow"
floating-point exception. This implementation does not.
PORTABILITY
C99, POSIX.
*/
#include <limits.h>
#include <math.h>
#include "fdlibm.h"
#ifndef _DOUBLE_IS_32BITS
#ifdef __STDC__
double remquo(double x, double y, int *quo) /* wrapper remquo */
#else
double remquo(x,y,quo) /* wrapper remquo */
double x,y;
int *quo;
/* For quotient, return either all 31 bits that can from calculation (using
* int32_t), or as many as can fit into an int that is smaller than 32 bits. */
#if INT_MAX > 0x7FFFFFFFL
#define QUO_MASK 0x7FFFFFFF
# else
#define QUO_MASK INT_MAX
#endif
static const double Zero[] = {0.0, -0.0,};
/*
* Return the IEEE remainder and set *quo to the last n bits of the
* quotient, rounded to the nearest integer. We choose n=31--if that many fit--
* because we wind up computing all the integer bits of the quotient anyway as
* a side-effect of computing the remainder by the shift and subtract
* method. In practice, this is far more bits than are needed to use
* remquo in reduction algorithms.
*/
double
remquo(double x, double y, int *quo)
{
int signx, signy, signres;
int mswx;
int mswy;
double x_over_y;
__int32_t n,hx,hy,hz,ix,iy,sx,i;
__uint32_t lx,ly,lz,q,sxy;
GET_HIGH_WORD(mswx, x);
GET_HIGH_WORD(mswy, y);
EXTRACT_WORDS(hx,lx,x);
EXTRACT_WORDS(hy,ly,y);
sxy = (hx ^ hy) & 0x80000000;
sx = hx&0x80000000; /* sign of x */
hx ^=sx; /* |x| */
hy &= 0x7fffffff; /* |y| */
signx = (mswx & 0x80000000) >> 31;
signy = (mswy & 0x80000000) >> 31;
/* purge off exception values */
if((hy|ly)==0||(hx>=0x7ff00000)|| /* y=0,or x not finite */
((hy|((ly|-ly)>>31))>0x7ff00000)) { /* or y is NaN */
*quo = 0; /* Not necessary, but return consistent value */
return (x*y)/(x*y);
}
if(hx<=hy) {
if((hx<hy)||(lx<ly)) {
q = 0;
goto fixup; /* |x|<|y| return x or x-y */
}
if(lx==ly) {
*quo = 1;
return Zero[(__uint32_t)sx>>31]; /* |x|=|y| return x*0*/
}
}
signres = (signx ^ signy) ? -1 : 1;
/* determine ix = ilogb(x) */
if(hx<0x00100000) { /* subnormal x */
if(hx==0) {
for (ix = -1043, i=lx; i>0; i<<=1) ix -=1;
} else {
for (ix = -1022,i=(hx<<11); i>0; i<<=1) ix -=1;
}
} else ix = (hx>>20)-1023;
x_over_y = fabs(x / y);
/* determine iy = ilogb(y) */
if(hy<0x00100000) { /* subnormal y */
if(hy==0) {
for (iy = -1043, i=ly; i>0; i<<=1) iy -=1;
} else {
for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
}
} else iy = (hy>>20)-1023;
*quo = signres * (lrint(x_over_y) & 0x7f);
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -1022)
hx = 0x00100000|(0x000fffff&hx);
else { /* subnormal x, shift x to normal */
n = -1022-ix;
if(n<=31) {
hx = (hx<<n)|(lx>>(32-n));
lx <<= n;
} else {
hx = lx<<(n-32);
lx = 0;
}
}
if(iy >= -1022)
hy = 0x00100000|(0x000fffff&hy);
else { /* subnormal y, shift y to normal */
n = -1022-iy;
if(n<=31) {
hy = (hy<<n)|(ly>>(32-n));
ly <<= n;
} else {
hy = ly<<(n-32);
ly = 0;
}
}
return remainder(x,y);
/* fix point fmod */
n = ix - iy;
q = 0;
while(n--) {
hz=hx-hy;lz=lx-ly; if(lx<ly) hz -= 1;
if(hz<0){hx = hx+hx+(lx>>31); lx = lx+lx;}
else {hx = hz+hz+(lz>>31); lx = lz+lz; q++;}
q <<= 1;
}
hz=hx-hy;lz=lx-ly; if(lx<ly) hz -= 1;
if(hz>=0) {hx=hz;lx=lz;q++;}
/* convert back to floating value and restore the sign */
if((hx|lx)==0) { /* return sign(x)*0 */
q &= QUO_MASK;
*quo = (sxy ? -q : q);
return Zero[(__uint32_t)sx>>31];
}
while(hx<0x00100000) { /* normalize x */
hx = hx+hx+(lx>>31); lx = lx+lx;
iy -= 1;
}
if(iy>= -1022) { /* normalize output */
hx = ((hx-0x00100000)|((iy+1023)<<20));
} else { /* subnormal output */
n = -1022 - iy;
if(n<=20) {
lx = (lx>>n)|((__uint32_t)hx<<(32-n));
hx >>= n;
} else if (n<=31) {
lx = (hx<<(32-n))|(lx>>n); hx = sx;
} else {
lx = hx>>(n-32); hx = sx;
}
}
fixup:
INSERT_WORDS(x,hx,lx);
y = fabs(y);
if (y < 0x1p-1021) {
if (x+x>y || (x+x==y && (q & 1))) {
q++;
x-=y;
}
} else if (x>0.5*y || (x==0.5*y && (q & 1))) {
q++;
x-=y;
}
GET_HIGH_WORD(hx,x);
SET_HIGH_WORD(x,hx^sx);
q &= QUO_MASK;
*quo = (sxy ? -q : q);
return x;
}
#endif /* defined(_DOUBLE_IS_32BITS) */

34
newlib/libm/common/s_rint.c
View File

@@ -10,6 +10,40 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<rint>>, <<rintf>>--round to integer
INDEX
rint
INDEX
rintf
ANSI_SYNOPSIS
#include <math.h>
double rint(double <[x]>);
float rintf(float <[x]>);
DESCRIPTION
The <<rint>> functions round their argument to an integer value in
floating-point format, using the current rounding direction. They
raise the "inexact" floating-point exception if the result differs
in value from the argument. See the <<nearbyint>> functions for the
same function with the "inexact" floating-point exception never being
raised. Newlib does not directly support floating-point exceptions.
The <<rint>> functions are written so that the "inexact" exception is
raised in hardware implementations that support it, even though Newlib
does not provide access.
RETURNS
<[x]> rounded to an integral value, using the current rounding direction.
PORTABILITY
ANSI C, POSIX
SEEALSO
<<nearbyint>>, <<round>>
*/
/*
* rint(x)

32
newlib/libm/common/s_round.c
View File

@@ -8,6 +8,38 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<round>>, <<roundf>>--round to integer, to nearest
INDEX
round
INDEX
roundf
ANSI_SYNOPSIS
#include <math.h>
double round(double <[x]>);
float roundf(float <[x]>);
DESCRIPTION
The <<round>> functions round their argument to the nearest integer
value in floating-point format, rounding halfway cases away from zero,
regardless of the current rounding direction. (While the "inexact"
floating-point exception behavior is unspecified by the C standard, the
<<round>> functions are written so that "inexact" is not raised if the
result does not equal the argument, which behavior is as recommended by
IEEE 754 for its related functions.)
RETURNS
<[x]> rounded to an integral value.
PORTABILITY
ANSI C, POSIX
SEEALSO
<<nearbyint>>, <<rint>>
*/
#include "fdlibm.h"

44
newlib/libm/common/s_scalbn.c
View File

@@ -13,38 +13,44 @@
/*
FUNCTION
<<scalbn>>, <<scalbnf>>---scale by power of two
<<scalbn>>, <<scalbnf>>, <<scalbln>>, <<scalblnf>>--scale by power of FLT_RADIX (=2)
INDEX
scalbn
INDEX
scalbnf
INDEX
scalbln
INDEX
scalblnf
ANSI_SYNOPSIS
#include <math.h>
double scalbn(double <[x]>, int <[y]>);
float scalbnf(float <[x]>, int <[y]>);
TRAD_SYNOPSIS
#include <math.h>
double scalbn(<[x]>,<[y]>)
double <[x]>;
int <[y]>;
float scalbnf(<[x]>,<[y]>)
float <[x]>;
int <[y]>;
double scalbn(double <[x]>, int <[n]>);
float scalbnf(float <[x]>, int <[n]>);
double scalbln(double <[x]>, long int <[n]>);
float scalblnf(float <[x]>, long int <[n]>);
DESCRIPTION
<<scalbn>> and <<scalbnf>> scale <[x]> by <[n]>, returning <[x]> times
2 to the power <[n]>. The result is computed by manipulating the
exponent, rather than by actually performing an exponentiation or
multiplication.
The <<scalbn>> and <<scalbln>> functions compute
@ifnottex
<[x]> times FLT_RADIX to the power <[n]>.
@end ifnottex
@tex
$x \cdot FLT\_RADIX^n$.
@end tex
efficiently. The result is computed by manipulating the exponent, rather than
by actually performing an exponentiation or multiplication. In this
floating-point implementation FLT_RADIX=2, which makes the <<scalbn>>
functions equivalent to the <<ldexp>> functions.
RETURNS
<[x]> times 2 to the power <[n]>.
<[x]> times 2 to the power <[n]>. A range error may occur.
PORTABILITY
Neither <<scalbn>> nor <<scalbnf>> is required by ANSI C or by the System V
Interface Definition (Issue 2).
ANSI C, POSIX
SEEALSO
<<ldexp>>
*/

29
newlib/libm/common/s_signbit.c
View File

@@ -3,6 +3,35 @@
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
*/
/*
FUNCTION
<<signbit>>--Does floating-point number have negative sign?
INDEX
signbit
ANSI_SYNOPSIS
#include <math.h>
int signbit(real-floating <[x]>);
DESCRIPTION
The <<signbit>> macro determines whether the sign of its argument value is
negative. The macro reports the sign of all values, including infinities,
zeros, and NaNs. If zero is unsigned, it is treated as positive. As shown in
the synopsis, the argument is "real-floating," meaning that any of the real
floating-point types (float, double, etc.) may be given to it.
Note that because of the possibilities of signed 0 and NaNs, the expression
"<[x]> < 0.0" does not give the same result as <<signbit>> in all cases.
RETURNS
The <<signbit>> macro returns a nonzero value if and only if the sign of its
argument value is negative.
PORTABILITY
C99, POSIX.
*/
#include "fdlibm.h"

29
newlib/libm/common/s_trunc.c
View File

@@ -8,6 +8,35 @@
* is preserved.
* ====================================================
*/
/*
FUNCTION
<<trunc>>, <<truncf>>--round to integer, towards zero
INDEX
trunc
INDEX
truncf
ANSI_SYNOPSIS
#include <math.h>
double trunc(double <[x]>);
float truncf(float <[x]>);
DESCRIPTION
The <<trunc>> functions round their argument to the integer value, in
floating format, nearest to but no larger in magnitude than the
argument, regardless of the current rounding direction. (While the
"inexact" floating-point exception behavior is unspecified by the C
standard, the <<trunc>> functions are written so that "inexact" is not
raised if the result does not equal the argument, which behavior is as
recommended by IEEE 754 for its related functions.)
RETURNS
<[x]> truncated to an integral value.
PORTABILITY
ANSI C, POSIX
*/
#include "fdlibm.h"

6
newlib/libm/common/sf_fdim.c
View File

@@ -15,8 +15,10 @@
#endif
{
int c = __fpclassifyf(x);
if (c == FP_NAN || c == FP_INFINITE)
return HUGE_VAL;
if (c == FP_NAN) return(x);
if (__fpclassifyf(y) == FP_NAN) return(y);
if (c == FP_INFINITE)
return HUGE_VALF;
return x > y ? x - y : 0.0;
}

10
newlib/libm/common/sf_fma.c
View File

@@ -15,8 +15,14 @@
float z;
#endif
{
/* Let the implementation handle this. */
return (x * y) + z;
/* NOTE: The floating-point exception behavior of this is not as
* required. But since the basic function is not really done properly,
* it is not worth bothering to get the exceptions right, either. */
/* Let the implementation handle this. */ /* <= NONSENSE! */
/* In floating-point implementations in which double is larger than float,
* computing as double should provide the desired function. Otherwise,
* the behavior will not be as specified in the standards. */
return (float) (((double) x * (double) y) + (double) z);
}
#ifdef _DOUBLE_IS_32BITS

55
newlib/libm/common/sf_llround.c Normal file
View File

@@ -0,0 +1,55 @@
/* lroundf adapted to be llroundf for Newlib, 2009 by Craig Howland. */
/*
* ====================================================
* 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.
* ====================================================
*/
#include "fdlibm.h"
long long int
llroundf(float x)
{
__int32_t exponent_less_127;
__uint32_t w;
long long int result;
__int32_t sign;
GET_FLOAT_WORD (w, x);
exponent_less_127 = ((w & 0x7f800000) >> 23) - 127;
sign = (w & 0x80000000) != 0 ? -1 : 1;
w &= 0x7fffff;
w |= 0x800000;
if (exponent_less_127 < (int)((8 * sizeof (long long int)) - 1))
{
if (exponent_less_127 < 0)
return exponent_less_127 < -1 ? 0 : sign;
else if (exponent_less_127 >= 23)
result = (long long int) w << (exponent_less_127 - 23);
else
{
w += 0x400000 >> exponent_less_127;
result = w >> (23 - exponent_less_127);
}
}
else
return (long long int) x;
return sign * result;
}
#ifdef _DOUBLE_IS_32BITS
long long int
llround(double x)
{
return llroundf((float) x);
}
#endif /* defined(_DOUBLE_IS_32BITS) */

46
newlib/libm/common/sf_logb.c
View File

@@ -1,7 +1,4 @@
/* sf_logb.c -- float version of s_logb.c.
* Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
*/
/* 2009 for Newlib: Sun's sf_ilogb.c converted to be sf_logb.c. */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
@@ -13,24 +10,41 @@
* ====================================================
*/
/* float logb(float x)
* return the binary exponent of non-zero x
* logbf(0) = -inf, raise divide-by-zero floating point exception
* logbf(+inf|-inf) = +inf (no signal is raised)
* logbf(NaN) = NaN (no signal is raised)
* Per C99 recommendation, a NaN argument is returned unchanged.
*/
#include "fdlibm.h"
float
#ifdef __STDC__
float logbf(float x)
logbf(float x)
#else
float logbf(x)
float x;
logbf(x)
float x;
#endif
{
__int32_t ix;
GET_FLOAT_WORD(ix,x);
ix &= 0x7fffffff; /* high |x| */
if(FLT_UWORD_IS_ZERO(ix)) return (float)-1.0/fabsf(x);
if(!FLT_UWORD_IS_FINITE(ix)) return x*x;
if((ix>>=23)==0) /* IEEE 754 logb */
return -126.0;
else
return (float) (ix-127);
__int32_t hx,ix;
GET_FLOAT_WORD(hx,x);
hx &= 0x7fffffff;
if(FLT_UWORD_IS_ZERO(hx)) {
float xx;
/* arg==0: return -inf and raise divide-by-zero exception */
SET_FLOAT_WORD(xx,hx); /* +0.0 */
return -1./xx; /* logbf(0) = -inf */
}
if(FLT_UWORD_IS_SUBNORMAL(hx)) {
for (ix = -126,hx<<=8; hx>0; hx<<=1) ix -=1;
return (float) ix;
}
else if (FLT_UWORD_IS_INFINITE(hx)) return HUGE_VALF; /* x==+|-inf */
else if (FLT_UWORD_IS_NAN(hx)) return x;
else return (float) ((hx>>23)-127);
}
#ifdef _DOUBLE_IS_32BITS

2
newlib/libm/common/sf_lrint.c
View File

@@ -95,7 +95,7 @@ TWO23[2]={
double x;
#endif
{
return (double) lrintf((float) x);
return lrintf((float) x);
}
#endif /* defined(_DOUBLE_IS_32BITS) */

2
newlib/libm/common/sf_lround.c
View File

@@ -56,7 +56,7 @@
double x;
#endif
{
return (double) lroundf((float) x);
return lroundf((float) x);
}
#endif /* defined(_DOUBLE_IS_32BITS) */

152
newlib/libm/common/sf_remquo.c
View File

@@ -1,50 +1,130 @@
/* Copyright (C) 2002 by Red Hat, Incorporated. All rights reserved.
/* Adapted for Newlib, 2009. (Allow for int < 32 bits; return *quo=0 during
* errors to make test scripts easier.) */
/* @(#)e_fmod.c 1.3 95/01/18 */
/*-
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this software
* is freely granted, provided that this notice is preserved.
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
#include <math.h>
#include "fdlibm.h"
#ifdef __STDC__
float remquof(float x, float y, int *quo) /* wrapper remquof */
#else
float remquof(x,y,quo) /* wrapper remquof */
float x,y;
int *quo;
/* For quotient, return either all 31 bits that can from calculation (using
* int32_t), or as many as can fit into an int that is smaller than 32 bits. */
#if INT_MAX > 0x7FFFFFFFL
#define QUO_MASK 0x7FFFFFFF
# else
#define QUO_MASK INT_MAX
#endif
static const float Zero[] = {0.0, -0.0,};
/*
* Return the IEEE remainder and set *quo to the last n bits of the
* quotient, rounded to the nearest integer. We choose n=31--if that many fit--
* we wind up computing all the integer bits of the quotient anyway as
* a side-effect of computing the remainder by the shift and subtract
* method. In practice, this is far more bits than are needed to use
* remquo in reduction algorithms.
*/
float
remquof(float x, float y, int *quo)
{
int signx, signy, signres;
int wx;
int wy;
float x_over_y;
__int32_t n,hx,hy,hz,ix,iy,sx,i;
__uint32_t q,sxy;
GET_FLOAT_WORD(wx, x);
GET_FLOAT_WORD(wy, y);
GET_FLOAT_WORD(hx,x);
GET_FLOAT_WORD(hy,y);
sxy = (hx ^ hy) & 0x80000000;
sx = hx&0x80000000; /* sign of x */
hx ^=sx; /* |x| */
hy &= 0x7fffffff; /* |y| */
signx = (wx & 0x80000000) >> 31;
signy = (wy & 0x80000000) >> 31;
/* purge off exception values */
if(hy==0||hx>=0x7f800000||hy>0x7f800000) { /* y=0,NaN;or x not finite */
*quo = 0; /* Not necessary, but return consistent value */
return (x*y)/(x*y);
}
if(hx<hy) {
q = 0;
goto fixup; /* |x|<|y| return x or x-y */
} else if(hx==hy) {
*quo = 1;
return Zero[(__uint32_t)sx>>31]; /* |x|=|y| return x*0*/
}
signres = (signx ^ signy) ? -1 : 1;
/* determine ix = ilogb(x) */
if(hx<0x00800000) { /* subnormal x */
for (ix = -126,i=(hx<<8); i>0; i<<=1) ix -=1;
} else ix = (hx>>23)-127;
x_over_y = fabsf(x / y);
/* determine iy = ilogb(y) */
if(hy<0x00800000) { /* subnormal y */
for (iy = -126,i=(hy<<8); i>0; i<<=1) iy -=1;
} else iy = (hy>>23)-127;
*quo = signres * (lrintf(x_over_y) & 0x7f);
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -126)
hx = 0x00800000|(0x007fffff&hx);
else { /* subnormal x, shift x to normal */
n = -126-ix;
hx <<= n;
}
if(iy >= -126)
hy = 0x00800000|(0x007fffff&hy);
else { /* subnormal y, shift y to normal */
n = -126-iy;
hy <<= n;
}
return remainderf(x,y);
/* fix point fmod */
n = ix - iy;
q = 0;
while(n--) {
hz=hx-hy;
if(hz<0) hx = hx << 1;
else {hx = hz << 1; q++;}
q <<= 1;
}
hz=hx-hy;
if(hz>=0) {hx=hz;q++;}
/* convert back to floating value and restore the sign */
if(hx==0) { /* return sign(x)*0 */
*quo = (sxy ? -q : q);
return Zero[(__uint32_t)sx>>31];
}
while(hx<0x00800000) { /* normalize x */
hx <<= 1;
iy -= 1;
}
if(iy>= -126) { /* normalize output */
hx = ((hx-0x00800000)|((iy+127)<<23));
} else { /* subnormal output */
n = -126 - iy;
hx >>= n;
}
fixup:
SET_FLOAT_WORD(x,hx);
y = fabsf(y);
if (y < 0x1p-125f) {
if (x+x>y || (x+x==y && (q & 1))) {
q++;
x-=y;
}
} else if (x>0.5f*y || (x==0.5f*y && (q & 1))) {
q++;
x-=y;
}
GET_FLOAT_WORD(hx,x);
SET_FLOAT_WORD(x,hx^sx);
q &= 0x7fffffff;
*quo = (sxy ? -q : q);
return x;
}
#ifdef _DOUBLE_IS_32BITS
#ifdef __STDC__
double remquo(double x, double y, int *quo) /* wrapper remquof */
#else
double remquo(x,y,quo) /* wrapper remquof */
double x,y;
int *quo;
#endif
{
return (double) remquof((float) x, (float) y, quo);
}
#endif /* defined(_DOUBLE_IS_32BITS) */

104
newlib/libm/math/Makefile.am
View File

@@ -65,14 +65,14 @@ endif # USE_LIBTOOL
include $(srcdir)/../../Makefile.shared
chobj = wacos.def wacosh.def wasin.def sasinh.def \
satan.def watan2.def watanh.def wj0.def \
wcosh.def serf.def wexp.def \
sfabs.def sfloor.def wfmod.def sfrexp.def \
wgamma.def whypot.def sldexp.def wlog.def \
wlog10.def \
wpow.def wremainder.def ssin.def wsinh.def \
wsqrt.def stan.def stanh.def
chobj = w_acos.def w_acosh.def w_asin.def s_asinh.def \
s_atan.def w_atan2.def w_atanh.def w_j0.def \
w_cosh.def s_erf.def w_exp.def w_exp2.def \
s_fabs.def s_floor.def w_fmod.def s_frexp.def \
w_gamma.def w_hypot.def s_ldexp.def w_log.def \
w_log10.def \
w_pow.def w_remainder.def s_sin.def w_sinh.def \
w_sqrt.def s_tan.def s_tanh.def
SUFFIXES = .def
@@ -89,94 +89,6 @@ doc: $(chobj)
CLEANFILES = $(chobj) *.ref
# Texinfo does not appear to support underscores in file names, so we
# name the .def files without underscores.
wacos.def: w_acos.c
$(CHEW) < $(srcdir)/w_acos.c >$@ 2>/dev/null
touch stmp-def
wacosh.def: w_acosh.c
$(CHEW) < $(srcdir)/w_acosh.c >$@ 2>/dev/null
touch stmp-def
wasin.def: w_asin.c
$(CHEW) < $(srcdir)/w_asin.c >$@ 2>/dev/null
touch stmp-def
sasinh.def: s_asinh.c
$(CHEW) < $(srcdir)/s_asinh.c >$@ 2>/dev/null
touch stmp-def
satan.def: s_atan.c
$(CHEW) < $(srcdir)/s_atan.c >$@ 2>/dev/null
touch stmp-def
watan2.def: w_atan2.c
$(CHEW) < $(srcdir)/w_atan2.c >$@ 2>/dev/null
touch stmp-def
watanh.def: w_atanh.c
$(CHEW) < $(srcdir)/w_atanh.c >$@ 2>/dev/null
touch stmp-def
wj0.def: w_j0.c
$(CHEW) < $(srcdir)/w_j0.c >$@ 2>/dev/null
touch stmp-def
scopysign.def: s_copysign.c
$(CHEW) < $(srcdir)/../common/s_copysign.c >$@ 2>/dev/null
touch stmp-def
wcosh.def: w_cosh.c
$(CHEW) < $(srcdir)/w_cosh.c >$@ 2>/dev/null
touch stmp-def
serf.def: s_erf.c
$(CHEW) < $(srcdir)/s_erf.c >$@ 2>/dev/null
touch stmp-def
wexp.def: w_exp.c
$(CHEW) < $(srcdir)/w_exp.c >$@ 2>/dev/null
touch stmp-def
sfabs.def: s_fabs.c
$(CHEW) < $(srcdir)/s_fabs.c >$@ 2>/dev/null
touch stmp-def
sfloor.def: s_floor.c
$(CHEW) < $(srcdir)/s_floor.c >$@ 2>/dev/null
touch stmp-def
wfmod.def: w_fmod.c
$(CHEW) < $(srcdir)/w_fmod.c >$@ 2>/dev/null
touch stmp-def
sfrexp.def: s_frexp.c
$(CHEW) < $(srcdir)/s_frexp.c >$@ 2>/dev/null
touch stmp-def
wgamma.def: w_gamma.c
$(CHEW) < $(srcdir)/w_gamma.c >$@ 2>/dev/null
touch stmp-def
whypot.def: w_hypot.c
$(CHEW) < $(srcdir)/w_hypot.c >$@ 2>/dev/null
touch stmp-def
sldexp.def: s_ldexp.c
$(CHEW) < $(srcdir)/s_ldexp.c >$@ 2>/dev/null
touch stmp-def
wlog.def: w_log.c
$(CHEW) < $(srcdir)/w_log.c >$@ 2>/dev/null
touch stmp-def
wlog10.def: w_log10.c
$(CHEW) < $(srcdir)/w_log10.c >$@ 2>/dev/null
touch stmp-def
wpow.def: w_pow.c
$(CHEW) < $(srcdir)/w_pow.c >$@ 2>/dev/null
touch stmp-def
wremainder.def: w_remainder.c
$(CHEW) < $(srcdir)/w_remainder.c >$@ 2>/dev/null
touch stmp-def
ssin.def: s_sin.c
$(CHEW) < $(srcdir)/s_sin.c >$@ 2>/dev/null
touch stmp-def
wsinh.def: w_sinh.c
$(CHEW) < $(srcdir)/w_sinh.c >$@ 2>/dev/null
touch stmp-def
wsqrt.def: w_sqrt.c
$(CHEW) < $(srcdir)/w_sqrt.c >$@ 2>/dev/null
touch stmp-def
stan.def: s_tan.c
$(CHEW) < $(srcdir)/s_tan.c >$@ 2>/dev/null
touch stmp-def
stanh.def: s_tanh.c
$(CHEW) < $(srcdir)/s_tanh.c >$@ 2>/dev/null
touch stmp-def
# A partial dependency list.
$(lib_a_OBJECTS): $(srcdir)/../../libc/include/math.h $(srcdir)/../common/fdlibm.h

104
newlib/libm/math/Makefile.in
View File

@@ -373,14 +373,14 @@ libmath_la_LDFLAGS = -Xcompiler -nostdlib
@USE_LIBTOOL_FALSE@noinst_LIBRARIES = lib.a
@USE_LIBTOOL_FALSE@lib_a_SOURCES = $(src) $(fsrc)
@USE_LIBTOOL_FALSE@lib_a_CFLAGS = $(AM_CFLAGS)
chobj = wacos.def wacosh.def wasin.def sasinh.def \
satan.def watan2.def watanh.def wj0.def \
wcosh.def serf.def wexp.def \
sfabs.def sfloor.def wfmod.def sfrexp.def \
wgamma.def whypot.def sldexp.def wlog.def \
wlog10.def \
wpow.def wremainder.def ssin.def wsinh.def \
wsqrt.def stan.def stanh.def
chobj = w_acos.def w_acosh.def w_asin.def s_asinh.def \
s_atan.def w_atan2.def w_atanh.def w_j0.def \
w_cosh.def s_erf.def w_exp.def w_exp2.def \
s_fabs.def s_floor.def w_fmod.def s_frexp.def \
w_gamma.def w_hypot.def s_ldexp.def w_log.def \
w_log10.def \
w_pow.def w_remainder.def s_sin.def w_sinh.def \
w_sqrt.def s_tan.def s_tanh.def
SUFFIXES = .def
CHEW = ../../doc/makedoc -f $(srcdir)/../../doc/doc.str
@@ -1424,94 +1424,6 @@ objectlist.awk.in: $(noinst_LTLIBRARIES)
doc: $(chobj)
cat $(srcdir)/math.tex >> $(TARGETDOC)
# Texinfo does not appear to support underscores in file names, so we
# name the .def files without underscores.
wacos.def: w_acos.c
$(CHEW) < $(srcdir)/w_acos.c >$@ 2>/dev/null
touch stmp-def
wacosh.def: w_acosh.c
$(CHEW) < $(srcdir)/w_acosh.c >$@ 2>/dev/null
touch stmp-def
wasin.def: w_asin.c
$(CHEW) < $(srcdir)/w_asin.c >$@ 2>/dev/null
touch stmp-def
sasinh.def: s_asinh.c
$(CHEW) < $(srcdir)/s_asinh.c >$@ 2>/dev/null
touch stmp-def
satan.def: s_atan.c
$(CHEW) < $(srcdir)/s_atan.c >$@ 2>/dev/null
touch stmp-def
watan2.def: w_atan2.c
$(CHEW) < $(srcdir)/w_atan2.c >$@ 2>/dev/null
touch stmp-def
watanh.def: w_atanh.c
$(CHEW) < $(srcdir)/w_atanh.c >$@ 2>/dev/null
touch stmp-def
wj0.def: w_j0.c
$(CHEW) < $(srcdir)/w_j0.c >$@ 2>/dev/null
touch stmp-def
scopysign.def: s_copysign.c
$(CHEW) < $(srcdir)/../common/s_copysign.c >$@ 2>/dev/null
touch stmp-def
wcosh.def: w_cosh.c
$(CHEW) < $(srcdir)/w_cosh.c >$@ 2>/dev/null
touch stmp-def
serf.def: s_erf.c
$(CHEW) < $(srcdir)/s_erf.c >$@ 2>/dev/null
touch stmp-def
wexp.def: w_exp.c
$(CHEW) < $(srcdir)/w_exp.c >$@ 2>/dev/null
touch stmp-def
sfabs.def: s_fabs.c
$(CHEW) < $(srcdir)/s_fabs.c >$@ 2>/dev/null
touch stmp-def
sfloor.def: s_floor.c
$(CHEW) < $(srcdir)/s_floor.c >$@ 2>/dev/null
touch stmp-def
wfmod.def: w_fmod.c
$(CHEW) < $(srcdir)/w_fmod.c >$@ 2>/dev/null
touch stmp-def
sfrexp.def: s_frexp.c
$(CHEW) < $(srcdir)/s_frexp.c >$@ 2>/dev/null
touch stmp-def
wgamma.def: w_gamma.c
$(CHEW) < $(srcdir)/w_gamma.c >$@ 2>/dev/null
touch stmp-def
whypot.def: w_hypot.c
$(CHEW) < $(srcdir)/w_hypot.c >$@ 2>/dev/null
touch stmp-def
sldexp.def: s_ldexp.c
$(CHEW) < $(srcdir)/s_ldexp.c >$@ 2>/dev/null
touch stmp-def
wlog.def: w_log.c
$(CHEW) < $(srcdir)/w_log.c >$@ 2>/dev/null
touch stmp-def
wlog10.def: w_log10.c
$(CHEW) < $(srcdir)/w_log10.c >$@ 2>/dev/null
touch stmp-def
wpow.def: w_pow.c
$(CHEW) < $(srcdir)/w_pow.c >$@ 2>/dev/null
touch stmp-def
wremainder.def: w_remainder.c
$(CHEW) < $(srcdir)/w_remainder.c >$@ 2>/dev/null
touch stmp-def
ssin.def: s_sin.c
$(CHEW) < $(srcdir)/s_sin.c >$@ 2>/dev/null
touch stmp-def
wsinh.def: w_sinh.c
$(CHEW) < $(srcdir)/w_sinh.c >$@ 2>/dev/null
touch stmp-def
wsqrt.def: w_sqrt.c
$(CHEW) < $(srcdir)/w_sqrt.c >$@ 2>/dev/null
touch stmp-def
stan.def: s_tan.c
$(CHEW) < $(srcdir)/s_tan.c >$@ 2>/dev/null
touch stmp-def
stanh.def: s_tanh.c
$(CHEW) < $(srcdir)/s_tanh.c >$@ 2>/dev/null
touch stmp-def
# A partial dependency list.
$(lib_a_OBJECTS): $(srcdir)/../../libc/include/math.h $(srcdir)/../common/fdlibm.h

229
newlib/libm/math/math.tex
View File

@@ -9,6 +9,7 @@ Two definitions from @file{math.h} are of particular interest.
@item
The representation of infinity as a @code{double} is defined as
@code{HUGE_VAL}; this number is returned on overflow by many functions.
The macro @code{HUGE_VALF} is a corresponding value for @code{float}.
@item
The structure @code{exception} is used when you write customized error
@@ -44,43 +45,59 @@ machines---are available when you include @file{fastmath.h} instead of
* atan:: Arctangent
* atan2:: Arctangent of y/x
* atanh:: Inverse hyperbolic tangent
* jN:: Bessel functions (jN, yN)
* jN:: Bessel functions (jN, yN)
* cbrt:: Cube root
* copysign:: Sign of Y, magnitude of X
* cosh:: Hyperbolic cosine
* erf:: Error function (erf, erfc)
* exp:: Exponential
* expm1:: Exponential of x, - 1
* erf:: Error function (erf, erfc)
* exp:: Exponential, base e
* exp2:: Exponential, base 2
* expm1:: Exponential, base e, of x - 1
* fabs:: Absolute value (magnitude)
* floor:: Floor and ceiling (floor, ceil)
* fdim:: Positive difference
* floor:: Floor and ceiling (floor, ceil)
* fma:: Floating multiply add
* fmax:: Maximum
* fmin:: Minimum
* fmod:: Floating-point remainder (modulo)
* fpclassify:: Floating-point classification macro
* frexp:: Split floating-point number
* gamma:: Logarithmic gamma function
* hypot:: Distance from origin
* ilogb:: Get exponent
* infinity:: Floating infinity
* isnan:: Check type of number
* ldexp:: Load exponent
* isgreater:: Comparison macros
* ldexp:: Scale by a power of 2
* log:: Natural logarithms
* log10:: Base 10 logarithms
* log1p:: Log of 1 + X
* log2:: Base 2 logarithms
* logb:: Get exponent
* lrint:: Round to integer
* lround:: Round to integer, away from zero (lround, llround)
* matherr:: Modifiable math error handler
* modf:: Split fractional and integer parts
* nan:: Floating Not a Number
* nearbyint:: Round to integer
* nextafter:: Get next representable number
* pow:: X to the power Y
* remainder:: remainder of X divided by Y
* scalbn:: scalbn
* remquo:: Remainder and part of quotient
* rint:: Round to integer
* round:: Round to integer, away from zero
* scalbn:: Scale by a power of FLT_RADIX (2)
* signbit:: Does floating-point number have negative sign?
* sin:: Sine or cosine (sin, cos)
* sinh:: Hyperbolic sine
* sqrt:: Positive square root
* tan:: Tangent
* tanh:: Hyperbolic tangent
* trunc:: Round to integer, towards zero
@end menu
@page
@node version
@section Version of library
@section Error Handling
There are four different versions of the math library routines: IEEE,
POSIX, X/Open, or SVID. The version may be selected at runtime by
@@ -114,119 +131,89 @@ log: DOMAIN error
The library is set to X/Open mode by default.
@page
@include math/wacos.def
The aforementioned error reporting is the supported Newlib libm error
handling method. However, the majority of the functions are written
so as to produce the floating-point exceptions (e.g. "invalid",
"divide-by-zero") as required by the C and POSIX standards, for
floating-point implementations that support them. Newlib does not provide
the floating-point exception access routines defined in the standards
for fenv.h, though, which is why they are considered unsupported. It is
mentioned in case you have separately-provided access routines so that
you are aware that they can be caused.
@page
@include math/wacosh.def
@section Standards Compliance And Portability
Most of the individual function descriptions describe the standards to which
each function complies. However, these descriptions are mostly out of date,
having been written before C99 was released. One of these days we'll get
around to updating the rest of them. (If you'd like to help, please let us
know.)
@page
@include math/wasin.def
``C99'' refers to ISO/IEC 9899:1999, ``Programming languages--C''.
``POSIX'' refers to IEEE Standard 1003.1. POSIX@registeredsymbol{} is a
registered trademark of The IEEE.
@page
@include math/sasinh.def
@c To sort the include list easily, keep the indentation right because want to
@c skip the s_|w_ at the start of most--but not all--of the file names.
@c (e.g., isgreater.def does not have a leading s nor w.) Then, sort
@c based on the column. For example: "sort -t@ -k3.17"
@c A few hand-edits might be appropriate after a sort, although not necessary
@c and are a nuisance as ought to be kept in sync with menu list above:
@c atan2 after atan, exp2 after exp, log first in log list, and w_j0 to place
@c to reflect function name of Bessel (as opposed to j; e.g. after atanh,
@c before cbrt).
@page
@include math/satan.def
@page
@include math/watan2.def
@page
@include math/watanh.def
@page
@include math/wj0.def
@page
@include common/scbrt.def
@page
@include common/scopysign.def
@page
@include math/wcosh.def
@page
@include math/serf.def
@page
@include math/wexp.def
@page
@include common/sexpm1.def
@page
@include math/sfabs.def
@page
@include math/sfloor.def
@page
@include math/wfmod.def
@page
@include math/sfrexp.def
@page
@include math/wgamma.def
@page
@include math/whypot.def
@page
@include common/silogb.def
@page
@include common/sinfinity.def
@page
@include common/sisnan.def
@page
@include math/sldexp.def
@page
@include math/wlog.def
@page
@include math/wlog10.def
@page
@include common/slog1p.def
@page
@include common/smatherr.def
@page
@include common/smodf.def
@page
@include common/snan.def
@page
@include common/snextafter.def
@page
@include math/wpow.def
@page
@include math/wremainder.def
@page
@include common/sscalbn.def
@page
@include math/wsqrt.def
@page
@include math/ssin.def
@page
@include math/wsinh.def
@page
@include math/stan.def
@page
@include math/stanh.def
@page @include math/w_acos.def
@page @include math/w_acosh.def
@page @include math/w_asin.def
@page @include math/s_asinh.def
@page @include math/s_atan.def
@page @include math/w_atan2.def
@page @include math/w_atanh.def
@page @include math/w_j0.def
@page @include common/s_cbrt.def
@page @include common/s_copysign.def
@page @include math/w_cosh.def
@page @include math/s_erf.def
@page @include math/w_exp.def
@page @include math/w_exp2.def
@page @include common/s_expm1.def
@page @include math/s_fabs.def
@page @include common/s_fdim.def
@page @include math/s_floor.def
@page @include common/s_fma.def
@page @include common/s_fmax.def
@page @include common/s_fmin.def
@page @include math/w_fmod.def
@page @include math/s_frexp.def
@page @include math/w_gamma.def
@page @include math/w_hypot.def
@page @include common/s_ilogb.def
@page @include common/s_infinity.def
@page @include common/isgreater.def
@page @include common/s_isnan.def
@page @include math/s_ldexp.def
@page @include math/w_log.def
@page @include math/w_log10.def
@page @include common/s_log1p.def
@page @include common/s_log2.def
@page @include common/s_logb.def
@page @include common/s_lrint.def
@page @include common/s_lround.def
@page @include common/s_matherr.def
@page @include common/s_modf.def
@page @include common/s_nan.def
@page @include common/s_nearbyint.def
@page @include common/s_nextafter.def
@page @include math/w_pow.def
@page @include math/w_remainder.def
@page @include common/s_remquo.def
@page @include common/s_rint.def
@page @include common/s_round.def
@page @include common/s_scalbn.def
@page @include common/s_signbit.def
@page @include math/s_sin.def
@page @include math/w_sinh.def
@page @include math/w_sqrt.def
@page @include math/s_tan.def
@page @include math/s_tanh.def
@page @include common/s_trunc.def

13
newlib/libm/math/w_exp2.c
View File

@@ -13,7 +13,7 @@
/*
FUNCTION
<<exp2>>, <<exp2f>>---exponential
<<exp2>>, <<exp2f>>--exponential, base 2
INDEX
exp2
INDEX
@@ -24,14 +24,6 @@ ANSI_SYNOPSIS
double exp2(double <[x]>);
float exp2f(float <[x]>);
TRAD_SYNOPSIS
#include <math.h>
double exp2(<[x]>);
double <[x]>;
float exp2f(<[x]>);
float <[x]>;
DESCRIPTION
<<exp2>> and <<exp2f>> calculate 2 ^ <[x]>, that is,
@ifnottex
@@ -50,6 +42,9 @@ RETURNS
result overflows, the returned value is <<HUGE_VAL>>. In
either case, <<errno>> is set to <<ERANGE>>.
PORTABILITY
ANSI C, POSIX.
*/
/*

53
newlib/libm/math/w_gamma.c
View File

@@ -12,11 +12,20 @@
*
*/
/* BUG: FIXME?
According to Linux man pages for tgamma, lgamma, and gamma, the gamma
function was originally defined in BSD as implemented here--the log of the gamma
function. BSD 4.3 changed the name to lgamma, apparently removing gamma. BSD
4.4 re-introduced the gamma name with the more intuitive, without logarithm,
plain gamma function. The C99 standard apparently wanted to avoid a problem
with the poorly-named earlier gamma and used tgamma when adding a plain
gamma function.
So the current gamma is matching an old, bad definition, and not
matching a newer, better definition. */
/*
FUNCTION
<<gamma>>, <<gammaf>>, <<lgamma>>, <<lgammaf>>, <<gamma_r>>,
<<gammaf_r>>, <<lgamma_r>>, <<lgammaf_r>>---logarithmic gamma
function
<<gamma>>, <<gammaf>>, <<lgamma>>, <<lgammaf>>, <<gamma_r>>, <<gammaf_r>>, <<lgamma_r>>, <<lgammaf_r>>, <<tgamma>>, and <<tgammaf>>--logarithmic and plain gamma functions
INDEX
gamma
INDEX
@@ -33,6 +42,10 @@ INDEX
lgamma_r
INDEX
lgammaf_r
INDEX
tgamma
INDEX
tgammaf
ANSI_SYNOPSIS
#include <math.h>
@@ -44,6 +57,8 @@ double gamma_r(double <[x]>, int *<[signgamp]>);
float gammaf_r(float <[x]>, int *<[signgamp]>);
double lgamma_r(double <[x]>, int *<[signgamp]>);
float lgammaf_r(float <[x]>, int *<[signgamp]>);
double tgamma(double <[x]>);
float tgammaf(float <[x]>);
TRAD_SYNOPSIS
#include <math.h>
@@ -67,15 +82,19 @@ int <[signgamp]>;
float lgammaf_r(<[x]>, <[signgamp]>)
float <[x]>;
int <[signgamp]>;
double tgamma(<[x]>)
double <[x]>;
float tgammaf(<[x]>)
float <[x]>;
DESCRIPTION
<<gamma>> calculates
<<gamma>> calculates
@tex
$\mit ln\bigl(\Gamma(x)\bigr)$,
$\mit ln\bigl(\Gamma(x)\bigr)$,
@end tex
the natural logarithm of the gamma function of <[x]>. The gamma function
(<<exp(gamma(<[x]>))>>) is a generalization of factorial, and retains
the property that
the property that
@ifnottex
<<exp(gamma(N))>> is equivalent to <<N*exp(gamma(N-1))>>.
@end ifnottex
@@ -83,13 +102,13 @@ the property that
$\mit \Gamma(N)\equiv N\times\Gamma(N-1)$.
@end tex
Accordingly, the results of the gamma function itself grow very
quickly. <<gamma>> is defined as
quickly. <<gamma>> is defined as
@tex
$\mit ln\bigl(\Gamma(x)\bigr)$ rather than simply $\mit \Gamma(x)$
@end tex
@ifnottex
the natural log of the gamma function, rather than the gamma function
itself,
itself,
@end ifnottex
to extend the useful range of results representable.
@@ -113,8 +132,17 @@ variable <<signgam>> is not used. These functions may be used for
reentrant calls (but they will still set the global variable <<errno>>
if an error occurs).
<<tgamma>> and <<tgammaf>> are the "true gamma" functions, returning
@tex
$\mit \Gamma(x)$,
@end tex
the gamma function of <[x]>--without a logarithm.
(They are apparently so named because of the prior existence of the old,
poorly-named <<gamma>> functions which returned the log of gamma up
through BSD 4.2.)
RETURNS
Normally, the computed result is returned.
Normally, the computed result is returned.
When <[x]> is a nonpositive integer, <<gamma>> returns <<HUGE_VAL>>
and <<errno>> is set to <<EDOM>>. If the result overflows, <<gamma>>
@@ -123,7 +151,12 @@ returns <<HUGE_VAL>> and <<errno>> is set to <<ERANGE>>.
You can modify this error treatment using <<matherr>>.
PORTABILITY
Neither <<gamma>> nor <<gammaf>> is ANSI C. */
Neither <<gamma>> nor <<gammaf>> is ANSI C. It is better not to use either
of these; use <<lgamma>> or <<tgamma>> instead.@*
<<lgamma>>, <<lgammaf>>, <<tgamma>>, and <<tgammaf>> are nominally C standard
in terms of the base return values, although the <<matherr>> error-handling
is not standard, nor is the <[signgam]> global for <<lgamma>>.
*/
/* double gamma(double x)
* Return the logarithm of the Gamma function of x.