16221: implement Struve H/L functions

Author: rwst

src/sage/functions/all.py

@@ -27,7 +27,8 @@
 from transcendental import (zeta, zetaderiv, zeta_symmetric, hurwitz_zeta,
                             dickman_rho)
 
-from sage.functions.bessel import (bessel_I, bessel_J, bessel_K, bessel_Y, Bessel)
+from sage.functions.bessel import (bessel_I, bessel_J, bessel_K, bessel_Y,
+                                   Bessel, struve_H, struve_L)
 
 from special import (hypergeometric_U,
                      spherical_bessel_J, spherical_bessel_Y,

src/sage/functions/bessel.py

@@ -1666,3 +1666,205 @@ def plot(self,a,b):
             f = lambda z: _bessel_Y(nu,z,s)
             P = plot(f,a,b)
         return P
+
+class Function_Struve_H(BuiltinFunction):
+    r"""
+    The Struve functions, solutions to the non-homogeneous Bessel differential equation:
+    
+    .. math::
+    
+        x^2\frac{d^2y}{dx^2}+x\frac{dy}{dx}+(x^2-\alpha^2)y=\frac{4\bigl(\frac{x}{2}\bigr)^{\alpha+1}}{\sqrt\pi\Gamma(\alpha+\tfrac12)},
+
+    .. math::
+    
+        \mathrm{H}_\alpha(x) = y(x)
+    
+    EXAMPLES::
+    
+        sage: struve_H(-1/2,x)
+        sqrt(2)*sqrt(1/(pi*x))*sin(x)
+        sage: struve_H(2,x)
+        struve_H(2, x)
+        sage: struve_H(1/2.,pi)
+        0.900316316157106
+        
+    REFERENCES:
+    
+        - Abramowitz and Stegun: Handbook of Mathematical Functions,
+          http://www.math.sfu.ca/~cbm/aands/
+    
+        - http://en.wikipedia.org/wiki/Struve_function
+    """
+    def __init__(self):
+        r"""
+    ^    EXAMPLES::
+    
+            sage: maxima("struve_h(n,x);").sage()
+            struve_H(n, x)
+            sage: struve_H(7/5,1)._maxima_()
+            struve_h(7/5,1)
+        """
+        BuiltinFunction.__init__(self, 'struve_H', nargs=2,
+                                 conversions=dict(maple='StruveH',
+                                                  mathematica='StruveH',
+                                                  maxima='struve_h',
+                                                  sympy='struveh'))
+
+    def _eval_(self, a, z):
+        """
+        EXAMPLES::
+        
+            sage: struve_H(0,0)
+            0
+            sage: struve_H(-1/2,x)
+            sqrt(2)*sqrt(1/(pi*x))*sin(x)
+            sage: struve_H(1/2,-1)
+            -sqrt(2)*sqrt(-1/pi)*(cos(1) - 1)
+            sage: struve_H(1/2,pi)
+            2*sqrt(2)/pi
+            sage: struve_H(2,x)
+            struve_H(2, x)
+            sage: struve_H(1/2.,pi)
+            0.900316316157106
+            sage: struve_H(1/2,pi).n(200)
+            0.9003163161571060695551991910...
+        """
+        if a.is_zero() and z.is_zero():
+            return 0
+        if a == -Integer(1)/2:
+            from sage.functions.trig import sin
+            return sqrt(2/(pi*z)) * sin(z)
+        if a == Integer(1)/2:
+            from sage.functions.trig import cos
+            return sqrt(2/(pi*z)) * (1-cos(z))
+        if isinstance(a, Expression) or isinstance(z, Expression):
+            return None
+        if is_inexact(a) or is_inexact(z):
+            a, z = get_coercion_model().canonical_coercion(a, z)
+            return self._evalf_(a, z, parent(z))
+        return None
+    
+    def _evalf_(self, a, z, parent=None, algorithm=None):
+        """
+        EXAMPLES::
+        
+            sage: struve_H(1/2.,pi)
+            0.900316316157106
+        """
+        import mpmath
+        return mpmath_utils.call(mpmath.struveh, a, z, parent=parent)
+    
+    def _derivative_(self, a, z, diff_param=None):
+        """
+        EXAMPLES::
+        
+        sage: diff(struve_H(3/2,x),x)
+        -1/2*sqrt(2)*sqrt(1/(pi*x))*(cos(x) - 1) + 1/16*sqrt(2)*x^(3/2)/sqrt(pi) - 1/2*struve_H(5/2, x)
+        """
+        if diff_param == 0:
+            raise ValueError("cannot differentiate struve_H in the first parameter")
+
+        from sage.functions.other import sqrt, gamma
+        return (z**a/(sqrt(pi)*2**a*gamma(a+Integer(3)/Integer(2)))-struve_H(a+1,z)+struve_H(a-1,z))/2
+    
+    def _print_latex_(self, a, z):
+        """
+            sage: latex(struve_H(2,x))
+            H_{{2}}({x})
+        """
+        return r"H_{{%s}}({%s})" % (a, z)
+    
+struve_H = Function_Struve_H()
+
+class Function_Struve_L(BuiltinFunction):
+    r"""
+    The modified Struve functions.
+
+    .. math::
+    
+        \mathrm{L}_\alpha(x) = -i\cdot e^{-\frac{i\alpha\pi}{2}}\cdot\mathrm{H}_\alpha(ix)
+    
+    EXAMPLES::
+    
+        sage: struve_L(2,x)
+        struve_L(2, x)
+        sage: struve_L(1/2.,pi)
+        4.76805417696286
+        sage: diff(struve_L(1,x),x)
+        1/3*x/pi - 1/2*struve_L(2, x) + 1/2*struve_L(0, x)
+        
+    REFERENCES:
+    
+        - Abramowitz and Stegun: Handbook of Mathematical Functions,
+          http://www.math.sfu.ca/~cbm/aands/
+    
+        - http://en.wikipedia.org/wiki/Struve_function
+    """
+    def __init__(self):
+        r"""
+    ^    EXAMPLES::
+    
+            sage: maxima("struve_l(n,x);").sage()
+            struve_L(n, x)
+            sage: struve_L(7/5,1)._maxima_()
+            struve_l(7/5,1)
+        """
+        BuiltinFunction.__init__(self, 'struve_L', nargs=2,
+                                 conversions=dict(maple='StruveL',
+                                                  mathematica='StruveL',
+                                                  maxima='struve_l',
+                                                  sympy='struvel'))
+
+    def _eval_(self, a, z):
+        """
+        EXAMPLES::
+        
+            sage: struve_L(0,0)
+            0
+            sage: struve_L(2,x)
+            struve_L(2, x)
+            sage: struve_L(1/2.,pi)
+            4.76805417696286
+            sage: struve_L(1/2,pi).n(200)
+            4.768054176962864289162484345...
+        """
+        if a.is_zero() and z.is_zero():
+            return 0
+        if isinstance(a, Expression) or isinstance(z, Expression):
+            return None
+        if is_inexact(a) or is_inexact(z):
+            a, z = get_coercion_model().canonical_coercion(a, z)
+            return self._evalf_(a, z, parent(z))
+        return None
+    
+    def _evalf_(self, a, z, parent=None, algorithm=None):
+        """
+        EXAMPLES::
+        
+            sage: struve_L(1/2.,pi)
+            4.76805417696286
+        """
+        import mpmath
+        return mpmath_utils.call(mpmath.struvel, a, z, parent=parent)
+    
+    def _derivative_(self, a, z, diff_param=None):
+        """
+        EXAMPLES::
+        
+            sage: diff(struve_L(1,x),x)
+            1/3*x/pi - 1/2*struve_L(2, x) + 1/2*struve_L(0, x)
+        """
+        if diff_param == 0:
+            raise ValueError("cannot differentiate struve_L in the first parameter")
+
+        from sage.functions.other import sqrt, gamma
+        return (z**a/(sqrt(pi)*2**a*gamma(a+Integer(3)/Integer(2)))-struve_L(a+1,z)+struve_L(a-1,z))/2
+    
+    def _print_latex_(self, a, z):
+        """
+            sage: latex(struve_L(2,x))
+            L_{{2}}({x})
+        """
+        return r"L_{{%s}}({%s})" % (a, z)
+    
+struve_L = Function_Struve_L()