Ticket #7748: trac_7748-incomplete_gamma_ver2.4.3.3.alpha1.patch

File trac_7748-incomplete_gamma_ver2.4.3.3.alpha1.patch, 7.2 KB (added by burcin, 12 years ago)

rebased to 4.3.3.alpha1

  • sage/functions/all.py 

    # HG changeset patch
# User Flavia Stan <flavia.stan@gmail.com>
# Date 1266276321 -3600
# Node ID a0ce93292e40d05415d90614464b257b5cce98ec
# Parent  a2f3ef9bff74d59da1629398f19da126d7beb5a9
Trac 7748: Make incomplete gamma function symbolic.

diff --git a/sage/functions/all.py b/sage/functions/all.py
    a b  
    1717
    1818from other import ( ceil, floor, gamma, psi, factorial,
    1919                    abs_symbolic, erf, sqrt,
     20                    gamma_inc, incomplete_gamma,
    2021                    real_part, real,
    2122                    imag_part, imag, imaginary, conjugate)
    2223
    … …  
    2425
    2526
    2627from transcendental import (exponential_integral_1,
    27                             gamma_inc, incomplete_gamma,
    2828                            zeta, zeta_symmetric,
    2929                            Li, Ei,
    3030                            dickman_rho)

  • sage/functions/other.py 

    diff --git a/sage/functions/other.py b/sage/functions/other.py
    a b  
    55from sage.symbolic.expression import Expression
    66from sage.libs.pari.gen import pari
    77from sage.symbolic.all import SR
    8 from sage.rings.all import Integer, Rational, RealField, CC, RR
     8from sage.rings.all import Integer, Rational, RealField, CC, RR, \
     9     is_ComplexNumber, ComplexField
    910from sage.misc.latex import latex
    1011import math
    1112
     13import sage.structure.element
     14coercion_model = sage.structure.element.get_coercion_model()
     15
     16from sage.structure.coerce import parent
     17from sage.symbolic.constants import pi
     18from sage.symbolic.function import is_inexact
     19from sage.functions.log import exp
     20from sage.functions.transcendental import Ei
     21
    1222one_half = ~SR(2)
    1323
    1424class Function_erf(BuiltinFunction):
    … …  
    539549
    540550gamma = Function_gamma()
    541551
     552class Function_gamma_inc(BuiltinFunction):
     553    def __init__(self):
     554        """
     555        The incomplete gamma function.
     556
     557        EXAMPLES::
     558
     559            sage: gamma_inc(CDF(0,1), 3)
     560            0.00320857499337 + 0.0124061858119*I
     561            sage: gamma_inc(RDF(1), 3)
     562            0.0497870683678639
     563            sage: gamma_inc(3,2)
     564            gamma(3, 2)
     565            sage: gamma_inc(x,0)
     566            gamma(x)
     567            sage: latex(gamma_inc(3,2))
     568            \Gamma\left(3, 2\right)
     569            sage: loads(dumps((gamma_inc(3,2))))
     570            gamma(3, 2)
     571            sage: i = ComplexField(30).0; gamma_inc(2, 1 + i)
     572            0.70709210 - 0.42035364*I
     573            sage: gamma_inc(2., 5)
     574            0.0404276819945128
     575        """
     576        BuiltinFunction.__init__(self, "gamma", nargs=2, latex_name=r"\Gamma")
     577
     578    def _eval_(self, x, y):
     579        """
     580        EXAMPLES::
     581
     582            sage: gamma_inc(2.,0)
     583            1.00000000000000
     584            sage: gamma_inc(2,0)
     585            1
     586            sage: gamma_inc(1/2,2)
     587            -(erf(sqrt(2)) - 1)*sqrt(pi)
     588            sage: gamma_inc(1/2,1)
     589            -(erf(1) - 1)*sqrt(pi)
     590            sage: gamma_inc(1/2,0)
     591            sqrt(pi)
     592            sage: gamma_inc(x,0)
     593            gamma(x)
     594            sage: gamma_inc(1,2)
     595            e^(-2)
     596            sage: gamma_inc(0,2)
     597            -Ei(-2)
     598        """
     599        if not isinstance(x, Expression) and not isinstance(y, Expression) and \
     600               (is_inexact(x) or is_inexact(y)):
     601            x, y = coercion_model.canonical_coercion(x, y)
     602            return self._evalf_(x, y, parent(x))
     603
     604        if y == 0:
     605            return gamma(x)
     606        if x == 1:
     607            return exp(-y)
     608        if x == 0:
     609            return -Ei(-y)
     610        if x == Rational(1)/2: #only for x>0
     611            return sqrt(pi)*(1-erf(sqrt(y)))
     612        return None
     613
     614    def _evalf_(self, x, y, parent=None):
     615        """
     616        EXAMPLES::
     617
     618            sage: gamma_inc(0,2)
     619            -Ei(-2)
     620            sage: gamma_inc(0,2.)
     621            0.0489005107080611
     622            sage: gamma_inc(3,2).n()
     623            1.35335283236613
     624        """
     625        try:
     626            return x.gamma_inc(y)
     627        except AttributeError:
     628            if not (is_ComplexNumber(x)):
     629                if is_ComplexNumber(y):
     630                    C = y.parent()
     631                else:
     632                    C = ComplexField()
     633                    x = C(x)
     634            return x.gamma_inc(y)
     635           
     636# synonym.
     637incomplete_gamma = gamma_inc=Function_gamma_inc()
     638
     639
    542640class Function_psi1(GinacFunction):
    543641    def __init__(self):
    544642        r"""

  • sage/functions/transcendental.py 

    diff --git a/sage/functions/transcendental.py b/sage/functions/transcendental.py
    a b  
    101101    else:
    102102        return [float(z) for z in pari(x).eint1(n)]
    103103
    104 def gamma_inc(s, t):
    105     """
    106     Incomplete Gamma function Gamma(s,t).
    107    
    108     EXAMPLES::
    109    
    110         sage: gamma_inc(CDF(0,1), 3)
    111         0.00320857499337 + 0.0124061858119*I
    112         sage: gamma_inc(3, 3)
    113         0.846380162253687
    114         sage: gamma_inc(RDF(1), 3)
    115         0.0497870683678639
    116     """
    117     try:
    118         return s.gamma_inc(t)
    119     except AttributeError:
    120         if not (is_ComplexNumber(s)):
    121             if is_ComplexNumber(t):
    122                 C = t.parent()
    123             else:
    124                 C = ComplexField()
    125             s = C(s)
    126         return s.gamma_inc(t)
    127 
    128 
    129 # synonym.
    130 incomplete_gamma = gamma_inc
    131104
    132105def zeta(s):
    133106    """

  • sage/rings/complex_number.pyx 

    diff --git a/sage/rings/complex_number.pyx b/sage/rings/complex_number.pyx
    a b  
    18301830            -0.0013781309 + 0.0065198200*I
    18311831            sage: C(2).gamma_inc(1 + i)
    18321832            0.70709210 - 0.42035364*I
    1833             sage: gamma_inc(2, 1 + i)   
    1834             0.70709210 - 0.42035364*I
    1835             sage: gamma_inc(2, 5)
     1833            sage: CC(2).gamma_inc(5)
    18361834            0.0404276819945128
    18371835        """
    18381836        return self._parent(self._pari_().incgam(t))

  • sage/symbolic/random_tests.py 

    diff --git a/sage/symbolic/random_tests.py b/sage/symbolic/random_tests.py
    a b  
    202202
    203203        sage: from sage.symbolic.random_tests import *
    204204        sage: random_expr(50, nvars=3, coeff_generator=CDF.random_element)
    205         sinh(sinh(-coth(v2)/csc((-0.615863165633 + 0.879368031485*I)*v1^2*v3) - gamma(pi) + erf((-0.708874026302 - 0.954135400334*I)*v3)))^coth(-cosh(-polylog((v2 + 0.913564344312 + 0.0898040160336*I)^((-0.723896589334 - 0.799038508886*I)*v2), -v1 - v3))/arcsin((-0.0263902659909 + 0.153261789843*I)*cosh(pi*catalan)))
     205        sinh(sinh((-0.314177274493 + 0.144437996366*I)/csc(-v1^2*e/v3) + erf((-0.708874026302 - 0.954135400334*I)*v3) + 0.0275857401668 - 0.479027260657*I))^(-cosh(-polylog((v2^2 + (0.067987275089 + 1.08529153495*I)*v3)^(-v1 - v3), (5.29385548262 + 2.57440711353*I)*e/cosh((-0.436810529675 + 0.736945423566*I)*arccot(pi)))))
    206206        sage: random_expr(5, verbose=True)
    207207        About to apply <built-in function add> to [v1, v1]
    208         About to apply <built-in function div> to [-1/3, 2*v1]
    209         -1/6/v1
     208        About to apply <built-in function mul> to [v1, -1/2]
     209        About to apply <built-in function mul> to [2*v1, -1/2*v1]
     210        -v1^2
    210211    """
    211212    vars = [(1.0, sage.calculus.calculus.var('v%d' % (n+1))) for n in range(nvars)]
    212213    if ncoeffs is None: