Ticket #11948: 11948.patch

File 11948.patch, 6.1 KB (added by jdemeyer, 7 years ago)
  • sage/functions/other.py

    # HG changeset patch
    # User Jeroen Demeyer <jdemeyer@cage.ugent.be>
    # Date 1319474587 -7200
    # Node ID 5da8e8e0d10e0abc38bcd88f5099044d5370cdfb
    # Parent  bda8f14b0bf36072c1a064b2acdf17a94c78d1ee
    Fix numeric evaluation of error function
    
    diff --git a/sage/functions/other.py b/sage/functions/other.py
    a b  
    2828    _eval_ = BuiltinFunction._eval_default
    2929    def __init__(self):
    3030        r"""
    31         The error function, defined as
     31        The error function, defined for real values as
    3232        `\text{erf}(x) = \frac{2}{\sqrt{\pi}} \int_0^x e^{-t^2} dt`.
     33        This function is also defined for complex values, via analytic
     34        continuation.
    3335
    34         Sage currently only implements the error function (via a call to
    35         PARI) when the input is real.
     36        Sage implements the error function via the ``erfc()`` function in PARI.
    3637
    3738        EXAMPLES::
    3839
     
    4344            sage: loads(dumps(erf))
    4445            erf
    4546
    46         The following fails because we haven't implemented
    47         erf yet for complex values::
    48        
    49             sage: complex(erf(3*I))
    50             Traceback (most recent call last):
    51             ...
    52             TypeError: unable to simplify to complex approximation
    53 
    5447        TESTS:
    5548
    5649        Check if conversion from maxima elements work::
     
    6154        """
    6255        BuiltinFunction.__init__(self, "erf", latex_name=r"\text{erf}")
    6356
    64     def _evalf_(self, x, parent=None):
     57    def _evalf_(self, x, parent):
    6558        """
    6659        EXAMPLES::
    6760
    6861            sage: erf(2).n()
    6962            0.995322265018953
    70             sage: erf(2).n(150)
    71             Traceback (most recent call last):
    72             ...
    73             NotImplementedError: erf not implemented for precision higher than 53
     63            sage: erf(2).n(200)
     64            0.99532226501895273416206925636725292861089179704006007673835
     65            sage: erf(pi - 1/2*I).n(100)
     66            1.0000111669099367825726058952 + 1.6332655417638522934072124548e-6*I
     67
     68        TESTS:
     69
     70        Check that PARI/GP through the GP interface gives the same answer::
     71
     72            sage: gp.set_real_precision(59)  # random
     73            38
     74            sage: print gp.eval("1 - erfc(1)"); print erf(1).n(200);
     75            0.84270079294971486934122063508260925929606699796630290845994
     76            0.84270079294971486934122063508260925929606699796630290845994
    7477        """
    7578        try:
    7679            prec = parent.prec()
    7780        except AttributeError: # not a Sage parent
    7881            prec = 0
    79         if prec > 53:
    80             raise NotImplementedError, "erf not implemented for precision higher than 53"
    81         return parent(1 - pari(float(x)).erfc())
     82        return parent(1) - parent(pari(x).erfc(prec))
    8283       
    8384    def _derivative_(self, x, diff_param=None):
    8485        """
  • sage/symbolic/expression.pyx

    diff --git a/sage/symbolic/expression.pyx b/sage/symbolic/expression.pyx
    a b  
    988988        EXAMPLES::
    989989
    990990            sage: complex(I)
    991             1j
     991            1j
    992992            sage: complex(erf(3*I))
    993             Traceback (most recent call last):
    994             ...
    995             TypeError: unable to simplify to complex approximation
     993            (1.0000000000000002+1629.8673238578601j)
    996994        """
    997995        try:
    998996            return self._eval_self(complex)
     
    29272925            sage: taylor(a*log(z), z, 2, 3)
    29282926            1/24*(z - 2)^3*a - 1/8*(z - 2)^2*a + 1/2*(z - 2)*a + a*log(2)
    29292927
    2930         ::
     2928        ::
    29312929
    29322930            sage: taylor(sqrt (sin(x) + a*x + 1), x, 0, 3)
    29332931            1/48*(3*a^3 + 9*a^2 + 9*a - 1)*x^3 - 1/8*(a^2 + 2*a + 1)*x^2 + 1/2*(a + 1)*x + 1
    29342932
    2935         ::
     2933        ::
    29362934
    29372935            sage: taylor (sqrt (x + 1), x, 0, 5)
    29382936            7/256*x^5 - 5/128*x^4 + 1/16*x^3 - 1/8*x^2 + 1/2*x + 1
    29392937
    2940         ::
     2938        ::
    29412939
    29422940            sage: taylor (1/log (x + 1), x, 0, 3)
    29432941            -19/720*x^3 + 1/24*x^2 - 1/12*x + 1/x + 1/2
    29442942
    2945         ::
     2943        ::
    29462944
    29472945            sage: taylor (cos(x) - sec(x), x, 0, 5)
    29482946            -1/6*x^4 - x^2
    29492947
    2950         ::
     2948        ::
    29512949
    29522950            sage: taylor ((cos(x) - sec(x))^3, x, 0, 9)
    29532951            -1/2*x^8 - x^6
    29542952
    2955         ::
     2953        ::
    29562954
    29572955            sage: taylor (1/(cos(x) - sec(x))^3, x, 0, 5)
    29582956            -15377/7983360*x^4 - 6767/604800*x^2 + 11/120/x^2 + 1/2/x^4 - 1/x^6 - 347/15120
    29592957
    2960  
    2961         Ticket #7472 fixed (Taylor polynomial in more variables) ::
     2958        TESTS:
     2959
     2960        Check that ticket #7472 is fixed (Taylor polynomial in more variables)::
    29622961 
    29632962            sage: x,y=var('x y'); taylor(x*y^3,(x,1),(y,1),4)
    29642963            (y - 1)^3*(x - 1) + (y - 1)^3 + 3*(y - 1)^2*(x - 1) + 3*(y - 1)^2 + 3*(y - 1)*(x - 1) + x + 3*y - 3
     
    75487547        -  ``multiplicities`` - bool (default: False); if True,
    75497548           return corresponding multiplicities.  This keyword is
    75507549           incompatible with ``to_poly_solve=True`` and does not make
    7551            any sense when solving inequality.
     7550           any sense when solving inequality.
    75527551       
    75537552        -  ``solution_dict`` - bool (default: False); if True or non-zero,
    75547553           return a list of dictionaries containing solutions. Not used
     
    75627561           Maxima's ``to_poly_solver`` package to search for more possible
    75637562           solutions, but possibly encounter approximate solutions.
    75647563           This keyword is incompatible with ``multiplicities=True``
    7565            and is not used when solving inequality. Setting ``to_poly_solve``
     7564           and is not used when solving inequality. Setting ``to_poly_solve``
    75667565           to 'force' (string) omits Maxima's solve command (usefull when
    75677566           some solution of trigonometric equations are lost).
    75687567       
     
    76107609            sage: solve(Q*sqrt(Q^2 + 2) - 1, Q)
    76117610            [Q == 1/sqrt(Q^2 + 2)]
    76127611            sage: solve(Q*sqrt(Q^2 + 2) - 1, Q, to_poly_solve=True)
    7613                         [Q == 1/sqrt(-sqrt(2) + 1), Q == 1/sqrt(sqrt(2) + 1)]
     7612            [Q == 1/sqrt(-sqrt(2) + 1), Q == 1/sqrt(sqrt(2) + 1)]
    76147613
    76157614        In some cases there may be infinitely many solutions indexed
    76167615        by a dummy variable.  If it begins with ``z``, it is implicitly