1.\" Copyright (c) 1985 Regents of the University of California. 2.\" All rights reserved. 3.\" 4.\" Redistribution and use in source and binary forms, with or without 5.\" modification, are permitted provided that the following conditions 6.\" are met: 7.\" 1. Redistributions of source code must retain the above copyright 8.\" notice, this list of conditions and the following disclaimer. 9.\" 2. Redistributions in binary form must reproduce the above copyright 10.\" notice, this list of conditions and the following disclaimer in the 11.\" documentation and/or other materials provided with the distribution. 12.\" 3. All advertising materials mentioning features or use of this software 13.\" must display the following acknowledgement: 14.\" This product includes software developed by the University of 15.\" California, Berkeley and its contributors. 16.\" 4. Neither the name of the University nor the names of its contributors 17.\" may be used to endorse or promote products derived from this software 18.\" without specific prior written permission. 19.\" 20.\" THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND 21.\" ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 22.\" IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 23.\" ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE 24.\" FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 25.\" DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 26.\" OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 27.\" HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 28.\" LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 29.\" OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 30.\" SUCH DAMAGE. 31.\" 32.\" from: @(#)math.3 6.10 (Berkeley) 5/6/91 33.\" $FreeBSD$ 34.\" 35.Dd January 11, 2005 36.Dt MATH 3 37.Os 38.if n \{\ 39.char \[sr] "sqrt 40.\} 41.Sh NAME 42.Nm math 43.Nd "floating-point mathematical library" 44.Sh LIBRARY 45.Lb libm 46.Sh SYNOPSIS 47.In math.h 48.Sh DESCRIPTION 49These functions constitute the C math library. 50.Sh "LIST OF FUNCTIONS" 51Each of the following 52.Vt double 53functions has a 54.Vt float 55counterpart with an 56.Ql f 57appended to the name and a 58.Vt "long double" 59counterpart with an 60.Ql l 61appended. 62As an example, the 63.Vt float 64and 65.Vt "long double" 66counterparts of 67.Ft double 68.Fn acos "double x" 69are 70.Ft float 71.Fn acosf "float x" 72and 73.Ft "long double" 74.Fn acosl "long double x" , 75respectively. 76.Pp 77The programs are accurate to within the numbers 78of 79.Em ulp Ns s 80tabulated below; an 81.Em ulp 82is one 83.Em U Ns nit 84in the 85.Em L Ns ast 86.Em P Ns lace . 87.Bl -column "nexttoward" "remainder with partial quotient" 88.Em "Name Description Error Bound (ULPs)" 89.\" XXX Many of these error bounds are wrong for the current implementation! 90acos inverse trigonometric function ??? 91acosh inverse hyperbolic function ??? 92asin inverse trigonometric function ??? 93asinh inverse hyperbolic function ??? 94atan inverse trigonometric function ??? 95atanh inverse hyperbolic function ??? 96atan2 inverse trigonometric function ??? 97cbrt cube root 1 98ceil integer no less than 0 99copysign copy sign bit 0 100cos trigonometric function 1 101cosh hyperbolic function ??? 102erf error function 1 103erfc complementary error function 1 104exp exponential base e 1 105.\" exp2 exponential base 2 ??? 106expm1 exp(x)\-1 1 107fabs absolute value 0 108fdim positive difference 1 109floor integer no greater than 0 110.\" fma multiply-add ??? 111fmax maximum function 0 112fmin minimum function 0 113fmod remainder function ??? 114frexp extract mantissa and exponent 0 115hypot Euclidean distance 1 116ilogb exponent extraction 0 117j0 bessel function ??? 118j1 bessel function ??? 119jn bessel function ??? 120ldexp multiply by power of 2 0 121lgamma log gamma function 1 122llrint round to integer 0 123llround round to nearest integer 0 124log natural logarithm 1 125log10 logarithm to base 10 1 126log1p log(1+x) 1 127.\" log2 base 2 logarithm 0 128logb exponent extraction 0 129lrint round to integer 0 130lround round to nearest integer 0 131modf extract fractional part 0 132.\" nan return quiet \*(Na) 0 133nearbyint round to integer 0 134nextafter next representable value 0 135.\" nexttoward next representable value 0 136pow exponential x**y 60-500 137remainder remainder 0 138.\" remquo remainder with partial quotient ??? 139rint round to nearest integer 0 140round round to nearest integer 0 141scalbln exponent adjustment 0 142scalbn exponent adjustment 0 143sin trigonometric function 1 144sinh hyperbolic function ??? 145sqrt square root 1 146tan trigonometric function 1 147tanh hyperbolic function ??? 148tgamma gamma function 1 149trunc round towards zero 0 150y0 bessel function ??? 151y1 bessel function ??? 152yn bessel function ??? 153.El 154.Sh NOTES 155Virtually all modern floating-point units attempt to support 156IEEE Standard 754 for Binary Floating-Point Arithmetic. 157This standard does not cover particular routines in the math library 158except for the few documented in 159.Xr ieee 3 ; 160it primarily defines representations of numbers and abstract 161properties of arithmetic operations relating to precision, rounding, 162and exceptional cases, as described below. 163.Ss IEEE STANDARD 754 Floating-Point Arithmetic 164.\" XXX mention single- and extended-/quad- precisions 165Properties of IEEE 754 Double-Precision: 166.Bd -ragged -offset indent -compact 167Wordsize: 64 bits, 8 bytes. 168.Pp 169Radix: Binary. 170.Pp 171Precision: 53 significant bits, 172roughly like 16 significant decimals. 173.Bd -ragged -offset indent -compact 174If x and x' are consecutive positive Double-Precision 175numbers (they differ by 1 176.Em ulp ) , 177then 178.Bd -ragged -compact 1791.1e\-16 < 0.5**53 < (x'\-x)/x \(<= 0.5**52 < 2.3e\-16. 180.Ed 181.Ed 182.Pp 183.Bl -column "XXX" -compact 184Range: Overflow threshold = 2.0**1024 = 1.8e308 185 Underflow threshold = 0.5**1022 = 2.2e\-308 186.El 187.Bd -ragged -offset indent -compact 188Overflow goes by default to a signed \*(If. 189Underflow is 190.Em Gradual , 191rounding to the nearest 192integer multiple of 0.5**1074 = 4.9e\-324. 193.Ed 194.Pp 195Zero is represented ambiguously as +0 or \-0. 196.Bd -ragged -offset indent -compact 197Its sign transforms correctly through multiplication or 198division, and is preserved by addition of zeros 199with like signs; but x\-x yields +0 for every 200finite x. 201The only operations that reveal zero's 202sign are division by zero and 203.Fn copysign x \(+-0 . 204In particular, comparison (x > y, x \(>= y, etc.)\& 205cannot be affected by the sign of zero; but if 206finite x = y then \*(If = 1/(x\-y) \(!= \-1/(y\-x) = \-\*(If. 207.Ed 208.Pp 209\*(If is signed. 210.Bd -ragged -offset indent -compact 211It persists when added to itself 212or to any finite number. 213Its sign transforms 214correctly through multiplication and division, and 215(finite)/\(+-\*(If\0=\0\(+-0 216(nonzero)/0 = \(+-\*(If. 217But 218\*(If\-\*(If, \*(If\(**0 and \*(If/\*(If 219are, like 0/0 and sqrt(\-3), 220invalid operations that produce \*(Na. ... 221.Ed 222.Pp 223Reserved operands: 224.Bd -ragged -offset indent -compact 225there are 2**53\-2 of them, all 226called \*(Na 227.Em ( N Ns ot Em a N Ns umber ) . 228Some, called Signaling \*(Nas, trap any floating-point operation 229performed upon them; they are used to mark missing 230or uninitialized values, or nonexistent elements 231of arrays. 232The rest are Quiet \*(Nas; they are 233the default results of Invalid Operations, and 234propagate through subsequent arithmetic operations. 235If x \(!= x then x is \*(Na; every other predicate 236(x > y, x = y, x < y, ...) is FALSE if \*(Na is involved. 237.Pp 238NOTE: Trichotomy is violated by \*(Na. 239Besides being FALSE, predicates that entail ordered 240comparison, rather than mere (in)equality, 241signal Invalid Operation when \*(Na is involved. 242.Ed 243.Pp 244Rounding: 245.Bd -ragged -offset indent -compact 246Every algebraic operation (+, \-, \(**, /, 247\(sr) 248is rounded by default to within half an 249.Em ulp , 250and when the rounding error is exactly half an 251.Em ulp 252then 253the rounded value's least significant bit is zero. 254This kind of rounding is usually the best kind, 255sometimes provably so; for instance, for every 256x = 1.0, 2.0, 3.0, 4.0, ..., 2.0**52, we find 257(x/3.0)\(**3.0 == x and (x/10.0)\(**10.0 == x and ... 258despite that both the quotients and the products 259have been rounded. 260Only rounding like IEEE 754 can do that. 261But no single kind of rounding can be 262proved best for every circumstance, so IEEE 754 263provides rounding towards zero or towards 264+\*(If or towards \-\*(If 265at the programmer's option. 266And the 267same kinds of rounding are specified for 268Binary-Decimal Conversions, at least for magnitudes 269between roughly 1.0e\-10 and 1.0e37. 270.Ed 271.Pp 272Exceptions: 273.Bd -ragged -offset indent -compact 274IEEE 754 recognizes five kinds of floating-point exceptions, 275listed below in declining order of probable importance. 276.Bl -column -offset indent "Invalid Operation" "Gradual Underflow" 277.Em "Exception Default Result" 278Invalid Operation \*(Na, or FALSE 279Overflow \(+-\*(If 280Divide by Zero \(+-\*(If 281Underflow Gradual Underflow 282Inexact Rounded value 283.El 284.Pp 285NOTE: An Exception is not an Error unless handled 286badly. 287What makes a class of exceptions exceptional 288is that no single default response can be satisfactory 289in every instance. 290On the other hand, if a default 291response will serve most instances satisfactorily, 292the unsatisfactory instances cannot justify aborting 293computation every time the exception occurs. 294.Ed 295.Pp 296For each kind of floating-point exception, IEEE 754 297provides a Flag that is raised each time its exception 298is signaled, and stays raised until the program resets 299it. 300Programs may also test, save and restore a flag. 301Thus, IEEE 754 provides three ways by which programs 302may cope with exceptions for which the default result 303might be unsatisfactory: 304.Bl -enum 305.It 306Test for a condition that might cause an exception 307later, and branch to avoid the exception. 308.It 309Test a flag to see whether an exception has occurred 310since the program last reset its flag. 311.It 312Test a result to see whether it is a value that only 313an exception could have produced. 314.Pp 315CAUTION: The only reliable ways to discover 316whether Underflow has occurred are to test whether 317products or quotients lie closer to zero than the 318underflow threshold, or to test the Underflow 319flag. 320(Sums and differences cannot underflow in 321IEEE 754; if x \(!= y then x\-y is correct to 322full precision and certainly nonzero regardless of 323how tiny it may be.) 324Products and quotients that 325underflow gradually can lose accuracy gradually 326without vanishing, so comparing them with zero 327(as one might on a VAX) will not reveal the loss. 328Fortunately, if a gradually underflowed value is 329destined to be added to something bigger than the 330underflow threshold, as is almost always the case, 331digits lost to gradual underflow will not be missed 332because they would have been rounded off anyway. 333So gradual underflows are usually 334.Em provably 335ignorable. 336The same cannot be said of underflows flushed to 0. 337.El 338.Pp 339At the option of an implementor conforming to IEEE 754, 340other ways to cope with exceptions may be provided: 341.Bl -enum 342.It 343ABORT. 344This mechanism classifies an exception in 345advance as an incident to be handled by means 346traditionally associated with error-handling 347statements like "ON ERROR GO TO ...". 348Different 349languages offer different forms of this statement, 350but most share the following characteristics: 351.Bl -dash 352.It 353No means is provided to substitute a value for 354the offending operation's result and resume 355computation from what may be the middle of an 356expression. 357An exceptional result is abandoned. 358.It 359In a subprogram that lacks an error-handling 360statement, an exception causes the subprogram to 361abort within whatever program called it, and so 362on back up the chain of calling subprograms until 363an error-handling statement is encountered or the 364whole task is aborted and memory is dumped. 365.El 366.It 367STOP. 368This mechanism, requiring an interactive 369debugging environment, is more for the programmer 370than the program. 371It classifies an exception in 372advance as a symptom of a programmer's error; the 373exception suspends execution as near as it can to 374the offending operation so that the programmer can 375look around to see how it happened. 376Quite often 377the first several exceptions turn out to be quite 378unexceptionable, so the programmer ought ideally 379to be able to resume execution after each one as if 380execution had not been stopped. 381.It 382\&... Other ways lie beyond the scope of this document. 383.El 384.Ed 385.Pp 386Ideally, each 387elementary function should act as if it were indivisible, or 388atomic, in the sense that ... 389.Bl -enum 390.It 391No exception should be signaled that is not deserved by 392the data supplied to that function. 393.It 394Any exception signaled should be identified with that 395function rather than with one of its subroutines. 396.It 397The internal behavior of an atomic function should not 398be disrupted when a calling program changes from 399one to another of the five or so ways of handling 400exceptions listed above, although the definition 401of the function may be correlated intentionally 402with exception handling. 403.El 404.Pp 405The functions in 406.Nm libm 407are only approximately atomic. 408They signal no inappropriate exception except possibly ... 409.Bl -tag -width indent -offset indent -compact 410.It Xo 411Over/Underflow 412.Xc 413when a result, if properly computed, might have lain barely within range, and 414.It Xo 415Inexact in 416.Fn cabs , 417.Fn cbrt , 418.Fn hypot , 419.Fn log10 420and 421.Fn pow 422.Xc 423when it happens to be exact, thanks to fortuitous cancellation of errors. 424.El 425Otherwise, ... 426.Bl -tag -width indent -offset indent -compact 427.It Xo 428Invalid Operation is signaled only when 429.Xc 430any result but \*(Na would probably be misleading. 431.It Xo 432Overflow is signaled only when 433.Xc 434the exact result would be finite but beyond the overflow threshold. 435.It Xo 436Divide-by-Zero is signaled only when 437.Xc 438a function takes exactly infinite values at finite operands. 439.It Xo 440Underflow is signaled only when 441.Xc 442the exact result would be nonzero but tinier than the underflow threshold. 443.It Xo 444Inexact is signaled only when 445.Xc 446greater range or precision would be needed to represent the exact result. 447.El 448.Sh BUGS 449Several functions required by 450.St -isoC-99 451are missing, and many functions are not available in their 452.Vt "long double" 453variants. 454.Pp 455On some architectures, trigonometric argument reduction is not 456performed accurately, resulting in errors greater than 1 457.Em ulp 458for large arguments to 459.Fn cos , 460.Fn sin , 461and 462.Fn tan . 463.Sh SEE ALSO 464.Xr fenv 3 , 465.Xr ieee 3 466.Pp 467An explanation of IEEE 754 and its proposed extension p854 468was published in the IEEE magazine MICRO in August 1984 under 469the title "A Proposed Radix- and Word-length-independent 470Standard for Floating-point Arithmetic" by 471.An "W. J. Cody" 472et al. 473The manuals for Pascal, C and BASIC on the Apple Macintosh 474document the features of IEEE 754 pretty well. 475Articles in the IEEE magazine COMPUTER vol.\& 14 no.\& 3 (Mar.\& 4761981), and in the ACM SIGNUM Newsletter Special Issue of 477Oct.\& 1979, may be helpful although they pertain to 478superseded drafts of the standard. 479.Sh HISTORY 480A math library with many of the present functions appeared in 481.At v7 . 482The library was substantially rewritten for 483.Bx 4.3 484to provide 485better accuracy and speed on machines supporting either VAX 486or IEEE 754 floating-point. 487Most of this library was replaced with FDLIBM, developed at Sun 488Microsystems, in 489.Fx 1.1.5 . 490