time-to-botec

Benchmark sampling in different programming languages
Log | Files | Refs | README

j0.js (3066B)


      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 original C++ code and copyright notice are from the [Boost library]{@link https://github.com/boostorg/math/blob/develop/include/boost/math/special_functions/detail/bessel_j0.hpp}. The implementation has been modified for JavaScript.
     22 *
     23 * ```text
     24 * Copyright Xiaogang Zhang, 2006.
     25 *
     26 * Use, modification and distribution are subject to the
     27 * Boost Software License, Version 1.0. (See accompanying file
     28 * LICENSE or copy at http://www.boost.org/LICENSE_1_0.txt)
     29 * ```
     30 */
     31 
     32 'use strict';
     33 
     34 // MODULES //
     35 
     36 var sqrt = require( './../../../../base/special/sqrt' );
     37 var sincos = require( './../../../../base/special/sincos' );
     38 var PINF = require( '@stdlib/constants/float64/pinf' );
     39 var poly1 = require( './rational_p1q1.js' );
     40 var poly2 = require( './rational_p2q2.js' );
     41 var polyC = require( './rational_pcqc.js' );
     42 var polyS = require( './rational_psqs.js' );
     43 
     44 
     45 // VARIABLES //
     46 
     47 var ONE_DIV_SQRT_PI = 0.5641895835477563;
     48 var x1 = 2.4048255576957727686e+00;
     49 var x2 = 5.5200781102863106496e+00;
     50 var x11 = 6.160e+02;
     51 var x12 = -1.42444230422723137837e-03;
     52 var x21 = 1.4130e+03;
     53 var x22 = 5.46860286310649596604e-04;
     54 
     55 // `sincos` workspace:
     56 var sc = [ 0.0, 0.0 ]; // WARNING: not thread safe
     57 
     58 
     59 // MAIN //
     60 
     61 /**
     62 * Computes the Bessel function of the first kind of order zero.
     63 *
     64 * @param {number} x - input value
     65 * @returns {number} evaluated Bessel function
     66 *
     67 * @example
     68 * var v = j0( 0.0 );
     69 * // returns 1.0
     70 *
     71 * v = j0( 1.0 );
     72 * // returns ~0.765
     73 *
     74 * v = j0( Infinity );
     75 * // returns 0.0
     76 *
     77 * v = j0( -Infinity );
     78 * // returns 0.0
     79 *
     80 * v = j0( NaN );
     81 * // returns NaN
     82 */
     83 function j0( x ) {
     84 	var rc;
     85 	var rs;
     86 	var y2;
     87 	var r;
     88 	var y;
     89 	var f;
     90 
     91 	if ( x < 0 ) {
     92 		x = -x;
     93 	}
     94 	if ( x === PINF ) {
     95 		return 0.0;
     96 	}
     97 	if ( x === 0 ) {
     98 		return 1.0;
     99 	}
    100 	if ( x <= 4.0 ) {
    101 		y = x * x;
    102 		r = poly1( y );
    103 		f = ( x+x1 ) * ( (x - (x11/256.0)) - x12 );
    104 		return f * r;
    105 	}
    106 	if ( x <= 8.0 ) {
    107 		y = 1.0 - ( ( x*x )/64.0 );
    108 		r = poly2( y );
    109 		f = ( x+x2 ) * ( (x - (x21/256.0)) - x22 );
    110 		return f * r;
    111 	}
    112 	y = 8.0 / x;
    113 	y2 = y * y;
    114 	rc = polyC( y2 );
    115 	rs = polyS( y2 );
    116 	f = ONE_DIV_SQRT_PI / sqrt(x);
    117 
    118 	/*
    119 	* What follows is really just:
    120 	*
    121 	* ```
    122 	* var z = x - pi/4;
    123 	* return f * (rc * cos(z) - y * rs * sin(z));
    124 	* ```
    125 	*
    126 	* But using the addition formulae for sin and cos, plus the special values for sin/cos of `π/4`.
    127 	*/
    128 	sincos( sc, x );
    129 	return f * ( ( rc * (sc[1]+sc[0]) ) - ( (y*rs) * (sc[0]-sc[1]) ) );
    130 }
    131 
    132 
    133 // EXPORTS //
    134 
    135 module.exports = j0;