Special case of besselY

paulmasson · May 28, 2022 · d994dda · d994dda
1 parent ec1346b
commit d994dda
Show file tree

Hide file tree

Showing 2 changed files with 16 additions and 0 deletions.
diff --git a/build/math.js b/build/math.js
@@ -589,6 +589,11 @@ function besselY( n, x ) {
 
   if ( isComplex(n) || isComplex(x) ) {
 
+    if ( isZero(x) ) {
+      if ( !isComplex(n) ) n = complex(n);
+      if ( fractionalPart(n.re) === -.5 && n.im === 0 ) return complex(0);
+    }
+
     // dlmf.nist.gov/10.2.3
     if ( isInteger(n) )
       return div( add( diff( n => besselJ(n,x), n ),
@@ -601,6 +606,9 @@ function besselY( n, x ) {
 
   if ( x < 0 ) return besselY( n, complex(x) );
 
+  if ( x === 0 )
+    if ( fractionalPart(n) === -.5 ) return 0;
+
   // dlmf.nist.gov/10.2.3
   if ( Number.isInteger(n) )
     return ( diff( n => besselJ(n,x), n ) + (-1)**n * diff( n => besselJ(n,x), -n ) ) / pi;

diff --git a/src/functions/bessel.js b/src/functions/bessel.js
@@ -58,6 +58,11 @@ function besselY( n, x ) {
 
   if ( isComplex(n) || isComplex(x) ) {
 
+    if ( isZero(x) ) {
+      if ( !isComplex(n) ) n = complex(n);
+      if ( fractionalPart(n.re) === -.5 && n.im === 0 ) return complex(0);
+    }
+
     // dlmf.nist.gov/10.2.3
     if ( isInteger(n) )
       return div( add( diff( n => besselJ(n,x), n ),
@@ -70,6 +75,9 @@ function besselY( n, x ) {
 
   if ( x < 0 ) return besselY( n, complex(x) );
 
+  if ( x === 0 )
+    if ( fractionalPart(n) === -.5 ) return 0;
+
   // dlmf.nist.gov/10.2.3
   if ( Number.isInteger(n) )
     return ( diff( n => besselJ(n,x), n ) + (-1)**n * diff( n => besselJ(n,x), -n ) ) / pi;