ln.js (4547B)
1 /** 2 * @license Apache-2.0 3 * 4 * Copyright (c) 2018 The Stdlib Authors. 5 * 6 * Licensed under the Apache License, Version 2.0 (the "License"); 7 * you may not use this file except in compliance with the License. 8 * You may obtain a copy of the License at 9 * 10 * http://www.apache.org/licenses/LICENSE-2.0 11 * 12 * Unless required by applicable law or agreed to in writing, software 13 * distributed under the License is distributed on an "AS IS" BASIS, 14 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 15 * See the License for the specific language governing permissions and 16 * limitations under the License. 17 * 18 * 19 * ## Notice 20 * 21 * The following copyright and license were part of the original implementation available as part of [FreeBSD]{@link https://svnweb.freebsd.org/base/release/9.3.0/lib/msun/src/e_log.c}. The implementation follows the original, but has been modified for JavaScript. 22 * 23 * ```text 24 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 25 * 26 * Developed at SunPro, a Sun Microsystems, Inc. business. 27 * Permission to use, copy, modify, and distribute this 28 * software is freely granted, provided that this notice 29 * is preserved. 30 * ``` 31 */ 32 33 'use strict'; 34 35 // MODULES // 36 37 var getHighWord = require( '@stdlib/number/float64/base/get-high-word' ); 38 var setHighWord = require( '@stdlib/number/float64/base/set-high-word' ); 39 var isnan = require( './../../../../base/assert/is-nan' ); 40 var BIAS = require( '@stdlib/constants/float64/exponent-bias' ); 41 var NINF = require( '@stdlib/constants/float64/ninf' ); 42 var polyvalP = require( './polyval_p.js' ); 43 var polyvalQ = require( './polyval_q.js' ); 44 45 46 // VARIABLES // 47 48 var LN2_HI = 6.93147180369123816490e-01; // 3FE62E42 FEE00000 49 var LN2_LO = 1.90821492927058770002e-10; // 3DEA39EF 35793C76 50 var TWO54 = 1.80143985094819840000e+16; // 0x43500000, 0x00000000 51 var ONE_THIRD = 0.33333333333333333; 52 53 // 0x000fffff = 1048575 => 0 00000000000 11111111111111111111 54 var HIGH_SIGNIFICAND_MASK = 0x000fffff|0; // asm type annotation 55 56 // 0x7ff00000 = 2146435072 => 0 11111111111 00000000000000000000 => biased exponent: 2047 = 1023+1023 => 2^1023 57 var HIGH_MAX_NORMAL_EXP = 0x7ff00000|0; // asm type annotation 58 59 // 0x00100000 = 1048576 => 0 00000000001 00000000000000000000 => biased exponent: 1 = -1022+1023 => 2^-1022 60 var HIGH_MIN_NORMAL_EXP = 0x00100000|0; // asm type annotation 61 62 // 0x3ff00000 = 1072693248 => 0 01111111111 00000000000000000000 => biased exponent: 1023 = 0+1023 => 2^0 = 1 63 var HIGH_BIASED_EXP_0 = 0x3ff00000|0; // asm type annotation 64 65 66 // MAIN // 67 68 /** 69 * Evaluates the natural logarithm. 70 * 71 * @param {NonNegativeNumber} x - input value 72 * @returns {number} function value 73 * 74 * @example 75 * var v = ln( 4.0 ); 76 * // returns ~1.386 77 * 78 * @example 79 * var v = ln( 0.0 ); 80 * // returns -Infinity 81 * 82 * @example 83 * var v = ln( Infinity ); 84 * // returns Infinity 85 * 86 * @example 87 * var v = ln( NaN ); 88 * // returns NaN 89 * 90 * @example 91 * var v = ln( -4.0 ); 92 * // returns NaN 93 */ 94 function ln( x ) { 95 var hfsq; 96 var hx; 97 var t2; 98 var t1; 99 var k; 100 var R; 101 var f; 102 var i; 103 var j; 104 var s; 105 var w; 106 var z; 107 108 if ( x === 0.0 ) { 109 return NINF; 110 } 111 if ( isnan( x ) || x < 0.0 ) { 112 return NaN; 113 } 114 hx = getHighWord( x ); 115 k = 0|0; // asm type annotation 116 if ( hx < HIGH_MIN_NORMAL_EXP ) { 117 // Case: 0 < x < 2**-1022 118 k -= 54|0; // asm type annotation 119 120 // Subnormal number, scale up `x`: 121 x *= TWO54; 122 hx = getHighWord( x ); 123 } 124 if ( hx >= HIGH_MAX_NORMAL_EXP ) { 125 return x + x; 126 } 127 k += ( ( hx>>20 ) - BIAS )|0; // asm type annotation 128 hx &= HIGH_SIGNIFICAND_MASK; 129 i = ( (hx+0x95f64) & 0x100000 )|0; // asm type annotation 130 131 // Normalize `x` or `x/2`... 132 x = setHighWord( x, hx|(i^HIGH_BIASED_EXP_0) ); 133 k += ( i>>20 )|0; // asm type annotation 134 f = x - 1.0; 135 if ( (HIGH_SIGNIFICAND_MASK&(2+hx)) < 3 ) { 136 // Case: -2**-20 <= f < 2**-20 137 if ( f === 0.0 ) { 138 if ( k === 0 ) { 139 return 0.0; 140 } 141 return (k * LN2_HI) + (k * LN2_LO); 142 } 143 R = f * f * ( 0.5 - (ONE_THIRD*f) ); 144 if ( k === 0 ) { 145 return f - R; 146 } 147 return (k * LN2_HI) - ( (R-(k*LN2_LO)) - f ); 148 } 149 s = f / (2.0 + f); 150 z = s * s; 151 i = ( hx - 0x6147a )|0; // asm type annotation 152 w = z * z; 153 j = ( 0x6b851 - hx )|0; // asm type annotation 154 t1 = w * polyvalP( w ); 155 t2 = z * polyvalQ( w ); 156 i |= j; 157 R = t2 + t1; 158 if ( i > 0 ) { 159 hfsq = 0.5 * f * f; 160 if ( k === 0 ) { 161 return f - ( hfsq - (s * (hfsq+R)) ); 162 } 163 return (k * LN2_HI) - ( hfsq - ((s*(hfsq+R))+(k*LN2_LO)) - f ); 164 } 165 if ( k === 0 ) { 166 return f - (s*(f-R)); 167 } 168 return (k * LN2_HI) - ( ( (s*(f-R)) - (k*LN2_LO) ) - f ); 169 } 170 171 172 // EXPORTS // 173 174 module.exports = ln;