In computer science, a math library (or maths library) is a component of a programming language's standard library containing functions (or subroutines) for the most common mathematical functions, such as trigonometry and exponentiation. Bit-twiddling and control functionalities related to floating point numbers may also be included (such as in C).
Examples include:
- the C standard library math functions,[1]
- Java maths library[2]
- 'Prelude.Math' in haskell.[3]
In some languages (such as haskell) parts of the standard library (including maths) are imported by default.[4]
More advanced functionality such as linear algebra is usually provided in 3rd party libraries, such as a linear algebra library or vector maths library.
Implementation outline
editBasic operations
editIn a math library, it is frequently useful to use type punning to interpret a floating-point number as an unsigned integer of the same size. This allows for faster inspection of certain numeral properties (positive or not) and number comparison. In more advanced cases, bit twiddling may be used to modify a number in a specific way.
For more exact operation, a double double or even triple double format may be used. In this case, a high-precision number is expressed as the sum of two or three floating-point numbers.[5]
Transcendental functions
editTranscendental functions such as log, exponential, and trig functions make up the backbone of any math library. These functions are generally implemented by a polynomial fit, usually a Taylor polynomial or a Chebyshev polynomial derived by the Remez algorithm (having the benefit of an improved error bound), but the pre-processing steps are equally important.[5]
Trignometry
editRange reduction (also argument reduction, domain-splitting) is the first step for any function, after checks for unusual values (infinity and NaN) are performed. The goal here is to reduce the domain of the argument for the polynomial to process, using the function's symmetry and periodicity (if any), setting flags to indicate e.g. whether to negate the result in the end (if needed). It is worth noting that periodic functions require higher-than-input precision when reducing, with the prototypical method being the Payne–Hanek–Corbett algorithm.[6] After range reduction, near-zero values may be subject to a "fast path": for example, a tiny input x can be returned directly in sin, while 1 ? |x| may be used for cos.
The next step is the evaluation of the polynomial, with a conventional Horner's method used. After that, the sign of the result can be flipped according to the information from the range-reduction routine before being returned.
Logarithm and exponential
editLogarithm in base 2 is relatively straightforward, as the integer part k is already in the floating-point exponent; a preliminary range reduction is accordingly performed, yielding k. The mantissa x (where log2(x) is between -1/2 and 1/2) is then compared to a table and intervals for further reduction into a z with known log2 and an in-range x/z, along with polynomial coefficients used for the interval x/z is in. The result is then log(z) + log(x/z) + k.[7] Logarithm in other bases follow a similar approach, with the differences being (a) a different table is used; (b) k needs to be multiplied by log2(new-base).[7]
Exponential in base 2 is similarly the "base case" due to floating-point structure. The procedure is simply a combination of range reduction (usually by lookup) and a polynomial over the remaining mantissa.[7] The natural exponent may be implemented with a separate table for precision, while exp10 can simply be exp(x × log2(10)) when in-range.[7] Finally, the any-base exponent function pow() is built around exp() and log() in the general case.[7]
See also
editReferences
edit- ^ "C math library".
- ^ "java maths library".
- ^ "haskell prelude maths".
- ^ "haskell prelude".
- ^ a b Daramy-Loirat, Catherine; Defour, David; Dinechin, Florent de; Gallet, Matthieu; Gast, Nicolas; Lauter, Christoph; Muller, Jean-Michel (December 2006). CR-LIBM A library of correctly rounded elementary functions in double-precision (Report).
- ^ Brisebarre, N.; Defour, D.; Kornerup, P.; Muller, J.-M.; Revol, N. (March 2005). "A new range-reduction algorithm". IEEE Transactions on Computers. 54 (3): 331–339. doi:10.1109/TC.2005.36.
- ^ a b c d e musl v1.2.2 math directory, log1p.c, log2.c, log10.c, exp2.c