Ticket #9706: trac_9706_chebyshev.patch

File trac_9706_chebyshev.patch, 25.5 KB (added by jdemeyer, 7 years ago)

cheby changes

  • sage/functions/orthogonal_polys.py

    # HG changeset patch
    # User Stefan Reiterer <domors@gmx.net>
    # Date 1384035568 -3600
    # Node ID 2adf09336035d579faf00ead8ebadc209f79065c
    # Parent  f3185fb2a16984a691df7be3d086d79286c58854
    Implement symbolic Chebyshev polynomials
    
    diff --git a/sage/functions/orthogonal_polys.py b/sage/functions/orthogonal_polys.py
    a b  
    261261   (it loads "specfun", but I'm not sure if that is the real reason).
    262262   The next call is usually faster but not always.
    263263
    264 TODO: Implement associated Legendre polynomials and Zernike
    265 polynomials. (Neither is in Maxima.)
    266 http://en.wikipedia.org/wiki/Associated_Legendre_polynomials
    267 http://en.wikipedia.org/wiki/Zernike_polynomials
     264.. TODO::
     265
     266    Implement associated Legendre polynomials and Zernike
     267    polynomials. (Neither is in Maxima.)
     268    :wikipedia:`Associated_Legendre_polynomials`
     269    :wikipedia:`Zernike_polynomials`
    268270
    269271REFERENCES:
    270272
    271 Abramowitz and Stegun: Handbook of Mathematical Functions,
     273[ASHandbook] Abramowitz and Stegun: Handbook of Mathematical Functions,
    272274   http://www.math.sfu.ca/ cbm/aands/
    273275
    274 http://en.wikipedia.org/wiki/Chebyshev_polynomials
     276:wikipedia:`Chebyshev_polynomials`
    275277
    276 http://en.wikipedia.org/wiki/Legendre_polynomials
     278:wikipedia:`Legendre_polynomials`
    277279
    278 http://en.wikipedia.org/wiki/Hermite_polynomials
     280:wikipedia:`Hermite_polynomials`
    279281
    280282-  http://mathworld.wolfram.com/GegenbauerPolynomial.html
    281283
    282 http://en.wikipedia.org/wiki/Jacobi_polynomials
     284:wikipedia:`Jacobi_polynomials`
    283285
    284 http://en.wikipedia.org/wiki/Laguerre_polynomia
     286:wikipedia:`Laguerre_polynomia`
    285287
    286 http://en.wikipedia.org/wiki/Associated_Legendre_polynomials
     288:wikipedia:`Associated_Legendre_polynomials`
    287289
     290- [EffCheby] Wolfram Koepf: Effcient Computation of Chebyshev Polynomials
     291  in Computer Algebra
     292  Computer Algebra Systems: A Practical Guide.
     293  John Wiley, Chichester (1999): 79-99.
    288294
    289295AUTHORS:
    290296
    291297- David Joyner (2006-06)
     298- Stefan Reiterer (2010-)
    292299"""
    293300
    294301#*****************************************************************************
    295302#       Copyright (C) 2006 William Stein <wstein@gmail.com>
    296303#                     2006 David Joyner <wdj@usna.edu>
     304#                     2010 Stefan Reiterer <maldun.finsterschreck@gmail.com>
    297305#
    298306#  Distributed under the terms of the GNU General Public License (GPL)
    299307#
     
    311319from sage.rings.all import ZZ
    312320from sage.calculus.calculus import maxima
    313321
     322
     323from sage.symbolic.function import BuiltinFunction, GinacFunction, is_inexact
     324from sage.symbolic.expression import is_Expression
     325import sage.functions.special
     326from sage.functions.special import MaximaFunction, meval
     327from sage.functions.other import floor, gamma, factorial, abs, binomial
     328from sage.functions.other import sqrt, conjugate
     329from sage.functions.trig import sin, cos
     330from sage.functions.log import ln
     331import sage.symbolic.expression as expression
     332from sage.structure.parent import Parent
     333from sage.structure.coerce import parent
     334
    314335_done = False
    315336def _init():
    316337    """
     
    327348
    328349    Then after using one of these functions, it changes::
    329350
    330         sage: from sage.functions.orthogonal_polys import chebyshev_T
    331         sage: chebyshev_T(2,x)
    332         2*(x - 1)^2 + 4*x - 3
     351        sage: from sage.functions.orthogonal_polys import legendre_P
     352        sage: legendre_P(2,x)
     353        3/2*(x - 1)^2 + 3*x - 2
    333354        sage: sage.functions.orthogonal_polys._done
    334355        True
    335356
     357
    336358    Note that because here we use a Pynac variable ``x``,
    337359    the representation of the function is different from
    338360    its actual doctest, where a polynomial indeterminate
     
    348370    _done = True
    349371
    350372
    351 def chebyshev_T(n,x):
     373
     374class OrthogonalPolynomial(BuiltinFunction):
    352375    """
    353     Returns the Chebyshev function of the first kind for integers
    354     `n>-1`.
     376    Base class for orthogonal polynomials.
     377
     378    This class is an abstract base class for all orthogonal polynomials since
     379    they share similar properties. The evaluation as a polynomial
     380    is either done via maxima, or with pynac.
     381
     382    Convention: The first argument is always the order of the polynomial,
     383    the last one is always the value `x` where the polynomial is evaluated.
     384    """
     385    def __init__(self, name, nargs=2, latex_name=None, conversions={}):
     386        """
     387        :class:`OrthogonalPolynomial` class needs the same input parameter as
     388        it's parent class.
     389       
     390        EXAMPLES::
     391
     392            sage: from sage.functions.orthogonal_polys import OrthogonalPolynomial
     393            sage: new = OrthogonalPolynomial('testo_P')
     394            sage: new
     395            testo_P
     396        """
     397        try:
     398            self._maxima_name = conversions['maxima']
     399        except KeyError:
     400            self._maxima_name = None
     401
     402        BuiltinFunction.__init__(self, name=name, nargs=nargs,
     403                                 latex_name=latex_name, conversions=conversions)
     404
     405    def _maxima_init_evaled_(self, *args):
     406        r"""
     407        Return a string which represents this function evaluated at
     408        ``*args`` in Maxima.
     409
     410        In fact these are thought to be the old wrappers for the orthogonal
     411        polynomials. They are used when the other evaluation methods fail,
     412        or are not fast enough. It appears that direct computation
     413        with pynac is in most cases faster than maxima. Maxima comes into
     414        play when all other methods fail.
     415
     416        EXAMPLES::
     417
     418            sage: chebyshev_T(3,x)
     419            4*x^3 - 3*x
     420        """
     421        return None
     422
     423    def _clenshaw_method_(self, *args):
     424        """
     425        The Clenshaw method uses the three term recursion of the polynomial,
     426        or explicit formulas instead of maxima to evaluate the polynomial
     427        efficiently, if the `x` argument is not a symbolic expression.
     428        The name comes from the Clenshaw algorithm for fast evaluation of
     429        polynomials in chebyshev form.
     430
     431        EXAMPLES::
     432
     433            sage: from sage.functions.orthogonal_polys import OrthogonalPolynomial
     434            sage: new = OrthogonalPolynomial('testo_P')
     435            sage: new._clenshaw_method_(1,2.0)
     436            Traceback (most recent call last):
     437            ...
     438            NotImplementedError: no recursive calculation of values implemented
     439        """
     440        raise NotImplementedError("no recursive calculation of values implemented")
     441
     442    def _eval_special_values_(self,*args):
     443        """
     444        Evaluate the polynomial explicitly for special values.
     445
     446        EXAMPLES::
     447           
     448            sage: var('n')
     449            n
     450            sage: chebyshev_T(n,-1)
     451            (-1)^n
     452        """
     453        raise ValueError("no special values known")
     454
     455    def _eval_(self, *args):
     456        """
     457        The symbolic evaluation is either done with maxima, or with direct
     458        computaion in pynac.
     459
     460        For the fast numerical evaluation, another method should be used.
     461        We use Clenshaw's algorithm, which uses the recursion.
     462        The function also checks for special values, and if
     463        the order is an integer and in range.
     464       
     465        EXAMPLES::
     466
     467            sage: chebyshev_T(5,x)
     468            16*x^5 - 20*x^3 + 5*x 
     469            sage: var('n')
     470            n
     471            sage: chebyshev_T(n,-1)
     472            (-1)^n
     473            sage: chebyshev_T(-7,x)
     474            64*x^7 - 112*x^5 + 56*x^3 - 7*x
     475            sage: chebyshev_T(3/2,x)
     476            chebyshev_T(3/2, x)
     477        """
     478        if not is_Expression(args[0]):
     479            if not is_Expression(args[-1]) and is_inexact(args[-1]):
     480                try:
     481                    import sage.libs.mpmath.all as mpmath
     482                    return self._evalf_(*args)
     483                except AttributeError:
     484                    pass
     485                except mpmath.NoConvergence:
     486                    print "Warning: mpmath returns NoConvergence!"
     487                    print "Switching to clenshaw_method, but it may not be stable!"
     488                except ValueError:
     489                    pass
     490
     491            # A faster check would be nice...
     492            if args[0] != floor(args[0]):
     493                if not is_Expression(args[-1]):
     494                    try:
     495                        return self._evalf_(*args)
     496                    except AttributeError:
     497                        pass
     498                else:
     499                    return None
     500        try:
     501            return self._eval_special_values_(*args)
     502        except ValueError:
     503            pass
     504
     505        #if negative indices are not specified
     506        #in _eval_special_values only return symbolic
     507        #value
     508        if args[0] < 0 and args[0] in ZZ:
     509                return None
     510
     511        if not is_Expression(args[0]):
     512            try:
     513                return self._clenshaw_method_(*args)
     514            except NotImplementedError:
     515                pass
     516
     517        if self._maxima_name is None:
     518            return None
     519
     520        _init()
     521        try: # use maxima as last resort
     522            s = maxima(self._maxima_init_evaled_(*args))
     523            return self._maxima_init_evaled_(*args)
     524        except TypeError:
     525            return None
     526       
     527        if self._maxima_name in repr(s):
     528            return None
     529
     530        return s.sage()
     531
     532    def __call__(self,n,x):
     533        """
     534        This overides the call method from SageObject to avoid problems with coercions,
     535        since the _eval_ method is able to handle for more data types than symbolic functions
     536        would normally allow.
     537        Thus we have the distinction between algebraic objects (if n is an integer),
     538        and else as symbolic function.
     539       
     540        EXAMPLES::
     541       
     542            sage: K.<a> = NumberField(x^3-x-1)
     543            sage: chebyshev_T(5, a)
     544            16*a^2 + a - 4
     545            sage: chebyshev_T(5,MatrixSpace(ZZ, 2)([1, 2, -4, 7]))
     546            [-40799  44162]
     547            [-88324  91687]
     548            sage: R.<x> = QQ[]
     549            sage: parent(chebyshev_T(5, x))
     550            Univariate Polynomial Ring in x over Rational Field
     551        """
     552       
     553        if n in ZZ: # check if n is integer -> consider polynomial as algebraic structure
     554            return self._eval_(n,x) # Let eval methode decide which is best
     555        else: # Consider OrthogonalPolynomial as symbol
     556            return super(OrthogonalPolynomial,self).__call__(n,x)
     557
     558
     559class Func_chebyshev_T(OrthogonalPolynomial):
     560    """
     561    Chebyshev polynomials of the first kind.
     562
     563    REFERENCE:
     564
     565    - [ASHandbook]_ 22.5.31 page 778 and 6.1.22 page 256.
     566
     567    EXAMPLES::
     568
     569       sage: chebyshev_T(3,x)
     570       4*x^3 - 3*x
     571       sage: chebyshev_T(5,x)
     572       16*x^5 - 20*x^3 + 5*x
     573       sage: var('k')
     574       k
     575       sage: test = chebyshev_T(k,x)
     576       sage: test
     577       chebyshev_T(k, x)
     578    """
     579    def __init__(self):
     580        """
     581        Init method for the chebyshev polynomials of the first kind.
     582
     583        EXAMPLES::
     584
     585            sage: from sage.functions.orthogonal_polys import Func_chebyshev_T
     586            sage: chebyshev_T2 = Func_chebyshev_T()
     587            sage: chebyshev_T2(1,x)
     588            x
     589        """
     590        OrthogonalPolynomial.__init__(self, "chebyshev_T", nargs=2,
     591                                      conversions=dict(maxima='chebyshev_t',
     592                                                       mathematica='ChebyshevT'))
     593   
     594    def _eval_special_values_(self,*args):
     595        """
     596        Values known for special values of x.
     597        For details see [ASHandbook]_ 22.4 (p. 777)
     598
     599        EXAMPLES:
     600
     601            sage: var('n')
     602            n
     603            sage: chebyshev_T(n,1)
     604            1
     605            sage: chebyshev_T(n,-1)
     606            (-1)^n
     607        """
     608        if args[-1] == 1:
     609            return 1
     610       
     611        if args[-1] == -1:
     612            return (-1)**args[0]
     613
     614        if (args[-1] == 0):
     615            return (1+(-1)**args[0])*(-1)**(args[0]/2)/2
     616
     617        if args[0] < 0 and args[0] in ZZ:
     618            return self._eval_(-args[0],args[-1])
     619
     620        raise ValueError("value not found")
     621   
     622    def _evalf_(self, *args,**kwds):
     623        """
     624        Evaluates :class:`chebyshev_T` numerically with mpmath.
     625
     626        EXAMPLES::
     627
     628            sage: chebyshev_T(10,3).n(75)
     629            2.261953700000000000000e7
     630            sage: chebyshev_T(10,I).n()
     631            -3363.00000000000
     632            sage: chebyshev_T(5,0.3).n()
     633            0.998880000000000
     634        """
     635        try:
     636            step_parent = kwds['parent']
     637        except KeyError:
     638            step_parent = parent(args[-1])
     639
     640        try:
     641            precision = step_parent.prec()
     642        except AttributeError:
     643            precision = RR.prec()
     644
     645        from sage.libs.mpmath.all import call as mpcall
     646        from sage.libs.mpmath.all import chebyt as mpchebyt
     647
     648        return mpcall(mpchebyt,args[0],args[-1],prec = precision)
     649
     650    def _maxima_init_evaled_(self, *args):
     651        """
     652        Evaluate the Chebyshev polynomial ``self`` with maxima.
     653
     654        EXAMPLES::
     655
     656            sage: chebyshev_T._maxima_init_evaled_(1,x)
     657            x
     658        """
     659        n = args[0]
     660        x = args[1]
     661        return sage_eval(maxima.eval('chebyshev_t(%s,x)'%ZZ(n)), locals={'x':x})
     662
     663    def _clenshaw_method_(self,*args):
     664        """
     665        Clenshaw method for :class:`chebyshev_T` (means use recursions in this
     666        case). This is much faster for numerical evaluation than maxima!
     667        See [ASHandbook]_ 227 (p. 782) for details for the recurions.
     668        See also [EffCheby] for fast evaluation techniques.
     669
     670        EXAMPLES::
     671
     672            sage: chebyshev_T._clenshaw_method_(2,0.1) == chebyshev_T._evalf_(2,0.1)
     673            True
     674            sage: chebyshev_T(51,x)
     675            2*(2*(2*(2*(2*(2*x^2 - 1)^2 - 1)*(2*(2*x^2 - 1)*x - x) - x)*(2*(2*(2*x^2 - 1)*x - x)^2 - 1) - x)^2 - 1)*(2*(2*(2*(2*(2*x^2 - 1)^2 - 1)*(2*(2*x^2 - 1)*x - x) - x)*(2*(2*(2*x^2 - 1)*x - x)^2 - 1) - x)*(2*(2*(2*(2*x^2 - 1)*x - x)^2 - 1)^2 - 1) - x) - x
     676            sage: chebyshev_T._clenshaw_method_(10,x)
     677            512*x^10 - 1280*x^8 + 1120*x^6 - 400*x^4 + 50*x^2 - 1
     678
     679        """
     680        k = args[0]
     681        x = args[1]
     682
     683        if k == 0:
     684            return 1
     685        if k == 1:
     686            return x
     687
     688        # TODO: When evaluation of Symbolic Expressions works better
     689        # use these explicit formulas instead!
     690        #if -1 <= x <= 1:
     691        #    return cos(k*acos(x))
     692        #elif 1 < x:
     693        #    return cosh(k*acosh(x))
     694        #else: # x < -1
     695        #    return (-1)**(k%2)*cosh(k*acosh(-x))
     696       
     697        help1 = 1
     698        help2 = x
     699        if is_Expression(x) and k <= 50:
     700            # Recursion gives more compact representations for large k
     701            help3 = 0
     702            for j in xrange(0,floor(k/2)+1):
     703                f = factorial(k-j-1) / factorial(j) / factorial(k-2*j)
     704                help3 = help3 + (-1)**j * (2*x)**(k-2*j) * f
     705            help3 = help3 * k / 2
     706            return help3
     707        else:
     708            return self._cheb_recur_(k,x)[0]
     709
     710    def _cheb_recur_(self,n, x, both=False):
     711        """
     712        Generalized recursion formula for Chebyshev polynomials.
     713        Implementation suggested by Frederik Johansson.
     714        returns (T(n,x), T(n-1,x)), or (T(n,x), _) if both=False
     715       
     716        EXAMPLES::
     717       
     718            sage: chebyshev_T._cheb_recur_(5,x)
     719            (2*(2*(2*x^2 - 1)*x - x)*(2*x^2 - 1) - x, False)
     720        """
     721        if n == 0:
     722            return 1, x
     723        if n == 1:
     724            return x, 1
     725        a, b = self._cheb_recur_((n+1)//2, x, both or n % 2)
     726        if n % 2 == 0:
     727            return 2*a**2 - 1, both and 2*a*b - x
     728        else:
     729            return 2*a*b - x, both and 2*b**2 - 1
     730
     731
     732    def _eval_numpy_(self, *args):
     733        """
     734        Evaluate ``self`` using numpy.
     735
     736        EXAMPLES::
     737
     738            sage: import numpy
     739            sage: z = numpy.array([1,2])
     740            sage: z2 = numpy.array([[1,2],[1,2]])
     741            sage: z3 = numpy.array([1,2,3.])
     742            sage: chebyshev_T(1,z)
     743            array([1, 2])
     744            sage: chebyshev_T(1,z2)
     745            array([[1, 2],
     746                   [1, 2]])
     747            sage: chebyshev_T(1,z3)
     748            array([ 1.,  2.,  3.])
     749            sage: chebyshev_T(z,0.1)
     750            array([ 0.1 , -0.98])
     751        """
     752        from scipy.special import eval_chebyt
     753        return eval_chebyt(args[0],args[-1])
     754
     755    def _derivative_(self, *args, **kwds):
     756        """
     757        Return the derivative of :class:`chebyshev_T` in form of the Chebyshev
     758        polynomial of the second kind :class:`chebyshev_U`.
     759
     760        EXAMPLES::
     761
     762            sage: var('k')
     763            k
     764            sage: derivative(chebyshev_T(k,x),x)
     765            k*chebyshev_U(k - 1, x)
     766            sage: derivative(chebyshev_T(3,x),x)
     767            12*x^2 - 3
     768            sage: derivative(chebyshev_T(k,x),k)
     769            Traceback (most recent call last):
     770            ...
     771            NotImplementedError: derivative w.r.t. to the index is not supported yet
     772        """
     773        diff_param = kwds['diff_param']
     774        if diff_param == 0:
     775            raise NotImplementedError("derivative w.r.t. to the index is not supported yet")
     776
     777        return args[0]*chebyshev_U(args[0]-1,args[1])
     778   
     779chebyshev_T = Func_chebyshev_T()
     780
     781class Func_chebyshev_U(OrthogonalPolynomial):
     782    """
     783    Class for the Chebyshev polynomial of the second kind.
    355784   
    356785    REFERENCE:
    357786
    358     - AS 22.5.31 page 778 and AS 6.1.22 page 256.
    359    
    360     EXAMPLES::
    361    
    362         sage: x = PolynomialRing(QQ, 'x').gen()
    363         sage: chebyshev_T(2,x)
    364         2*x^2 - 1
    365     """
    366     _init()
    367     return sage_eval(maxima.eval('chebyshev_t(%s,x)'%ZZ(n)), locals={'x':x})
    368 
    369 def chebyshev_U(n,x):
    370     """
    371     Returns the Chebyshev function of the second kind for integers `n>-1`.
    372    
    373     REFERENCE:
    374 
    375     - AS, 22.8.3 page 783 and AS 6.1.22 page 256.
     787    - [ASHandbook]_ 22.8.3 page 783 and 6.1.22 page 256.
    376788   
    377789    EXAMPLES::
    378790   
    379791        sage: x = PolynomialRing(QQ, 'x').gen()
    380792        sage: chebyshev_U(2,x)
    381793        4*x^2 - 1
     794        sage: chebyshev_U(3,x)
     795        8*x^3 - 4*x
    382796    """
    383     _init()
    384     return sage_eval(maxima.eval('chebyshev_u(%s,x)'%ZZ(n)), locals={'x':x})
     797    def __init__(self):
     798        """
     799        Init method for the chebyshev polynomials of the second kind.
     800
     801        EXAMPLES::
     802       
     803            sage: from sage.functions.orthogonal_polys import Func_chebyshev_U
     804            sage: chebyshev_U2 = Func_chebyshev_U()
     805            sage: chebyshev_U2(1,x)
     806            2*x
     807        """
     808        OrthogonalPolynomial.__init__(self, "chebyshev_U", nargs=2,
     809                                      conversions=dict(maxima='chebyshev_u',
     810                                                       mathematica='ChebyshevU'))
     811       
     812    def _clenshaw_method_(self,*args):
     813        """
     814        Clenshaw method for :class:`chebyshev_U` (means we use the recursion)
     815        This is much faster for numerical evaluation than maxima.
     816        See [ASHandbook]_ 227 (p. 782) for details on the recurions.
     817        See also [
     818
     819        EXAMPLES::
     820
     821            sage: chebyshev_U._clenshaw_method_(2,0.1) == chebyshev_U._evalf_(2,0.1)
     822            True
     823        """
     824        k = args[0]
     825        x = args[1]
     826           
     827        if k == 0:
     828            return 1
     829        if k == 1:
     830            return 2*x
     831
     832        help1 = 1
     833        help2 = 2*x
     834        if is_Expression(x) and k <= 50:
     835            # Recursion gives more compact representations for large k
     836            help3 = 0
     837            for j in xrange(0,floor(k/2)+1):
     838                f = factorial(k-j) / factorial(j) / factorial(k-2*j) # Change to a binomial?
     839                help3 = help3 + (-1)**j * (2*x)**(k-2*j) * f
     840            return help3
     841           
     842        else:
     843            return self._cheb_recur_(k,x)[0]
     844
     845 
     846    def _cheb_recur_(self,n, x, both=False):
     847        """
     848        Generalized recursion formula for Chebyshev polynomials.
     849        Implementation suggested by Frederik Johansson.
     850        returns (U(n,x), U(n-1,x)), or (U(n,x), _) if both=False
     851       
     852        EXAMPLES::
     853       
     854            sage: chebyshev_U._cheb_recur_(3,x)
     855            (4*(2*x^2 - 1)*x, False)
     856            sage: chebyshev_U._cheb_recur_(5,x)[0]
     857            -2*((2*x + 1)*(2*x - 1)*x - 4*(2*x^2 - 1)*x)*(2*x + 1)*(2*x - 1)
     858            sage: abs(pari('polchebyshev(5, 2, 0.1)') - chebyshev_U(5,0.1)) < 1e-10
     859            True
     860
     861
     862        """
     863       
     864        if n == 0:
     865            return 1, both and 2*x
     866        if n == 1:
     867            return 2*x, both and 4*x**2-1
     868           
     869        a, b = self._cheb_recur_((n-1)//2, x, True)
     870        if n % 2 == 0:
     871            return (b+a)*(b-a), both and 2*b*(x*b-a)
     872        else:
     873            return 2*a*(b-x*a), both and (b+a)*(b-a)
     874
     875    def _maxima_init_evaled_(self, *args):
     876        """
     877        Uses maxima to evaluate ``self``.
     878
     879        EXAMPLES::
     880
     881            sage: x = PolynomialRing(QQ, 'x').gen()
     882            sage: chebyshev_T(2,x)
     883            2*x^2 - 1
     884        """
     885        n = args[0]
     886        x = args[1]
     887        return sage_eval(maxima.eval('chebyshev_u(%s,x)'%ZZ(n)), locals={'x':x})
     888
     889    def _evalf_(self, *args,**kwds):
     890        """
     891        Evaluate :class:`chebyshev_U` numerically with mpmath.
     892
     893        EXAMPLES::
     894
     895            sage: chebyshev_U(5,-4+3.*I)
     896            98280.0000000000 - 11310.0000000000*I
     897            sage: chebyshev_U(10,3).n(75)
     898            4.661117900000000000000e7
     899        """
     900        try:
     901            step_parent = kwds['parent']
     902        except KeyError:
     903            step_parent = parent(args[-1])
     904
     905        try:
     906            precision = step_parent.prec()
     907        except AttributeError:
     908            precision = RR.prec()
     909
     910        from sage.libs.mpmath.all import call as mpcall
     911        from sage.libs.mpmath.all import chebyu as mpchebyu
     912
     913        return mpcall(mpchebyu,args[0],args[-1],prec = precision)
     914
     915    def _eval_special_values_(self,*args):
     916        """
     917        Special values that known. [ASHandbook]_ 22.4 (p.777).
     918
     919        EXAMPLES::
     920       
     921            sage: var('n')
     922            n
     923            sage: chebyshev_U(n,1)
     924            n + 1
     925            sage: chebyshev_U(n,-1)
     926            (-1)^n*(n + 1)
     927        """
     928        if args[-1] == 1:
     929            return (args[0]+1)
     930       
     931        if args[-1] == -1:
     932            return (-1)**args[0]*(args[0]+1)
     933
     934        if (args[-1] == 0):
     935            return (1+(-1)**args[0])*(-1)**(args[0]/2)/2
     936
     937        if args[0] < 0 and args[0] in ZZ:
     938            return -self._eval_(-args[0]-2,args[-1])
     939
     940        raise ValueError("value not found")
     941
     942    def _eval_numpy_(self, *args):
     943        """
     944        Evaluate ``self`` using numpy.
     945
     946        EXAMPLES::
     947
     948            sage: import numpy
     949            sage: z = numpy.array([1,2])
     950            sage: z2 = numpy.array([[1,2],[1,2]])
     951            sage: z3 = numpy.array([1,2,3.])
     952            sage: chebyshev_U(1,z)
     953            array([2, 4])
     954            sage: chebyshev_U(1,z2)
     955            array([[2, 4],
     956                   [2, 4]])
     957            sage: chebyshev_U(1,z3)
     958            array([ 2.,  4.,  6.])
     959            sage: chebyshev_U(z,0.1)
     960            array([ 0.2 , -0.96])
     961        """
     962        from scipy.special import eval_chebyu
     963        return eval_chebyu(args[0],args[1])
     964
     965
     966    def _derivative_(self, *args, **kwds):
     967        """
     968        Return the derivative of :class:`chebyshev_U` in form of the Chebyshev
     969        polynomials of the first and second kind.
     970       
     971        EXAMPLES::
     972
     973            sage: var('k')
     974            k
     975            sage: derivative(chebyshev_U(k,x),x)
     976            ((k + 1)*chebyshev_T(k + 1, x) - x*chebyshev_U(k, x))/(x^2 - 1)
     977            sage: derivative(chebyshev_U(3,x),x)
     978            24*x^2 - 4
     979            sage: derivative(chebyshev_U(k,x),k)
     980            Traceback (most recent call last):
     981            ...
     982            NotImplementedError: derivative w.r.t. to the index is not supported yet
     983        """
     984        diff_param = kwds['diff_param']
     985        if diff_param == 0:
     986            raise NotImplementedError("derivative w.r.t. to the index is not supported yet")
     987
     988        return ((args[0]+1)*chebyshev_T(args[0]+1,args[1])-args[1]*
     989                chebyshev_U(args[0],args[1]))/(args[1]**2-1)
     990
     991chebyshev_U = Func_chebyshev_U()
     992
    385993
    386994def gen_laguerre(n,a,x):
    387995    """
    388996    Returns the generalized Laguerre polynomial for integers `n > -1`.
    389     Typically, a = 1/2 or a = -1/2.
     997    Typically, `a = 1/2` or `a = -1/2`.
    390998   
    391999    REFERENCE:
    3921000
    393     - table on page 789 in AS.
     1001    - Table on page 789 in [ASHandbook]_.
    3941002   
    3951003    EXAMPLES::
    3961004   
     
    4931101   
    4941102    REFERENCE:
    4951103
    496     - AS 22.5.40 and 22.5.41, page 779.
     1104    - [ASHandbook]_ 22.5.40 and 22.5.41, page 779.
    4971105   
    4981106    EXAMPLES::
    4991107   
     
    5281136   
    5291137    REFERENCE:
    5301138
    531     - table on page 789 in AS.
     1139    - Table on page 789 in [ASHandbook]_.
    5321140   
    5331141    EXAMPLES::
    5341142   
     
    5431151
    5441152def laguerre(n,x):
    5451153    """
    546     Returns the Laguerre polynomial for integers `n > -1`.
     1154    Return the Laguerre polynomial for integers `n > -1`.
    5471155   
    5481156    REFERENCE:
    5491157
    550     - AS 22.5.16, page 778 and AS page 789.
     1158    - [ASHandbook]_ 22.5.16, page 778 and page 789.
    5511159   
    5521160    EXAMPLES::
    5531161   
     
    5691177   
    5701178    REFERENCE:
    5711179
    572     - AS 22.5.35 page 779.
     1180    - [ASHandbook]_ 22.5.35 page 779.
    5731181   
    5741182    EXAMPLES::
    5751183   
     
    5921200def legendre_Q(n,x):
    5931201    """
    5941202    Returns the Legendre function of the second kind for integers
    595     `n>-1`.
     1203    `n > -1`.
    5961204   
    5971205    Computed using Maxima.
    5981206   
     
    6201228   
    6211229    REFERENCE:
    6221230
    623     - AS 22.5.27
     1231    - [ASHandbook]_ 22.5.27
    6241232   
    6251233    EXAMPLES::
    6261234