version bump 0.3.0: cleanup
- eliminated array reduce (performance improvement) - explicit treatment of divergent values - better handling of non-integral order (fixes #3, h/t @vron) - reorganized source, removed voc dependency (fixes #2, h/t @hmalphettes) - new browser demo - more detailed test suite and coverage tests - updated travis versions for test - miscellaneous adjustments to tooling
This commit is contained in:
parent
31237637de
commit
73ee33e402
|
|
@ -0,0 +1,19 @@
|
|||
[ignore]
|
||||
.*/node_modules/.*
|
||||
.*/dist/.*
|
||||
.*/test.js
|
||||
|
||||
.*/bits/.*
|
||||
.*/ctest/.*
|
||||
.*/misc/.*
|
||||
.*/perf/.*
|
||||
|
||||
.*/bessel.js
|
||||
|
||||
[include]
|
||||
bessel.flow.js
|
||||
|
||||
[libs]
|
||||
bits/10_types.js
|
||||
|
||||
[options]
|
||||
|
|
@ -1 +1,2 @@
|
|||
node_modules
|
||||
misc/coverage.html
|
||||
|
|
|
|||
|
|
@ -0,0 +1,6 @@
|
|||
{
|
||||
"requireCommaBeforeLineBreak": true,
|
||||
"disallowTrailingWhitespace": true,
|
||||
"disallowTrailingComma": true
|
||||
}
|
||||
|
||||
12
.travis.yml
12
.travis.yml
|
|
@ -1,7 +1,15 @@
|
|||
language: node_js
|
||||
node_js:
|
||||
- "0.11"
|
||||
- "6"
|
||||
- "5"
|
||||
- "4.2"
|
||||
- "0.12"
|
||||
- "0.10"
|
||||
- "0.8"
|
||||
before_install:
|
||||
- "npm install -g mocha"
|
||||
- "npm install -g npm@next"
|
||||
- "npm install -g mocha voc"
|
||||
- "npm install blanket"
|
||||
- "npm install coveralls mocha-lcov-reporter"
|
||||
after_success:
|
||||
- "make coveralls-spin"
|
||||
|
|
|
|||
69
Makefile
69
Makefile
|
|
@ -1,7 +1,66 @@
|
|||
LIBRARY=bessel
|
||||
LIB=bessel
|
||||
REQS=
|
||||
ADDONS=
|
||||
AUXTARGETS=
|
||||
HTMLLINT=index.html
|
||||
|
||||
$(LIBRARY).js: bessel.md
|
||||
node_modules/.bin/voc $^ > $@
|
||||
ULIB=$(shell echo $(LIB) | tr a-z A-Z)
|
||||
DEPS=$(sort $(wildcard bits/*.js))
|
||||
TARGET=$(LIB).js
|
||||
FLOWTARGET=$(LIB).flow.js
|
||||
|
||||
test mocha:
|
||||
mocha -R spec
|
||||
## Main Targets
|
||||
|
||||
.PHONY: all
|
||||
all: $(TARGET) $(AUXTARGETS) ## Build library and auxiliary scripts
|
||||
|
||||
$(TARGET) $(AUXTARGETS): %.js : %.flow.js
|
||||
node -e 'process.stdout.write(require("fs").readFileSync("$<","utf8").replace(/^[ \t]*\/\*[:#][^*]*\*\/\s*(\n)?/gm,"").replace(/\/\*[:#][^*]*\*\//gm,""))' > $@
|
||||
|
||||
$(FLOWTARGET): $(DEPS)
|
||||
cat $^ | tr -d '\15\32' > $@
|
||||
|
||||
bits/01_version.js: package.json
|
||||
echo "BESSEL.version = '"`grep version package.json | awk '{gsub(/[^0-9a-z\.-]/,"",$$2); print $$2}'`"';" > $@
|
||||
|
||||
.PHONY: clean
|
||||
clean: ## Remove targets and build artifacts
|
||||
rm -f $(TARGET) $(FLOWTARGET)
|
||||
|
||||
## Testing
|
||||
|
||||
.PHONY: test mocha
|
||||
test mocha: test.js $(TARGET) ## Run test suite
|
||||
mocha -R spec -t 20000
|
||||
|
||||
## Code Checking
|
||||
|
||||
.PHONY: lint
|
||||
lint: $(TARGET) $(AUXTARGETS) ## Run jshint and jscs checks
|
||||
@jshint --show-non-errors $(TARGET) $(AUXTARGETS)
|
||||
@jshint --show-non-errors package.json
|
||||
@jshint --show-non-errors --extract=always $(HTMLLINT)
|
||||
@jscs $(TARGET) $(AUXTARGETS)
|
||||
|
||||
.PHONY: flow
|
||||
flow: lint ## Run flow checker
|
||||
@flow check --all --show-all-errors
|
||||
|
||||
.PHONY: cov
|
||||
cov: misc/coverage.html ## Run coverage test
|
||||
|
||||
misc/coverage.html: $(TARGET) test.js
|
||||
mocha --require blanket -R html-cov -t 20000 > $@
|
||||
|
||||
.PHONY: coveralls
|
||||
coveralls: ## Coverage Test + Send to coveralls.io
|
||||
mocha --require blanket --reporter mocha-lcov-reporter -t 20000 | node ./node_modules/coveralls/bin/coveralls.js
|
||||
|
||||
|
||||
.PHONY: help
|
||||
help:
|
||||
@grep -hE '(^[a-zA-Z_-][ a-zA-Z_-]*:.*?|^#[#*])' $(MAKEFILE_LIST) | bash misc/help.sh
|
||||
|
||||
#* To show a spinner, append "-spin" to any target e.g. cov-spin
|
||||
%-spin:
|
||||
@make $* & bash misc/spin.sh $$!
|
||||
|
|
|
|||
94
README.md
94
README.md
|
|
@ -1,36 +1,92 @@
|
|||
# Bessel.JS
|
||||
# bessel
|
||||
|
||||
Pure-JS implementation of the Bessel functions (J,Y,I,K), for node and browser
|
||||
Pure-JS implementation of Bessel functions J,Y,I,K (for the browser and nodejs).
|
||||
Emphasis on correctness and performance for integer order.
|
||||
|
||||
The standard notation is used here:
|
||||
|
||||
- J is the Bessel function of the first kind
|
||||
- Y is the Bessel function of the second kind
|
||||
- I is the modified Bessel function of the first kind
|
||||
- K is the modified Bessel function of the first kind
|
||||
- `J` is the Bessel function of the first kind
|
||||
- `Y` is the Bessel function of the second kind
|
||||
- `I` is the modified Bessel function of the first kind
|
||||
- `K` is the modified Bessel function of the second kind
|
||||
|
||||
# Usage
|
||||
## Installation
|
||||
|
||||
The functions `besselj`, `bessely`, `besseli`, `besselk` are exposed when you include
|
||||
the script `bessel.js`:
|
||||
With [npm](https://www.npmjs.org/package/bessel):
|
||||
|
||||
```html>
|
||||
```bash
|
||||
$ npm install bessel
|
||||
```
|
||||
|
||||
In the browser:
|
||||
|
||||
```html
|
||||
<script src="bessel.js"></script>
|
||||
<script>console.log(besselj(1,2));</script>
|
||||
```
|
||||
|
||||
See `test.html` for an example
|
||||
The browser exposes a variable `BESSEL`
|
||||
|
||||
In node, those four functions are exported:
|
||||
The script will manipulate `module.exports` if available (e.g. in a CommonJS
|
||||
`require` context). This is not always desirable. To prevent the behavior,
|
||||
define `DO_NOT_EXPORT_BESSEL`
|
||||
|
||||
```js>
|
||||
var besselj01 = require('bessel').besselj(0,1);
|
||||
## Usage
|
||||
|
||||
In all cases, the relevant function takes two arguments (`value`, `order`).
|
||||
|
||||
The return value is a JS number. `NaN` signals an error in calculation.
|
||||
|
||||
- `BESSEL.besselj(x, n)` computes `J_{n}(x)`
|
||||
|
||||
- `BESSEL.bessely(x, n)` computes `Y_{n}(x)`
|
||||
|
||||
- `BESSEL.besseli(x, n)` computes `I_{n}(x)`
|
||||
|
||||
- `BESSEL.besselk(x, n)` computes `K_{n}(x)`
|
||||
|
||||
For example:
|
||||
|
||||
```js
|
||||
> // var BESSEL = require('bessel'); // uncomment this line if in node
|
||||
> BESSEL.besselj(1.5,0) // 0.5118276712499389
|
||||
> BESSEL.bessely(1.5,0) // 0.38244892476502895
|
||||
> BESSEL.besseli(1.5,0) // 1.6467232021476754
|
||||
> BESSEL.besselk(1.5,0) // 0.2138055693236539
|
||||
|
||||
> var Y = BESSEL.bessely
|
||||
> Y(Math.PI, 5) + Y(Math.PI, 3) - (2 * 4 / Math.PI) * Y(Math.PI, 4) // 0
|
||||
```
|
||||
|
||||
Each function follows Excel semantics `(value, function-order)`. For example,
|
||||
## Testing
|
||||
|
||||
```js>
|
||||
bessel.besselj(1.5, 1)
|
||||
`make test` will run the nodejs-based test.
|
||||
|
||||
To generate the `excel.tsv` test cases, make a 6-column Excel Sheet:
|
||||
|
||||
| x | n | `BESSELI` | `BESSELJ` | `BESSELK` | `BESSELY` |
|
||||
|---|---|:--------------:|:--------------:|:--------------:|:--------------:|
|
||||
| x | n |`BESSELI(A1,B1)`|`BESSELJ(A1,B1)`|`BESSELK(A1,B1)`|`BESSELY(A1,B1)`|
|
||||
|
||||
To generate the `mma.tsv` test cases, use the Mathematica Bessel Functions:
|
||||
|
||||
```mathematica
|
||||
(* Bessel_[value, order] *)
|
||||
F[x_,n_]:={x/2,n,BesselI[n,x/2], BesselJ[n,x/2], BesselK[n,x/2], BesselY[n,x/2]}
|
||||
Do[Print[ExportString[N[F[x,n],10],"csv"]], {n,1,3}, {x,1,20} ]
|
||||
```
|
||||
|
||||
is the value of the bessel function J1 at the point x=1.5
|
||||
Note: Each function follows Excel semantics `(value, order)`. Other platforms
|
||||
like Mathematica reverse the argument order.
|
||||
|
||||
## License
|
||||
|
||||
Please consult the attached LICENSE file for details. All rights not explicitly
|
||||
granted by the MIT License are reserved by the Original Author.
|
||||
|
||||
## Badges
|
||||
|
||||
[](https://travis-ci.org/SheetJS/bessel)
|
||||
|
||||
[](https://coveralls.io/r/SheetJS/bessel?branch=master)
|
||||
|
||||
[](https://github.com/SheetJS/bessel)
|
||||
|
|
|
|||
|
|
@ -0,0 +1,252 @@
|
|||
/* bessel.js (C) 2013-present SheetJS -- http://sheetjs.com */
|
||||
/* vim: set ts=2: */
|
||||
/*exported BESSEL */
|
||||
var BESSEL;
|
||||
/*:: declare var DO_NOT_EXPORT_BESSEL: any; */
|
||||
/*:: declare var define: any; */
|
||||
(function (factory) {
|
||||
/*jshint ignore:start */
|
||||
if(typeof DO_NOT_EXPORT_BESSEL === 'undefined') {
|
||||
if('object' === typeof exports) {
|
||||
factory(exports);
|
||||
} else if ('function' === typeof define && define.amd) {
|
||||
define(function () {
|
||||
var module = {};
|
||||
factory(module);
|
||||
return module;
|
||||
});
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
/*jshint ignore:end */
|
||||
}(function(BESSEL) {
|
||||
BESSEL.version = '0.3.0';
|
||||
var M = Math;
|
||||
|
||||
/*::
|
||||
type BesselN = (x:number) => number;
|
||||
type BesselF = (x:number, n:number) => number;
|
||||
*/
|
||||
function _horner(arr/*:Array<number>*/, v/*:number*/)/*:number*/ { for(var i = 0, z = 0; i < arr.length; ++i) z = v * z + arr[i]; return z; }
|
||||
function _bessel_iter(x/*:number*/, n/*:number*/, f0/*:number*/, f1/*:number*/, sign/*:number*/)/*:number*/ {
|
||||
if(n === 0) return f0;
|
||||
if(n === 1) return f1;
|
||||
var tdx = 2 / x, f2 = f1;
|
||||
for(var o = 1; o < n; ++o) {
|
||||
f2 = f1 * o * tdx + sign * f0;
|
||||
f0 = f1; f1 = f2;
|
||||
}
|
||||
return f2;
|
||||
}
|
||||
function _bessel_wrap(bessel0/*:BesselN*/, bessel1/*:BesselN*/, name/*:string*/, nonzero/*:number*/, sign/*:number*/)/*:BesselF*/ {
|
||||
return function bessel(x/*:number*/,n/*:number*/) {
|
||||
if(nonzero) {
|
||||
if(x === 0) return (nonzero == 1 ? -Infinity : Infinity);
|
||||
else if(x < 0) return NaN;
|
||||
}
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) return NaN;
|
||||
n|=0;
|
||||
var b0/*:number*/ = bessel0(x), b1/*:number*/ = bessel1(x);
|
||||
return _bessel_iter(x, n, b0, b1, sign);
|
||||
};
|
||||
}
|
||||
var besselj/*:BesselF*/ = (function() {
|
||||
var W = 0.636619772; // 2 / Math.PI
|
||||
|
||||
var b0_a1a = [57568490574.0, -13362590354.0, 651619640.7, -11214424.18, 77392.33017, -184.9052456].reverse();
|
||||
var b0_a2a = [57568490411.0, 1029532985.0, 9494680.718, 59272.64853, 267.8532712, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934935152e-7].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x * x;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
a = a1 / a2;
|
||||
} else {
|
||||
var xx = x - 0.785398164;
|
||||
y = 64 / y;
|
||||
a1 = _horner(b0_a1b, y);
|
||||
a2 = _horner(b0_a2b, y);
|
||||
a = M.sqrt(W/x)*(M.cos(xx)*a1-M.sin(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
var b1_a1a = [72362614232.0, -7895059235.0, 242396853.1, -2972611.439, 15704.48260, -30.16036606].reverse();
|
||||
var b1_a2a = [144725228442.0, 2300535178.0, 18583304.74, 99447.43394, 376.9991397, 1.0].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = M.abs(x) - 2.356194491;
|
||||
if(Math.abs(x)< 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
a = a1 / a2;
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1=_horner(b1_a1b, y);
|
||||
a2=_horner(b1_a2b, y);
|
||||
a=M.sqrt(W/M.abs(x))*(M.cos(xx)*a1-M.sin(xx)*a2*8/M.abs(x));
|
||||
if(x < 0) a = -a;
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
return function besselj(x/*:number*/, n/*:number*/)/*:number*/ {
|
||||
n = Math.round(n);
|
||||
if(!isFinite(x)) return isNaN(x) ? x : 0;
|
||||
if(n < 0) return ((n%2)?-1:1)*besselj(x, -n);
|
||||
if(x < 0) return ((n%2)?-1:1)*besselj(-x, n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(x === 0) return 0;
|
||||
|
||||
var ret=0.0;
|
||||
if(x > n) {
|
||||
ret = _bessel_iter(x, n, bessel0(x), bessel1(x),-1);
|
||||
} else {
|
||||
var m=2*M.floor((n+M.floor(M.sqrt(40*n)))/2);
|
||||
var jsum=false;
|
||||
var bjp=0.0, sum=0.0;
|
||||
var bj=1.0, bjm = 0.0;
|
||||
var tox = 2 / x;
|
||||
for (var j=m;j>0;j--) {
|
||||
bjm=j*tox*bj-bjp;
|
||||
bjp=bj;
|
||||
bj=bjm;
|
||||
if (M.abs(bj) > 1E10) {
|
||||
bj *= 1E-10;
|
||||
bjp *= 1E-10;
|
||||
ret *= 1E-10;
|
||||
sum *= 1E-10;
|
||||
}
|
||||
if (jsum) sum += bj;
|
||||
jsum=!jsum;
|
||||
if (j == n) ret=bjp;
|
||||
}
|
||||
sum=2.0*sum-bj;
|
||||
ret /= sum;
|
||||
}
|
||||
return ret;
|
||||
};
|
||||
})();
|
||||
var bessely/*:BesselF*/ = (function() {
|
||||
var W = 0.636619772;
|
||||
|
||||
var b0_a1a = [-2957821389.0, 7062834065.0, -512359803.6, 10879881.29, -86327.92757, 228.4622733].reverse();
|
||||
var b0_a2a = [40076544269.0, 745249964.8, 7189466.438, 47447.26470, 226.1030244, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934945152e-7].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x * x, xx = x - 0.785398164;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
a = a1/a2 + W * besselj(x,0) * M.log(x);
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1 = _horner(b0_a1b, y);
|
||||
a2 = _horner(b0_a2b, y);
|
||||
a = M.sqrt(W/x)*(M.sin(xx)*a1+M.cos(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
var b1_a1a = [-0.4900604943e13, 0.1275274390e13, -0.5153438139e11, 0.7349264551e9, -0.4237922726e7, 0.8511937935e4].reverse();
|
||||
var b1_a2a = [0.2499580570e14, 0.4244419664e12, 0.3733650367e10, 0.2245904002e8, 0.1020426050e6, 0.3549632885e3, 1].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = x - 2.356194491;
|
||||
if(x < 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
a = a1/a2 + W * (besselj(x,1) * M.log(x) - 1 / x);
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1=_horner(b1_a1b, y);
|
||||
a2=_horner(b1_a2b, y);
|
||||
a=M.sqrt(W/x)*(M.sin(xx)*a1+M.cos(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
return _bessel_wrap(bessel0, bessel1, 'BESSELY', 1, -1);
|
||||
})();
|
||||
var besseli/*:BesselF*/ = (function() {
|
||||
var b0_a = [1.0, 3.5156229, 3.0899424, 1.2067492, 0.2659732, 0.360768e-1, 0.45813e-2].reverse();
|
||||
var b0_b = [0.39894228, 0.1328592e-1, 0.225319e-2, -0.157565e-2, 0.916281e-2, -0.2057706e-1, 0.2635537e-1, -0.1647633e-1, 0.392377e-2].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
if(x <= 3.75) return _horner(b0_a, x*x/(3.75*3.75));
|
||||
return M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b0_b, 3.75/M.abs(x));
|
||||
}
|
||||
|
||||
var b1_a = [0.5, 0.87890594, 0.51498869, 0.15084934, 0.2658733e-1, 0.301532e-2, 0.32411e-3].reverse();
|
||||
var b1_b = [0.39894228, -0.3988024e-1, -0.362018e-2, 0.163801e-2, -0.1031555e-1, 0.2282967e-1, -0.2895312e-1, 0.1787654e-1, -0.420059e-2].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
if(x < 3.75) return x * _horner(b1_a, x*x/(3.75*3.75));
|
||||
return (x < 0 ? -1 : 1) * M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b1_b, 3.75/M.abs(x));
|
||||
}
|
||||
|
||||
return function besseli(x/*:number*/, n/*:number*/)/*:number*/ {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) return NaN;
|
||||
if(M.abs(x) === 0) return 0;
|
||||
if(x == Infinity) return Infinity;
|
||||
|
||||
var ret = 0.0, j, tox = 2 / M.abs(x), bip = 0.0, bi=1.0, bim=0.0;
|
||||
var m=2*M.round((n+M.round(M.sqrt(40*n)))/2);
|
||||
for (j=m;j>0;j--) {
|
||||
bim=j*tox*bi + bip;
|
||||
bip=bi; bi=bim;
|
||||
if (M.abs(bi) > 1E10) {
|
||||
bi *= 1E-10;
|
||||
bip *= 1E-10;
|
||||
ret *= 1E-10;
|
||||
}
|
||||
if(j == n) ret = bip;
|
||||
}
|
||||
ret *= besseli(x, 0) / bi;
|
||||
return x < 0 && (n%2) ? -ret : ret;
|
||||
};
|
||||
|
||||
})();
|
||||
|
||||
var besselk/*:BesselF*/ = (function() {
|
||||
var b0_a = [-0.57721566, 0.42278420, 0.23069756, 0.3488590e-1, 0.262698e-2, 0.10750e-3, 0.74e-5].reverse();
|
||||
var b0_b = [1.25331414, -0.7832358e-1, 0.2189568e-1, -0.1062446e-1, 0.587872e-2, -0.251540e-2, 0.53208e-3].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
if(x <= 2) return -M.log(x/2) * besseli(x,0) + _horner(b0_a, x*x/4);
|
||||
return M.exp(-x) / M.sqrt(x) * _horner(b0_b, 2/x);
|
||||
}
|
||||
|
||||
var b1_a = [1.0, 0.15443144, -0.67278579, -0.18156897, -0.1919402e-1, -0.110404e-2, -0.4686e-4].reverse();
|
||||
var b1_b = [1.25331414, 0.23498619, -0.3655620e-1, 0.1504268e-1, -0.780353e-2, 0.325614e-2, -0.68245e-3].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
if(x <= 2) return M.log(x/2) * besseli(x,1) + (1/x) * _horner(b1_a, x*x/4);
|
||||
return M.exp(-x)/M.sqrt(x)*_horner(b1_b, 2/x);
|
||||
}
|
||||
|
||||
return _bessel_wrap(bessel0, bessel1, 'BESSELK', 2, 1);
|
||||
})();
|
||||
BESSEL.besselj = besselj;
|
||||
BESSEL.bessely = bessely;
|
||||
BESSEL.besseli = besseli;
|
||||
BESSEL.besselk = besselk;
|
||||
}));
|
||||
141
bessel.js
141
bessel.js
|
|
@ -1,55 +1,85 @@
|
|||
/* bessel.js (C) 2013-present SheetJS -- http://sheetjs.com */
|
||||
/* vim: set ts=2: */
|
||||
/*exported BESSEL */
|
||||
var BESSEL;
|
||||
(function (factory) {
|
||||
/*jshint ignore:start */
|
||||
if(typeof DO_NOT_EXPORT_BESSEL === 'undefined') {
|
||||
if('object' === typeof exports) {
|
||||
factory(exports);
|
||||
} else if ('function' === typeof define && define.amd) {
|
||||
define(function () {
|
||||
var module = {};
|
||||
factory(module);
|
||||
return module;
|
||||
});
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
/*jshint ignore:end */
|
||||
}(function(BESSEL) {
|
||||
BESSEL.version = '0.3.0';
|
||||
var M = Math;
|
||||
function _horner(arr, v) { return arr.reduce(function(z,w){return v * z + w;},0); };
|
||||
|
||||
function _horner(arr, v) { for(var i = 0, z = 0; i < arr.length; ++i) z = v * z + arr[i]; return z; }
|
||||
function _bessel_iter(x, n, f0, f1, sign) {
|
||||
if(!sign) sign = -1;
|
||||
var tdx = 2 / x, f2;
|
||||
if(n === 0) return f0;
|
||||
if(n === 1) return f1;
|
||||
for(var o = 1; o != n; ++o) {
|
||||
var tdx = 2 / x, f2 = f1;
|
||||
for(var o = 1; o < n; ++o) {
|
||||
f2 = f1 * o * tdx + sign * f0;
|
||||
f0 = f1; f1 = f2;
|
||||
}
|
||||
return f1;
|
||||
return f2;
|
||||
}
|
||||
function _bessel_wrap(bessel0, bessel1, name, nonzero, sign) {
|
||||
return function bessel(x,n) {
|
||||
if(nonzero) {
|
||||
if(x === 0) return (nonzero == 1 ? -Infinity : Infinity);
|
||||
else if(x < 0) return NaN;
|
||||
}
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) throw name + ': Order (' + n + ') must be nonnegative';
|
||||
if(nonzero == 1 && x === 0) throw name + ': Undefined when x == 0';
|
||||
if(nonzero == 2 && x <= 0) throw name + ': Undefined when x <= 0';
|
||||
if(n < 0) return NaN;
|
||||
n|=0;
|
||||
var b0 = bessel0(x), b1 = bessel1(x);
|
||||
return _bessel_iter(x, n, b0, b1, sign);
|
||||
};
|
||||
}
|
||||
var besselj = (function() {
|
||||
var b0_a1a = [57568490574.0,-13362590354.0,651619640.7,-11214424.18,77392.33017,-184.9052456].reverse();
|
||||
var b0_a2a = [57568490411.0,1029532985.0,9494680.718,59272.64853,267.8532712,1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934935152e-7].reverse();
|
||||
var W = 0.636619772; // 2 / Math.PI
|
||||
|
||||
var b0_a1a = [57568490574.0, -13362590354.0, 651619640.7, -11214424.18, 77392.33017, -184.9052456].reverse();
|
||||
var b0_a2a = [57568490411.0, 1029532985.0, 9494680.718, 59272.64853, 267.8532712, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934935152e-7].reverse();
|
||||
|
||||
function bessel0(x) {
|
||||
var a, a1, a2, y = x * x, xx = M.abs(x) - 0.785398164;
|
||||
if(M.abs(x) < 8) {
|
||||
var a=0, a1=0, a2=0, y = x * x;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
a = a1/a2;
|
||||
}
|
||||
else {
|
||||
a = a1 / a2;
|
||||
} else {
|
||||
var xx = x - 0.785398164;
|
||||
y = 64 / y;
|
||||
a1 = _horner(b0_a1b, y);
|
||||
a2 = _horner(b0_a2b, y);
|
||||
a = M.sqrt(W/M.abs(x))*(M.cos(xx)*a1-M.sin(xx)*a2*8/M.abs(x));
|
||||
a = M.sqrt(W/x)*(M.cos(xx)*a1-M.sin(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
var b1_a1a = [72362614232.0,-7895059235.0,242396853.1,-2972611.439, 15704.48260, -30.16036606].reverse();
|
||||
|
||||
var b1_a1a = [72362614232.0, -7895059235.0, 242396853.1, -2972611.439, 15704.48260, -30.16036606].reverse();
|
||||
var b1_a2a = [144725228442.0, 2300535178.0, 18583304.74, 99447.43394, 376.9991397, 1.0].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x) {
|
||||
var a, a1, a2, y = x*x, xx = M.abs(x) - 2.356194491;
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = M.abs(x) - 2.356194491;
|
||||
if(Math.abs(x)< 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
|
|
@ -63,22 +93,26 @@ var besselj = (function() {
|
|||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
return function besselj(x, n) {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(M.abs(x));
|
||||
if(n === 1) return bessel1(M.abs(x));
|
||||
if(n < 0) throw 'BESSELJ: Order (' + n + ') must be nonnegative';
|
||||
if(M.abs(x) === 0) return 0;
|
||||
if(!isFinite(x)) return isNaN(x) ? x : 0;
|
||||
if(n < 0) return ((n%2)?-1:1)*besselj(x, -n);
|
||||
if(x < 0) return ((n%2)?-1:1)*besselj(-x, n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(x === 0) return 0;
|
||||
|
||||
var ret, j, tox = 2 / M.abs(x), m, jsum, sum, bjp, bj, bjm;
|
||||
if(M.abs(x) > n) {
|
||||
ret = _bessel_iter(x, n, bessel0(M.abs(x)), bessel1(M.abs(x)),-1);
|
||||
var ret=0.0;
|
||||
if(x > n) {
|
||||
ret = _bessel_iter(x, n, bessel0(x), bessel1(x),-1);
|
||||
} else {
|
||||
m=2*M.floor((n+M.floor(M.sqrt(40*n)))/2);
|
||||
jsum=0;
|
||||
bjp=ret=sum=0.0;
|
||||
bj=1.0;
|
||||
for (j=m;j>0;j--) {
|
||||
var m=2*M.floor((n+M.floor(M.sqrt(40*n)))/2);
|
||||
var jsum=false;
|
||||
var bjp=0.0, sum=0.0;
|
||||
var bj=1.0, bjm = 0.0;
|
||||
var tox = 2 / x;
|
||||
for (var j=m;j>0;j--) {
|
||||
bjm=j*tox*bj-bjp;
|
||||
bjp=bj;
|
||||
bj=bjm;
|
||||
|
|
@ -95,18 +129,19 @@ var besselj = (function() {
|
|||
sum=2.0*sum-bj;
|
||||
ret /= sum;
|
||||
}
|
||||
return x < 0 && (n%2) ? -ret : ret;
|
||||
return ret;
|
||||
};
|
||||
})();
|
||||
var bessely = (function() {
|
||||
var W = 0.636619772;
|
||||
|
||||
var b0_a1a = [-2957821389.0, 7062834065.0, -512359803.6, 10879881.29, -86327.92757, 228.4622733].reverse();
|
||||
var b0_a2a = [40076544269.0, 745249964.8, 7189466.438, 47447.26470, 226.1030244, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934945152e-7].reverse();
|
||||
|
||||
var W = 0.636619772;
|
||||
function bessel0(x) {
|
||||
var a, a1, a2, y = x * x, xx = x - 0.785398164;
|
||||
var a=0, a1=0, a2=0, y = x * x, xx = x - 0.785398164;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
|
|
@ -124,8 +159,9 @@ var bessely = (function() {
|
|||
var b1_a2a = [0.2499580570e14, 0.4244419664e12, 0.3733650367e10, 0.2245904002e8, 0.1020426050e6, 0.3549632885e3, 1].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x) {
|
||||
var a, a1, a2, y = x*x, xx = x - 2.356194491;
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = x - 2.356194491;
|
||||
if(x < 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
|
|
@ -144,6 +180,7 @@ var bessely = (function() {
|
|||
var besseli = (function() {
|
||||
var b0_a = [1.0, 3.5156229, 3.0899424, 1.2067492, 0.2659732, 0.360768e-1, 0.45813e-2].reverse();
|
||||
var b0_b = [0.39894228, 0.1328592e-1, 0.225319e-2, -0.157565e-2, 0.916281e-2, -0.2057706e-1, 0.2635537e-1, -0.1647633e-1, 0.392377e-2].reverse();
|
||||
|
||||
function bessel0(x) {
|
||||
if(x <= 3.75) return _horner(b0_a, x*x/(3.75*3.75));
|
||||
return M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b0_b, 3.75/M.abs(x));
|
||||
|
|
@ -151,6 +188,7 @@ var besseli = (function() {
|
|||
|
||||
var b1_a = [0.5, 0.87890594, 0.51498869, 0.15084934, 0.2658733e-1, 0.301532e-2, 0.32411e-3].reverse();
|
||||
var b1_b = [0.39894228, -0.3988024e-1, -0.362018e-2, 0.163801e-2, -0.1031555e-1, 0.2282967e-1, -0.2895312e-1, 0.1787654e-1, -0.420059e-2].reverse();
|
||||
|
||||
function bessel1(x) {
|
||||
if(x < 3.75) return x * _horner(b1_a, x*x/(3.75*3.75));
|
||||
return (x < 0 ? -1 : 1) * M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b1_b, 3.75/M.abs(x));
|
||||
|
|
@ -159,14 +197,13 @@ var besseli = (function() {
|
|||
return function besseli(x, n) {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n == 1) return bessel1(x);
|
||||
if(n < 0) throw 'BESSELI Order (' + n + ') must be nonnegative';
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) return NaN;
|
||||
if(M.abs(x) === 0) return 0;
|
||||
if(x == Infinity) return Infinity;
|
||||
|
||||
var ret, j, tox = 2 / M.abs(x), m, bip, bi, bim;
|
||||
m=2*M.round((n+M.round(M.sqrt(40*n)))/2);
|
||||
bip=ret=0.0;
|
||||
bi=1.0;
|
||||
var ret = 0.0, j, tox = 2 / M.abs(x), bip = 0.0, bi=1.0, bim=0.0;
|
||||
var m=2*M.round((n+M.round(M.sqrt(40*n)))/2);
|
||||
for (j=m;j>0;j--) {
|
||||
bim=j*tox*bi + bip;
|
||||
bip=bi; bi=bim;
|
||||
|
|
@ -186,24 +223,24 @@ var besseli = (function() {
|
|||
var besselk = (function() {
|
||||
var b0_a = [-0.57721566, 0.42278420, 0.23069756, 0.3488590e-1, 0.262698e-2, 0.10750e-3, 0.74e-5].reverse();
|
||||
var b0_b = [1.25331414, -0.7832358e-1, 0.2189568e-1, -0.1062446e-1, 0.587872e-2, -0.251540e-2, 0.53208e-3].reverse();
|
||||
|
||||
function bessel0(x) {
|
||||
if(x <= 2) return -M.log(x/2)*besseli(x,0) + _horner(b0_a, x*x/4);
|
||||
return M.exp(-x)/M.sqrt(x)*_horner(b0_b, 2/x);
|
||||
if(x <= 2) return -M.log(x/2) * besseli(x,0) + _horner(b0_a, x*x/4);
|
||||
return M.exp(-x) / M.sqrt(x) * _horner(b0_b, 2/x);
|
||||
}
|
||||
|
||||
var b1_a = [1.0, 0.15443144, -0.67278579, -0.18156897, -0.1919402e-1, -0.110404e-2, -0.4686e-4].reverse();
|
||||
var b1_b = [1.25331414, 0.23498619, -0.3655620e-1, 0.1504268e-1, -0.780353e-2, 0.325614e-2, -0.68245e-3].reverse();
|
||||
|
||||
function bessel1(x) {
|
||||
if(x <= 2) return M.log(x/2)*besseli(x,1) + (1/x)*_horner(b1_a, x*x/4);
|
||||
if(x <= 2) return M.log(x/2) * besseli(x,1) + (1/x) * _horner(b1_a, x*x/4);
|
||||
return M.exp(-x)/M.sqrt(x)*_horner(b1_b, 2/x);
|
||||
}
|
||||
|
||||
return _bessel_wrap(bessel0, bessel1, 'BESSELK', 2, 1);
|
||||
})();
|
||||
if(typeof exports !== "undefined") {
|
||||
exports.besselj = besselj;
|
||||
exports.bessely = bessely;
|
||||
exports.besseli = besseli;
|
||||
exports.besselk = besselk;
|
||||
}
|
||||
|
||||
BESSEL.besselj = besselj;
|
||||
BESSEL.bessely = bessely;
|
||||
BESSEL.besseli = besseli;
|
||||
BESSEL.besselk = besselk;
|
||||
}));
|
||||
|
|
|
|||
10
bessel.md
10
bessel.md
|
|
@ -14,7 +14,7 @@ var M = Math;
|
|||
The methods use an approximating polynomial and evaluate using Horner's method:
|
||||
|
||||
```
|
||||
function _horner(arr, v) { return arr.reduce(function(z,w){return v * z + w;},0); };
|
||||
function _horner(arr/*:Array<number>*/, v/*:number*/) { return arr.reduce(function(z,w){return v * z + w;},0); }
|
||||
```
|
||||
|
||||
## Recurrence
|
||||
|
|
@ -29,7 +29,7 @@ So rather than go back and try to find solution for each order, we will build
|
|||
solutions for `n=0` and `n=1` and then apply the recurrence. The helper:
|
||||
|
||||
```js
|
||||
function _bessel_iter(x, n, f0, f1, sign) {
|
||||
function _bessel_iter(x/*:number*/, n/*:number*/, f0/*:number*/, f1/*:number*/, sign/*:?number*/) {
|
||||
if(!sign) sign = -1;
|
||||
var tdx = 2 / x, f2;
|
||||
if(n === 0) return f0;
|
||||
|
|
@ -48,7 +48,7 @@ checks since `Y_n` is undefined at 0 and `K_n` is real only when `x>0`
|
|||
|
||||
```
|
||||
function _bessel_wrap(bessel0, bessel1, name, nonzero, sign) {
|
||||
return function bessel(x,n) {
|
||||
return function bessel(x/*:number*/,n/*:number*/) {
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) throw name + ': Order (' + n + ') must be nonnegative';
|
||||
|
|
@ -143,7 +143,7 @@ For large values of x, the aforementioned iteration is fine, but for small
|
|||
values the expressions quickly blow up. Hence a more careful iteration is used:
|
||||
|
||||
```
|
||||
return function besselj(x, n) {
|
||||
return function besselj(x/*:number*/, n/*:number*/) {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(M.abs(x));
|
||||
if(n === 1) return bessel1(M.abs(x));
|
||||
|
|
@ -246,7 +246,7 @@ var besseli = (function() {
|
|||
return (x < 0 ? -1 : 1) * M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b1_b, 3.75/M.abs(x));
|
||||
}
|
||||
|
||||
return function besseli(x, n) {
|
||||
return function besseli(x/*:number*/, n/*:number*/) {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n == 1) return bessel1(x);
|
||||
|
|
|
|||
|
|
@ -0,0 +1,25 @@
|
|||
/* bessel.js (C) 2013-present SheetJS -- http://sheetjs.com */
|
||||
/* vim: set ts=2: */
|
||||
/*exported BESSEL */
|
||||
var BESSEL;
|
||||
/*:: declare var DO_NOT_EXPORT_BESSEL: any; */
|
||||
/*:: declare var define: any; */
|
||||
(function (factory) {
|
||||
/*jshint ignore:start */
|
||||
if(typeof DO_NOT_EXPORT_BESSEL === 'undefined') {
|
||||
if('object' === typeof exports) {
|
||||
factory(exports);
|
||||
} else if ('function' === typeof define && define.amd) {
|
||||
define(function () {
|
||||
var module = {};
|
||||
factory(module);
|
||||
return module;
|
||||
});
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
} else {
|
||||
factory(BESSEL = {});
|
||||
}
|
||||
/*jshint ignore:end */
|
||||
}(function(BESSEL) {
|
||||
|
|
@ -0,0 +1 @@
|
|||
BESSEL.version = '0.3.0';
|
||||
|
|
@ -0,0 +1,2 @@
|
|||
var M = Math;
|
||||
|
||||
|
|
@ -0,0 +1,4 @@
|
|||
/*::
|
||||
type BesselN = (x:number) => number;
|
||||
type BesselF = (x:number, n:number) => number;
|
||||
*/
|
||||
|
|
@ -0,0 +1 @@
|
|||
function _horner(arr/*:Array<number>*/, v/*:number*/)/*:number*/ { for(var i = 0, z = 0; i < arr.length; ++i) z = v * z + arr[i]; return z; }
|
||||
|
|
@ -0,0 +1,10 @@
|
|||
function _bessel_iter(x/*:number*/, n/*:number*/, f0/*:number*/, f1/*:number*/, sign/*:number*/)/*:number*/ {
|
||||
if(n === 0) return f0;
|
||||
if(n === 1) return f1;
|
||||
var tdx = 2 / x, f2 = f1;
|
||||
for(var o = 1; o < n; ++o) {
|
||||
f2 = f1 * o * tdx + sign * f0;
|
||||
f0 = f1; f1 = f2;
|
||||
}
|
||||
return f2;
|
||||
}
|
||||
|
|
@ -0,0 +1,14 @@
|
|||
function _bessel_wrap(bessel0/*:BesselN*/, bessel1/*:BesselN*/, name/*:string*/, nonzero/*:number*/, sign/*:number*/)/*:BesselF*/ {
|
||||
return function bessel(x/*:number*/,n/*:number*/) {
|
||||
if(nonzero) {
|
||||
if(x === 0) return (nonzero == 1 ? -Infinity : Infinity);
|
||||
else if(x < 0) return NaN;
|
||||
}
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) return NaN;
|
||||
n|=0;
|
||||
var b0/*:number*/ = bessel0(x), b1/*:number*/ = bessel1(x);
|
||||
return _bessel_iter(x, n, b0, b1, sign);
|
||||
};
|
||||
}
|
||||
|
|
@ -0,0 +1,83 @@
|
|||
var besselj/*:BesselF*/ = (function() {
|
||||
var W = 0.636619772; // 2 / Math.PI
|
||||
|
||||
var b0_a1a = [57568490574.0, -13362590354.0, 651619640.7, -11214424.18, 77392.33017, -184.9052456].reverse();
|
||||
var b0_a2a = [57568490411.0, 1029532985.0, 9494680.718, 59272.64853, 267.8532712, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934935152e-7].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x * x;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
a = a1 / a2;
|
||||
} else {
|
||||
var xx = x - 0.785398164;
|
||||
y = 64 / y;
|
||||
a1 = _horner(b0_a1b, y);
|
||||
a2 = _horner(b0_a2b, y);
|
||||
a = M.sqrt(W/x)*(M.cos(xx)*a1-M.sin(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
var b1_a1a = [72362614232.0, -7895059235.0, 242396853.1, -2972611.439, 15704.48260, -30.16036606].reverse();
|
||||
var b1_a2a = [144725228442.0, 2300535178.0, 18583304.74, 99447.43394, 376.9991397, 1.0].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = M.abs(x) - 2.356194491;
|
||||
if(Math.abs(x)< 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
a = a1 / a2;
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1=_horner(b1_a1b, y);
|
||||
a2=_horner(b1_a2b, y);
|
||||
a=M.sqrt(W/M.abs(x))*(M.cos(xx)*a1-M.sin(xx)*a2*8/M.abs(x));
|
||||
if(x < 0) a = -a;
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
return function besselj(x/*:number*/, n/*:number*/)/*:number*/ {
|
||||
n = Math.round(n);
|
||||
if(!isFinite(x)) return isNaN(x) ? x : 0;
|
||||
if(n < 0) return ((n%2)?-1:1)*besselj(x, -n);
|
||||
if(x < 0) return ((n%2)?-1:1)*besselj(-x, n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(x === 0) return 0;
|
||||
|
||||
var ret=0.0;
|
||||
if(x > n) {
|
||||
ret = _bessel_iter(x, n, bessel0(x), bessel1(x),-1);
|
||||
} else {
|
||||
var m=2*M.floor((n+M.floor(M.sqrt(40*n)))/2);
|
||||
var jsum=false;
|
||||
var bjp=0.0, sum=0.0;
|
||||
var bj=1.0, bjm = 0.0;
|
||||
var tox = 2 / x;
|
||||
for (var j=m;j>0;j--) {
|
||||
bjm=j*tox*bj-bjp;
|
||||
bjp=bj;
|
||||
bj=bjm;
|
||||
if (M.abs(bj) > 1E10) {
|
||||
bj *= 1E-10;
|
||||
bjp *= 1E-10;
|
||||
ret *= 1E-10;
|
||||
sum *= 1E-10;
|
||||
}
|
||||
if (jsum) sum += bj;
|
||||
jsum=!jsum;
|
||||
if (j == n) ret=bjp;
|
||||
}
|
||||
sum=2.0*sum-bj;
|
||||
ret /= sum;
|
||||
}
|
||||
return ret;
|
||||
};
|
||||
})();
|
||||
|
|
@ -0,0 +1,45 @@
|
|||
var bessely/*:BesselF*/ = (function() {
|
||||
var W = 0.636619772;
|
||||
|
||||
var b0_a1a = [-2957821389.0, 7062834065.0, -512359803.6, 10879881.29, -86327.92757, 228.4622733].reverse();
|
||||
var b0_a2a = [40076544269.0, 745249964.8, 7189466.438, 47447.26470, 226.1030244, 1.0].reverse();
|
||||
var b0_a1b = [1.0, -0.1098628627e-2, 0.2734510407e-4, -0.2073370639e-5, 0.2093887211e-6].reverse();
|
||||
var b0_a2b = [-0.1562499995e-1, 0.1430488765e-3, -0.6911147651e-5, 0.7621095161e-6, -0.934945152e-7].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x * x, xx = x - 0.785398164;
|
||||
if(x < 8) {
|
||||
a1 = _horner(b0_a1a, y);
|
||||
a2 = _horner(b0_a2a, y);
|
||||
a = a1/a2 + W * besselj(x,0) * M.log(x);
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1 = _horner(b0_a1b, y);
|
||||
a2 = _horner(b0_a2b, y);
|
||||
a = M.sqrt(W/x)*(M.sin(xx)*a1+M.cos(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
var b1_a1a = [-0.4900604943e13, 0.1275274390e13, -0.5153438139e11, 0.7349264551e9, -0.4237922726e7, 0.8511937935e4].reverse();
|
||||
var b1_a2a = [0.2499580570e14, 0.4244419664e12, 0.3733650367e10, 0.2245904002e8, 0.1020426050e6, 0.3549632885e3, 1].reverse();
|
||||
var b1_a1b = [1.0, 0.183105e-2, -0.3516396496e-4, 0.2457520174e-5, -0.240337019e-6].reverse();
|
||||
var b1_a2b = [0.04687499995, -0.2002690873e-3, 0.8449199096e-5, -0.88228987e-6, 0.105787412e-6].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
var a=0, a1=0, a2=0, y = x*x, xx = x - 2.356194491;
|
||||
if(x < 8) {
|
||||
a1 = x*_horner(b1_a1a, y);
|
||||
a2 = _horner(b1_a2a, y);
|
||||
a = a1/a2 + W * (besselj(x,1) * M.log(x) - 1 / x);
|
||||
} else {
|
||||
y = 64 / y;
|
||||
a1=_horner(b1_a1b, y);
|
||||
a2=_horner(b1_a2b, y);
|
||||
a=M.sqrt(W/x)*(M.sin(xx)*a1+M.cos(xx)*a2*8/x);
|
||||
}
|
||||
return a;
|
||||
}
|
||||
|
||||
return _bessel_wrap(bessel0, bessel1, 'BESSELY', 1, -1);
|
||||
})();
|
||||
|
|
@ -0,0 +1,43 @@
|
|||
var besseli/*:BesselF*/ = (function() {
|
||||
var b0_a = [1.0, 3.5156229, 3.0899424, 1.2067492, 0.2659732, 0.360768e-1, 0.45813e-2].reverse();
|
||||
var b0_b = [0.39894228, 0.1328592e-1, 0.225319e-2, -0.157565e-2, 0.916281e-2, -0.2057706e-1, 0.2635537e-1, -0.1647633e-1, 0.392377e-2].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
if(x <= 3.75) return _horner(b0_a, x*x/(3.75*3.75));
|
||||
return M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b0_b, 3.75/M.abs(x));
|
||||
}
|
||||
|
||||
var b1_a = [0.5, 0.87890594, 0.51498869, 0.15084934, 0.2658733e-1, 0.301532e-2, 0.32411e-3].reverse();
|
||||
var b1_b = [0.39894228, -0.3988024e-1, -0.362018e-2, 0.163801e-2, -0.1031555e-1, 0.2282967e-1, -0.2895312e-1, 0.1787654e-1, -0.420059e-2].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
if(x < 3.75) return x * _horner(b1_a, x*x/(3.75*3.75));
|
||||
return (x < 0 ? -1 : 1) * M.exp(M.abs(x))/M.sqrt(M.abs(x))*_horner(b1_b, 3.75/M.abs(x));
|
||||
}
|
||||
|
||||
return function besseli(x/*:number*/, n/*:number*/)/*:number*/ {
|
||||
n = Math.round(n);
|
||||
if(n === 0) return bessel0(x);
|
||||
if(n === 1) return bessel1(x);
|
||||
if(n < 0) return NaN;
|
||||
if(M.abs(x) === 0) return 0;
|
||||
if(x == Infinity) return Infinity;
|
||||
|
||||
var ret = 0.0, j, tox = 2 / M.abs(x), bip = 0.0, bi=1.0, bim=0.0;
|
||||
var m=2*M.round((n+M.round(M.sqrt(40*n)))/2);
|
||||
for (j=m;j>0;j--) {
|
||||
bim=j*tox*bi + bip;
|
||||
bip=bi; bi=bim;
|
||||
if (M.abs(bi) > 1E10) {
|
||||
bi *= 1E-10;
|
||||
bip *= 1E-10;
|
||||
ret *= 1E-10;
|
||||
}
|
||||
if(j == n) ret = bip;
|
||||
}
|
||||
ret *= besseli(x, 0) / bi;
|
||||
return x < 0 && (n%2) ? -ret : ret;
|
||||
};
|
||||
|
||||
})();
|
||||
|
||||
|
|
@ -0,0 +1,19 @@
|
|||
var besselk/*:BesselF*/ = (function() {
|
||||
var b0_a = [-0.57721566, 0.42278420, 0.23069756, 0.3488590e-1, 0.262698e-2, 0.10750e-3, 0.74e-5].reverse();
|
||||
var b0_b = [1.25331414, -0.7832358e-1, 0.2189568e-1, -0.1062446e-1, 0.587872e-2, -0.251540e-2, 0.53208e-3].reverse();
|
||||
|
||||
function bessel0(x/*:number*/)/*:number*/ {
|
||||
if(x <= 2) return -M.log(x/2) * besseli(x,0) + _horner(b0_a, x*x/4);
|
||||
return M.exp(-x) / M.sqrt(x) * _horner(b0_b, 2/x);
|
||||
}
|
||||
|
||||
var b1_a = [1.0, 0.15443144, -0.67278579, -0.18156897, -0.1919402e-1, -0.110404e-2, -0.4686e-4].reverse();
|
||||
var b1_b = [1.25331414, 0.23498619, -0.3655620e-1, 0.1504268e-1, -0.780353e-2, 0.325614e-2, -0.68245e-3].reverse();
|
||||
|
||||
function bessel1(x/*:number*/)/*:number*/ {
|
||||
if(x <= 2) return M.log(x/2) * besseli(x,1) + (1/x) * _horner(b1_a, x*x/4);
|
||||
return M.exp(-x)/M.sqrt(x)*_horner(b1_b, 2/x);
|
||||
}
|
||||
|
||||
return _bessel_wrap(bessel0, bessel1, 'BESSELK', 2, 1);
|
||||
})();
|
||||
|
|
@ -0,0 +1,4 @@
|
|||
BESSEL.besselj = besselj;
|
||||
BESSEL.bessely = bessely;
|
||||
BESSEL.besseli = besseli;
|
||||
BESSEL.besselk = besselk;
|
||||
|
|
@ -0,0 +1 @@
|
|||
}));
|
||||
|
|
@ -0,0 +1,84 @@
|
|||
<!DOCTYPE html>
|
||||
<!-- bessel.js (C) 2013-present SheetJS http://sheetjs.com -->
|
||||
<!-- vim: set ts=2: -->
|
||||
<html>
|
||||
<head>
|
||||
<title>Bessel Live Demo</title>
|
||||
<script src="bessel.js"></script>
|
||||
<script src="https://cdn.rawgit.com/SheetJS/printj/master/printj.js"></script>
|
||||
</head>
|
||||
<body>
|
||||
<b>Bessel Live Demo</b><br />
|
||||
<a href="https://github.com/SheetJS/bessel">Source Code</a><br />
|
||||
<a href="https://github.com/SheetJS/bessel/issues">Issues? Something look weird? Click here and report an issue</a><br />
|
||||
<br />
|
||||
<table>
|
||||
<tr><td><b>Number:</b></td><td><input type="text" id="val" value="1.5"></td></tr>
|
||||
<tr><td><b>Order:</b></td><td><input type="text" id="ord" value="1"></td></tr>
|
||||
<tr><td colspan="2"> </td></tr>
|
||||
<tr><td colspan="2"><b><pre id="jt">J</pre></b></td></tr>
|
||||
<tr><td colspan="2"><b><pre id="yt">Y</pre></b></td></tr>
|
||||
<tr><td colspan="2"><b><pre id="it">I</pre></b></td></tr>
|
||||
<tr><td colspan="2"><b><pre id="kt">K</pre></b></td></tr>
|
||||
</table>
|
||||
</body>
|
||||
<script>
|
||||
/*jshint browser:true */
|
||||
var V = document.getElementById('val');
|
||||
var O = document.getElementById('ord');
|
||||
var J = document.getElementById('jt');
|
||||
var Y = document.getElementById('yt');
|
||||
var I = document.getElementById('it');
|
||||
var K = document.getElementById('kt');
|
||||
|
||||
function update() {
|
||||
var v = Number(V.value);
|
||||
var o = Number(O.value);
|
||||
if(v !== v) { J.innerHTML=Y.innerHTML=I.innerHTML=K.innerHTML = "Value |" + V.value + "| not valid"; return; }
|
||||
if(o !== o) { J.innerHTML=Y.innerHTML=I.innerHTML=K.innerHTML = "Order |" + O.value + "| not valid"; return; }
|
||||
o |= 0;
|
||||
try { J.innerHTML = PRINTJ.sprintf("J(%.11g, %d) = %+.11g", v, o, BESSEL.besselj(v, o)); } catch(e) { J.innerHTML = String(e); }
|
||||
try { Y.innerHTML = PRINTJ.sprintf("Y(%.11g, %d) = %+.11g", v, o, BESSEL.bessely(v, o)); } catch(e) { Y.innerHTML = String(e); }
|
||||
try { I.innerHTML = PRINTJ.sprintf("I(%.11g, %d) = %+.11g", v, o, BESSEL.besseli(v, o)); } catch(e) { I.innerHTML = String(e); }
|
||||
try { K.innerHTML = PRINTJ.sprintf("K(%.11g, %d) = %+.11g", v, o, BESSEL.besselk(v, o)); } catch(e) { K.innerHTML = String(e); }
|
||||
}
|
||||
|
||||
/* Bind to relevant events */
|
||||
if(O.addEventListener) {
|
||||
O.addEventListener('keyup', update);
|
||||
V.addEventListener('keyup', update);
|
||||
J.addEventListener('change', update);
|
||||
Y.addEventListener('change', update);
|
||||
I.addEventListener('change', update);
|
||||
K.addEventListener('change', update);
|
||||
} else if(O.attachEvent) {
|
||||
O.attachEvent('onkeyup', update);
|
||||
V.attachEvent('onkeyup', update);
|
||||
J.attachEvent('onchange', update);
|
||||
Y.attachEvent('onchange', update);
|
||||
I.attachEvent('onchange', update);
|
||||
K.attachEvent('onchange', update);
|
||||
} else {
|
||||
O.oninput = update;
|
||||
V.oninput = update;
|
||||
J.onchange = update;
|
||||
Y.onchange = update;
|
||||
I.onchange = update;
|
||||
K.onchange = update;
|
||||
}
|
||||
|
||||
window.onload = function() { update(); };
|
||||
</script>
|
||||
<script type="text/javascript">
|
||||
var _gaq = _gaq || [];
|
||||
_gaq.push(['_setAccount', 'UA-36810333-1']);
|
||||
_gaq.push(['_trackPageview']);
|
||||
|
||||
(function() {
|
||||
var ga = document.createElement('script'); ga.type = 'text/javascript'; ga.async = true;
|
||||
ga.src = ('https:' == document.location.protocol ? 'https://ssl' : 'http://www') + '.google-analytics.com/ga.js';
|
||||
var s = document.getElementsByTagName('script')[0]; s.parentNode.insertBefore(ga, s);
|
||||
})();
|
||||
</script>
|
||||
</html>
|
||||
|
||||
|
|
@ -0,0 +1,42 @@
|
|||
#!/bin/bash
|
||||
# make_help.sh -- process listing of targets and special items in Makefile
|
||||
# Copyright (C) 2016-present SheetJS
|
||||
#
|
||||
# usage in makefile: pipe the output of the following command:
|
||||
# @grep -hE '(^[a-zA-Z_-][ a-zA-Z_-]*:.*?|^#[#*])' $(MAKEFILE_LIST)
|
||||
#
|
||||
# lines starting with "## " are treated as subtitles
|
||||
# lines starting with "#* " are treated as plaintext comments
|
||||
# multiple targets with "## " after the ":" are rendered as separate targets
|
||||
# if the presumed default target is labeled, it will be assigned a unique color
|
||||
|
||||
awk '
|
||||
BEGIN{recipes=0;}
|
||||
!/#[#*] .*$/ {next;}
|
||||
{multi=0; isrecipe=0;}
|
||||
/^[^#]*:/ {isrecipe=1; ++recipes;}
|
||||
/^[^ :]* .*:/ {multi=1}
|
||||
multi==0 && isrecipe>0 { if(recipes > 1) print; else print $0, "[default]"; next}
|
||||
isrecipe == 0 {print; next}
|
||||
multi>0 {
|
||||
k=split($0, msg, "##"); m=split($0, a, ":"); n=split(a[1], b, " ");
|
||||
for(i=1; i<=n; ++i) print b[i] ":", "##" msg[2], (recipes==1 && i==1 ? "[default]" : "")
|
||||
}
|
||||
END {}
|
||||
' | if [[ -t 1 ]]; then
|
||||
awk '
|
||||
BEGIN {FS = ":.*?## "}
|
||||
{color=36}
|
||||
/\[default\]/ {color=35}
|
||||
NF==1 && /^##/ {color=34}
|
||||
NF==1 && /^#\*/ {color=20; $1 = substr($1, 4)}
|
||||
{printf "\033[" color "m%-20s\033[0m %s\n", $1, $2;}
|
||||
END{}' -
|
||||
else
|
||||
awk '
|
||||
BEGIN {FS = ":.*?## "}
|
||||
/^#\* / {$1 = substr($1, 4)}
|
||||
{printf "%-20s %s\n", $1, $2;}
|
||||
END{}' -
|
||||
fi
|
||||
|
||||
|
|
@ -0,0 +1,14 @@
|
|||
#!/bin/bash
|
||||
# spin.sh -- show a spinner (for coverage test)
|
||||
# Copyright (C) 2014-present SheetJS
|
||||
|
||||
wpid=$1
|
||||
delay=1
|
||||
str="|/-\\"
|
||||
while [ $(ps -a|awk '$1=='$wpid' {print $1}') ]; do
|
||||
t=${str#?}
|
||||
printf " [%c]" "$str"
|
||||
str=$t${str%"$t"}
|
||||
sleep $delay
|
||||
printf "\b\b\b\b"
|
||||
done
|
||||
29
package.json
29
package.json
|
|
@ -1,19 +1,26 @@
|
|||
{
|
||||
"name": "bessel",
|
||||
"version": "0.2.0",
|
||||
"version": "0.3.0",
|
||||
"author": "SheetJS",
|
||||
"description": "Bessel Functions in pure JS",
|
||||
"keywords": ["bessel", "math", "specfun"],
|
||||
"main": "bessel.js",
|
||||
"repository": { "type": "git", "url": "https://github.com/SheetJS/bessel" },
|
||||
"license": "MIT",
|
||||
"description": "Pure-JS Bessel Functions",
|
||||
"keywords": [ "bessel", "math", "specfun" ],
|
||||
"main": "./bessel",
|
||||
"dependencies": {
|
||||
},
|
||||
"devDependencies": {
|
||||
"mocha":"",
|
||||
"voc":""
|
||||
},
|
||||
"repository": { "type":"git", "url":"git://github.com/SheetJS/bessel.git" },
|
||||
"scripts": {
|
||||
"test": "make test"
|
||||
},
|
||||
"dependencies": {
|
||||
"voc":""
|
||||
"config": {
|
||||
"blanket": {
|
||||
"pattern": "bessel.js"
|
||||
}
|
||||
},
|
||||
"devDependencies": {
|
||||
"mocha":""
|
||||
}
|
||||
"bugs": { "url": "https://github.com/SheetJS/bessel/issues" },
|
||||
"license": "MIT",
|
||||
"engines": { "node": ">=0.8" }
|
||||
}
|
||||
|
|
|
|||
|
|
@ -1,4 +0,0 @@
|
|||
<script src="bessel.js"></script>
|
||||
<script>
|
||||
console.log(besselj(1,2));
|
||||
</script>
|
||||
142
test.js
142
test.js
|
|
@ -1,22 +1,124 @@
|
|||
var assert = require('assert');
|
||||
var bessel;
|
||||
describe('When using bessel functions', function() {
|
||||
before(function() {
|
||||
bessel = require('./bessel');
|
||||
});
|
||||
it('It must load besseli', function() {
|
||||
assert.equal(typeof bessel.besseli, 'function');
|
||||
});
|
||||
it('It must compute the Modified Bessel function at 1.5 with an order of 1 (0.981666)', function() {
|
||||
assert(Math.floor((bessel.besseli(1.5, 1) - 0.981666) * 10e6) < 10);
|
||||
});
|
||||
it('It must compute the Bessel function at 1.9 with an order of 2 (0.329926)', function() {
|
||||
assert(Math.floor((bessel.besselj(1.9, 2) - 0.329926) * 10e6) < 10);
|
||||
});
|
||||
it('It must compute the Modified Bessel function at 1.5 with an order of 1 (0.277388)', function() {
|
||||
assert(Math.floor((bessel.besselk(1.5, 1) - 0.277388) * 10e6) < 10);
|
||||
});
|
||||
it('It must compute Weber\'s Bessel function at 2.5 with an order of 1 (0.145918)', function() {
|
||||
assert(Math.floor((bessel.besselk(2.5, 1) - 0.981666) * 10e6) < 10);
|
||||
/* vim: set ts=2: */
|
||||
var X;
|
||||
if(typeof require !== 'undefined') {
|
||||
assert = require('assert');
|
||||
describe('source',function(){
|
||||
it('should load',function(){X=require('./');});
|
||||
it('all bessel functions', function() {
|
||||
assert.equal(typeof X.besseli, 'function');
|
||||
assert.equal(typeof X.besselj, 'function');
|
||||
assert.equal(typeof X.besselk, 'function');
|
||||
assert.equal(typeof X.bessely, 'function');
|
||||
});
|
||||
});
|
||||
} else { X = BESSEL; }
|
||||
|
||||
function approx(a, b, t) { return Math.abs((a - b) / (Math.abs(a) + Math.abs(b))) <= t; }
|
||||
|
||||
function chk(func, x, n, v, F, t) {
|
||||
var V = func(x,n);
|
||||
if(!approx(V, v, t) && v != V && !(isNaN(v) && isNaN(V))) throw new Error(F + "(" + x + "," + n + ") = " + V + " !~ " + v);
|
||||
}
|
||||
|
||||
var funcs = ['besseli', 'besselj', 'besselk', 'bessely'];
|
||||
|
||||
var tests = {
|
||||
'besseli': [
|
||||
[1.5, 1, 0.981666],
|
||||
[1.9, 2, 0.603272435],
|
||||
[2.5, 1, 2.516716242]
|
||||
],
|
||||
'besselj': [
|
||||
[1.5, 1, 0.557936508],
|
||||
[1.9, 2, 0.329925829],
|
||||
[2.5, 1, 0.497094103]
|
||||
],
|
||||
'besselk': [
|
||||
[1.5, 1, 0.277388],
|
||||
[1.9, 2, 0.296909301],
|
||||
[2.5, 1, 0.073890816]
|
||||
],
|
||||
'bessely': [
|
||||
[1.5, 1, -0.412308627],
|
||||
[1.9, 2, -0.669878674],
|
||||
[2.5, 1, 0.14591814]
|
||||
] };
|
||||
|
||||
describe('correctness', function() {
|
||||
funcs.forEach(function(F) { it(F, function() {
|
||||
tests[F].forEach(function(t) { chk(X[F], t[0], t[1], t[2], F, 1e-6); });
|
||||
});});
|
||||
|
||||
if(typeof require != 'undefined') {
|
||||
var fs = require('fs');
|
||||
it('excel', function() {
|
||||
var xl = fs.readFileSync("test_files/excel.tsv", 'ascii').split("\n").map(function(l) { return l.split("\t");});
|
||||
xl.forEach(function(l) { if(l.length < 6) return;
|
||||
var x=Number(l[0]), n=Number(l[1]);
|
||||
var i=Number(l[2]), j=Number(l[3]), k=Number(l[4]), y=Number(l[5]);
|
||||
chk(X.besseli, x, n, i, 'besseli', 1e-4);
|
||||
chk(X.besselj, x, n, j, 'besselj', 1e-4);
|
||||
chk(X.besselk, x, n, k, 'besselk', 1e-4);
|
||||
chk(X.bessely, x, n, y, 'bessely', 1e-4);
|
||||
});
|
||||
});
|
||||
it('mma', function() {
|
||||
var mma = fs.readFileSync("test_files/mma.tsv", 'ascii').split("\n").map(function(l) { return l.split("\t");});
|
||||
mma.forEach(function(l) { if(l.length < 6) return;
|
||||
var x=Number(l[0]), n=Number(l[1]);
|
||||
var i=Number(l[2]), j=Number(l[3]), k=Number(l[4]), y=Number(l[5]);
|
||||
chk(X.besseli, x, n, i, 'besseli', 1e-5);
|
||||
chk(X.besselj, x, n, j, 'besselj', 1e-5);
|
||||
chk(X.besselk, x, n, k, 'besselk', 1e-5);
|
||||
chk(X.bessely, x, n, y, 'bessely', 1e-5);
|
||||
});
|
||||
});
|
||||
}
|
||||
});
|
||||
|
||||
var rand = [
|
||||
0.3946086812775882, 3.946086812775882, 39.46086812775882,
|
||||
0.03119419917289612, 0.3119419917289612, 3.119419917289612,
|
||||
0.18311903629136062, 1.8311903629136062, 18.311903629136062,
|
||||
0.5317681126599452, 5.317681126599452, 53.17681126599452,
|
||||
0.37724541457032235, 3.7724541457032235, 37.724541457032235,
|
||||
0.8287802926845655, 8.287802926845655, 82.87802926845655,
|
||||
0.7197724866379658, 7.197724866379658, 71.97724866379658
|
||||
];
|
||||
|
||||
describe('properties', function() {
|
||||
var t = 1e-4;
|
||||
it('limiting behavior', function() {
|
||||
chk(X.besselj, 0, 0, 1, 'besselj', 1e-6);
|
||||
chk(X.besseli, 0, 0, 1, 'besseli', 1e-6);
|
||||
chk(X.besselk, 0, 0, Infinity, 'besselk', 1e-6);
|
||||
chk(X.bessely, 0, 0, -Infinity, 'bessely', 1e-6);
|
||||
for(var i = 1; i < 20; ++i) {
|
||||
chk(X.besselj, 0, i, 0, 'besselj', 1e-6);
|
||||
chk(X.besseli, 0, i, 0, 'besseli', 1e-6);
|
||||
}
|
||||
});
|
||||
it('besselj', function() {
|
||||
var F = "BESSELJ", f = X.besselj;
|
||||
rand.forEach(function(r) {
|
||||
for(var n = 0; n < 20; ++n) {
|
||||
var pp = f( r, n), pn = f( r, -n), np = f(-r, n), nn = f(-r, -n);
|
||||
if((n%2)) {
|
||||
/* besselj odd if n is odd */
|
||||
if(!approx(np, -pp, t)) throw new Error(F + " np[" + np + "] != -pp[" + (-pp) + "] (" + r + "," + n + ")");
|
||||
if(!approx(nn, -pn, t)) throw new Error(F + " nn[" + nn + "] != -pn[" + (-pn) + "] (" + r + "," + n + ")");
|
||||
/* asymmetric in n */
|
||||
if(!approx(pn, -pp, t)) throw new Error(F + " pn[" + pn + "] != -pp[" + (-pp) + "] (" + r + "," + n + ")");
|
||||
if(!approx(nn, -np, t)) throw new Error(F + " nn[" + nn + "] != -pn[" + (-pn) + "] (" + r + "," + n + ")");
|
||||
} else {
|
||||
/* besselj even if n is even */
|
||||
if(!approx(np, pp, t)) throw new Error(F + " np[" + np + "] != pp[" + (pp) + "] (" + r + "," + n + ")");
|
||||
if(!approx(nn, pn, t)) throw new Error(F + " nn[" + nn + "] != pn[" + (pn) + "] (" + r + "," + n + ")");
|
||||
/* symmetric in n */
|
||||
if(!approx(pn, pp, t)) throw new Error(F + " pn[" + pn + "] != pp[" + (pp) + "] (" + r + "," + n + ")");
|
||||
if(!approx(nn, np, t)) throw new Error(F + " nn[" + nn + "] != pn[" + (pn) + "] (" + r + "," + n + ")");
|
||||
}
|
||||
}
|
||||
});
|
||||
});
|
||||
});
|
||||
|
|
|
|||
|
|
@ -0,0 +1,60 @@
|
|||
0.5 1 0.257894303 0.242268458 1.656441128 -1.471472392
|
||||
1 1 0.565159098 0.440050586 0.601907232 -0.781212821
|
||||
1.5 1 0.981666428 0.557936508 0.277387804 -0.412308627
|
||||
2 1 1.590636857 0.576724808 0.13986588 -0.107032432
|
||||
2.5 1 2.516716242 0.497094103 0.073890816 0.145918138
|
||||
3 1 3.953370217 0.339058958 0.040156431 0.324674424
|
||||
3.5 1 6.205834932 0.137377527 0.022239393 0.410188417
|
||||
4 1 9.759465157 -0.066043328 0.012483499 0.397925709
|
||||
4.5 1 15.38922293 -0.231060432 0.007078095 0.300997321
|
||||
5 1 24.33564185 -0.327579139 0.004044613 0.14786314
|
||||
5.5 1 38.58816378 -0.341438216 0.002325569 -0.023758243
|
||||
6 1 61.34193694 -0.276683859 0.00134392 -0.175010349
|
||||
6.5 1 97.73501304 -0.153841304 0.000779894 -0.274091281
|
||||
7 1 156.0390965 -0.004682826 0.000454182 -0.302667244
|
||||
7.5 1 249.5843679 0.135248424 0.000265297 -0.25912852
|
||||
8 1 399.8731348 0.234628387 0.000155369 -0.158060462
|
||||
8.5 1 641.6198938 0.273115291 9.11973E-05 -0.02616868
|
||||
9 1 1030.914709 0.245307412 5.3637E-05 0.104314575
|
||||
9.5 1 1658.453068 0.161262274 3.1602E-05 0.203179899
|
||||
10 1 2670.988321 0.04347225 1.86488E-05 0.249015424
|
||||
0.5 2 0.031906148 0.030604023 7.550183547 -5.441370834
|
||||
1 2 0.135747667 0.114903485 1.624838884 -1.650682613
|
||||
1.5 2 0.337834621 0.232087679 0.583655974 -0.932193761
|
||||
2 2 0.688948449 0.352834208 0.25375976 -0.617408098
|
||||
2.5 2 1.276466159 0.446059058 0.121460207 -0.381335848
|
||||
3 2 2.245212431 0.486091263 0.061510459 -0.160400399
|
||||
3.5 2 3.832012072 0.458629185 0.032307122 0.045371436
|
||||
4 2 6.422189499 0.364128143 0.017401426 0.215903599
|
||||
4.5 2 10.64151731 0.217848984 0.009545677 0.328481596
|
||||
5 2 17.50561501 0.046565119 0.005308944 0.367662879
|
||||
5.5 2 28.66258571 -0.117315483 0.00298437 0.330841237
|
||||
6 2 46.7870949 -0.242873213 0.001691968 0.229857908
|
||||
6.5 2 76.22054669 -0.307430389 0.000965899 0.088906659
|
||||
7 2 124.01131 -0.301417224 0.000554562 -0.0605266
|
||||
7.5 2 201.6054825 -0.230273411 0.000319924 -0.186414216
|
||||
8 2 327.5958389 -0.112993711 0.000185313 -0.263036605
|
||||
8.5 2 532.1925545 0.022323169 0.000107716 -0.276362442
|
||||
9 2 864.4962206 0.144846369 6.28006E-05 -0.226755682
|
||||
9.5 2 1404.333008 0.2278787 3.67109E-05 -0.128435911
|
||||
10 2 2281.518997 0.254630214 2.15098E-05 -0.005868083
|
||||
0.5 3 0.002645112 0.00256373 62.0579095 -42.05949428
|
||||
1 3 0.022168424 0.019563354 7.101262769 -5.821517632
|
||||
1.5 3 0.080774114 0.060963951 1.833803735 -2.073541401
|
||||
2 3 0.21273996 0.12894325 0.6473854 -1.127783765
|
||||
2.5 3 0.474370413 0.216600403 0.268227146 -0.756055495
|
||||
3 3 0.959753625 0.309062865 0.122170376 -0.538541623
|
||||
3.5 3 1.826392593 0.386770113 0.059161818 -0.358335346
|
||||
4 3 3.337275843 0.430171471 0.029884924 -0.18202211
|
||||
4.5 3 5.930096274 0.424703973 0.015563141 -0.009013679
|
||||
5 3 10.3311502 0.364831234 0.008291768 0.146267163
|
||||
5.5 3 17.7426483 0.256117865 0.00449602 0.264370052
|
||||
6 3 30.15054042 0.114768384 0.002471898 0.328248955
|
||||
6.5 3 50.83005827 -0.035346628 0.001374294 0.328803071
|
||||
7 3 85.17548647 -0.167555588 0.000771075 0.268080616
|
||||
7.5 3 142.0614436 -0.25806091 0.000435923 0.159707605
|
||||
8 3 236.0752263 -0.291125242 0.000248026 0.026542159
|
||||
8.5 3 391.1763672 -0.26261027 0.000141887 -0.103884234
|
||||
9 3 646.6942118 -0.180931248 8.15484E-05 -0.205094878
|
||||
9.5 3 1067.155009 -0.065313348 4.70593E-05 -0.257258177
|
||||
10 3 1758.380739 0.058379836 2.72527E-05 -0.251362657
|
||||
|
|
|
@ -0,0 +1,60 @@
|
|||
0.5 1 0.2578943054 0.2422684577 1.65644112 -1.471472393
|
||||
1 1 0.565159104 0.4400505857 0.6019072302 -0.7812128213
|
||||
1.5 1 0.9816664286 0.5579365079 0.2773878005 -0.412308627
|
||||
2 1 1.590636855 0.5767248078 0.1398658818 -0.1070324315
|
||||
2.5 1 2.516716245 0.4970941025 0.07389081635 0.145918138
|
||||
3 1 3.953370217 0.3390589585 0.04015643113 0.3246744248
|
||||
3.5 1 6.205834922 0.1373775274 0.02223939293 0.4101884179
|
||||
4 1 9.759465154 -0.06604332802 0.01248349889 0.3979257106
|
||||
4.5 1 15.38922275 -0.2310604319 0.007078094909 0.3009973231
|
||||
5 1 24.33564214 -0.3275791376 0.004044613445 0.1478631434
|
||||
5.5 1 38.58816462 -0.3414382154 0.002325569009 -0.02375823896
|
||||
6 1 61.34193678 -0.2766838581 0.001343919718 -0.1750103443
|
||||
6.5 1 97.73501077 -0.1538413014 0.0007798943982 -0.274091274
|
||||
7 1 156.0390929 -0.004682823482 0.0004541824869 -0.302667237
|
||||
7.5 1 249.5843654 0.1352484276 0.0002652973901 -0.2591285105
|
||||
8 1 399.8731368 0.2346363469 0.0001553692118 -0.1580604617
|
||||
8.5 1 641.6199025 0.2731219637 0.00009119724775 -0.0261686794
|
||||
9 1 1030.914723 0.2453117866 0.00005363701638 0.1043145752
|
||||
9.5 1 1658.453078 0.1612644308 0.00003160203411 0.2031798994
|
||||
10 1 2670.988304 0.04347274617 0.00001864877345 0.2490154242
|
||||
0.5 2 0.03190614918 0.03060402346 7.550183551 -5.441370837
|
||||
1 2 0.1357476698 0.1149034849 1.624838899 -1.650682607
|
||||
1.5 2 0.3378346183 0.2320876721 0.5836559633 -0.9321937598
|
||||
2 2 0.6889484477 0.3528340286 0.2537597546 -0.6174081042
|
||||
2.5 2 1.276466148 0.4460590584 0.1214602063 -0.3813358492
|
||||
3 2 2.245212441 0.4860912606 0.06151045847 -0.1604003935
|
||||
3.5 2 3.832012048 0.4586291842 0.0323071217 0.04537143773
|
||||
4 2 6.422189375 0.3641281459 0.01740142553 0.2159035946
|
||||
4.5 2 10.6415173 0.2178489837 0.009545677203 0.3284815967
|
||||
5 2 17.50561497 0.04656511628 0.005308943712 0.3676628826
|
||||
5.5 2 28.66258529 -0.1173154816 0.002984370024 0.3308412333
|
||||
6 2 46.78709472 -0.24287321 0.001691967567 0.2298579025
|
||||
6.5 2 76.22054724 -0.3074303906 0.0009658992748 0.08890665832
|
||||
7 2 124.0113105 -0.3014172201 0.0005545621667 -0.06052660947
|
||||
7.5 2 201.6054807 -0.2302734105 0.0003199235871 -0.1864142223
|
||||
8 2 327.5958315 -0.1129917204 0.0001853130082 -0.2630366048
|
||||
8.5 2 532.1925382 0.02232473961 0.0001077157423 -0.2763624417
|
||||
9 2 864.496194 0.1448473415 0.00006280064993 -0.2267556816
|
||||
9.5 2 1404.332974 0.2278791542 0.00003671094477 -0.1284359105
|
||||
10 2 2281.518968 0.2546303137 0.00002150981701 -0.005868082442
|
||||
0.5 3 0.002645111969 0.002563729995 62.05790953 -42.0594943
|
||||
1 3 0.02216842492 0.01956335398 7.101262825 -5.821517606
|
||||
1.5 3 0.08077411302 0.06096395114 1.833803702 -2.073541399
|
||||
2 3 0.2127399592 0.1289432495 0.6473853909 -1.127783777
|
||||
2.5 3 0.4743704088 0.216600391 0.2682271464 -0.7560554968
|
||||
3 3 0.9597536295 0.3090627223 0.1221703758 -0.5385416161
|
||||
3.5 3 1.826392582 0.3867701117 0.05916181773 -0.3583353462
|
||||
4 3 3.337275778 0.4301714739 0.02988492442 -0.182022116
|
||||
4.5 3 5.930096266 0.424703973 0.01556314131 -0.009013681594
|
||||
5 3 10.33115017 0.3648312306 0.008291768415 0.1462671627
|
||||
5.5 3 17.74264804 0.2561178651 0.004496019935 0.264370045
|
||||
6 3 30.1505403 0.1147683848 0.002471898096 0.328248946
|
||||
6.5 3 50.83005863 -0.03534663129 0.001374293952 0.3288030637
|
||||
7 3 85.17548684 -0.167555588 0.0007710751536 0.268080603
|
||||
7.5 3 142.0614424 -0.2580609132 0.0004359233032 0.1597075919
|
||||
8 3 236.075221 -0.2911322071 0.0002480257159 0.02654215932
|
||||
8.5 3 391.1763552 -0.2626162039 0.0001418870088 -0.1038842343
|
||||
9 3 646.6941919 -0.1809351903 0.00008154841635 -0.2050948781
|
||||
9.5 3 1067.154983 -0.06531531322 0.00004705927401 -0.2572581775
|
||||
10 3 1758.380717 0.05837937931 0.00002725270026 -0.2513626572
|
||||
|
Can't render this file because it has a wrong number of fields in line 2.
|
Loading…
Reference in New Issue