Fresnel integrals

[rws]

The SymPy-Sage interface raises an error because there are no Fresnel integral functions in Sage:

sage: integrate(sin(x^2), x, algorithm='sympy')
...
AttributeError: 

see:

sage: import sympy
sage: x = sympy.Symbol('x')
sage: sympy.integrate(sympy.sin(x^2), x)
3*sqrt(2)*sqrt(pi)*fresnels(sqrt(2)*x/sqrt(pi))*gamma(3/4)/(8*gamma(7/4))

They are holonomic so they should be in Sage.

• ​https://en.wikipedia.org/wiki/Fresnel_integral
• ​http://mathworld.wolfram.com/FresnelIntegrals.html

[mforets]

This is a tentative implementation of sine Fresnel integral. The numerical evaluation relies on mpmath. I doubted about how to write the custom evaluation _eval_ method. In fact I have two failing doctests. I expect:

sage: fresnel_sin(0)
0
sage: fresnel_sin(x).subs(x==0)
0

whereas i get fresnel_sin(0) in both cases. Any hints?

[mforets]

Oh, and the beta(... --> beta... is because the missing parentheses were confusing my text editor, so that the entire file from that line on was shown as comments. I should remember to undo this change :)

[rws]

Returning None (or just passing) in _eval_ amounts to returning a held expression. _evalf_ is only used for inexact arguments or .n():

sage: fresnel_sin(0.)
0.000000000000000
sage: fresnel_sin(0).n()
0.000000000000000

For results to exact or symbolic arguments you need to implement them in _eval_.

[rws]

You would just need to handle 0 and oo in _eval_, and handle oddness before it, i.e. if arg.is_negative(): return -fresnel_sin(-arg).

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​f158660

symbolics in fresnel_sin eval

[mforets]

ok, that fixed the 0, oo and -oo evaluations. i added two new tests for the I*oo cases, but this is still wrong. since

sage: SR(I*oo).is_infinity()
True

probably i could use x.is_real() too.

[rws]

You mean unsigned_infinity?

[mforets]

i should relate unsigned_infinity to test if the argument is complex and infinite? the doc of unsigned_infinity is precarious:

sage: unsigned_infinity?
Type:           UnsignedInfinity
String form:    Infinity
File:           ~/Tools/sage-master/local/lib/python2.7/site-packages/sage/rings/infinity.py
Docstring:         Initialize "self".
Init docstring:    Initialize "self". 

i'm tempted to add a ticket to explain this object more.

[rws]

Replying to mforets:

i should relate unsigned_infinity to test if the argument is complex and infinite?

No, I was not sure what you had in mind. If you want to test if the argument is complex and infinite then do: if arg.is_infinity() and not arg.is_real(). Technically a False for arg.is_real() just means "I don't know" so you get True for SR(unsigned_infinity).is_real(), but the correct behaviour must wait until we have is_complex(). Just use not arg.is_real() for now.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​ada47e4

handle complex infinities

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​869020b

add cosine Fresnel integral

[mforets]

oh, there's no conversion back to sage:

sage: fresnel_sin(x)._sympy_() # ok
fresnels(x)
sage: fresnel_sin(x)._sympy_()._sage_() # not working!

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​e930137

add sympy -> sage conversion, and doctest

[rws]

See how easy it is to add conversions now.

[mforets]

See how easy it is to add conversions now.

awesome! it took me what, 15 minutes? and it was the first time i was doing that.

[rws]

Patchbot has relevant doctest errors.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​2a0e907

fix laplace test

[mforets]

Indeed, thanks for noticing it!

Some patchbots reported an error RuntimeError: cases argument not a sequence happening with some random tests. This seems to me unrelated to this ticket, but maybe i'm wrong.

[rws]

Indeed, the cases fail has nothing to do with this ticket. It's a bug that is rare and happens randomly (#24284).

[rws]

Docs are looking good. The patchbot errors are unrelated so I would be positive. However, the derivative is missing and there is a typo in the docs: write \operatorname{C} not \operatorname(C) in the formula. Could you please add that?

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​e6b884e

fix operatorname and add derivatives

[rws]

Thanks for your work.

[mforets]

well, thanks for the review and all your help.

[vbraun]

Merge conflict

[rws]

I merged with 8.2.beta0 + #24171 to avoid another conflict.

---

New commits:

​669ea22	24171: Formal set membership function
​34c9773	24171: make RealSet inherit from Set_generic
​af82717	24171: remove RealSet import; SR set coercion doctests
​592ce87	24171: Py3 fixes
​aa11a3e	Merge branch 'develop' into t/24171/formal_set_membership_function
​44c2bdb	Merge branch 'u/rws/formal_set_membership_function' of git://trac.sagemath.org/sage into t/24212/fresnel_integrals

[mforets]

well done, thanks.

[vbraun]

See patchbot

[rws]

The depended-on #24171 needs a patch so I'll squash it all, keeping authorship info correct.

[rws]

Actually the fails are due to the blockers #24284 and #24378 but they are not reviewed yet.

[rws]

I moved the code into error.py as the Fresnel integrals are related to error functions. This also makes it independent from other tickets.

---

New commits:

​e3cd31a

24212: Fresnel integrals

[vbraun]

File "src/sage/symbolic/random_tests.py", line 265, in sage.symbolic.random_tests.?
Failed example:
    random_expr(50, nvars=3, coeff_generator=CDF.random_element) # random
Exception raised:
    Traceback (most recent call last):
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 517, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/doctest/forker.py", line 920, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.symbolic.random_tests.?[2]>", line 1, in <module>
        random_expr(Integer(50), nvars=Integer(3), coeff_generator=CDF.random_element) # random
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 288, in random_expr
        return random_expr_helper(size, internal, leaves, verbose)
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 232, in random_expr_helper
        children = [random_expr_helper(n+1, internal, leaves, verbose) for n in nodes_per_child]
      File "/home/buildbot/slave/sage_git/build/local/lib/python2.7/site-packages/sage/symbolic/random_tests.py", line 235, in random_expr_helper
        return r[1](*children)
      File "sage/symbolic/expression.pyx", line 3537, in sage.symbolic.expression.Expression.__invert__ (build/cythonized/sage/symbolic/expression.cpp:25798)
        return 1/self
      File "sage/structure/element.pyx", line 1660, in sage.structure.element.Element.__div__ (build/cythonized/sage/structure/element.c:12835)
        return coercion_model.bin_op(left, right, div)
      File "sage/structure/coerce.pyx", line 1081, in sage.structure.coerce.CoercionModel_cache_maps.bin_op (build/cythonized/sage/structure/coerce.c:9515)
        return PyObject_CallObject(op, xy)
      File "sage/structure/element.pyx", line 1655, in sage.structure.element.Element.__div__ (build/cythonized/sage/structure/element.c:12758)
        return (<Element>left)._div_(right)
      File "sage/symbolic/expression.pyx", line 3518, in sage.symbolic.expression.Expression._div_ (build/cythonized/sage/symbolic/expression.cpp:25588)
        raise ZeroDivisionError("Symbolic division by zero")
    ZeroDivisionError: Symbolic division by zero
**********************************************************************
1 item had failures:
   1 of   5 in sage.symbolic.random_tests.?
    [35 tests, 1 failure, 22.02 s]
----------------------------------------------------------------------
sage -t --long src/sage/symbolic/random_tests.py  # 1 doctest failed

[rws]

random_expr is really interesting. It should be made easier to debug, and errors caught/not caught should be reviewed because TypeErrors from ECL represent real bugs and RuntimeErrors from Pynac should be rebranded in the Python layer. Also sometimes it takes really really long time. After that, someone should do a day-long run to catch any remaining problems.

[rws]

The reason why the doctest fails is that, despite fixing the random seed, the functions and symbols picked by random_expr are different. This, at some point in the buildup of the expression, leads to the sub-expression kronecker_delta(-0.6165326641917295 + 0.15117833554974136*I, e) which immediately evaluates to 0. The 0 progresses eventually into a denominator, causing the error. random_expr is not written to avoid such errors if they happen, nor does it avoid other exceptions thrown. The author of the doctest apparently just chose a random seed that happens to not throw errors.

Of course the actual bug is that random_expr changes behaviour when unrelated code is changed.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​a17755c	24425: Fix inherently failing random_expr doctest
​96521a4	Merge branch 't/24425/fix_inherently_failing_random_expr_doctest' into t/24212/24212-2

[rws]

Dependency needs review.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​e4f1592

Merge branch 'develop' into t/24212/24212-2

[mforets]

The dependency has been approved and merged. Thanks.

I've merged with the latest development version and now i'm checking if the branch applies.

[mforets]

New commits:

​e3cd31a	24212: Fresnel integrals
​96521a4	Merge branch 't/24425/fix_inherently_failing_random_expr_doctest' into t/24212/24212-2
​e4f1592	Merge branch 'develop' into t/24212/24212-2
​61b2f76	Merge branch 'develop' into t/24212/24212-2

[mforets]

$ ./sage -t src/sage/calculus/calculus.py src/sage/functions/*
...
Git branch: t/24212/24212-2
Using --optional=mpir,python2,sage
...
----------------------------------------------------------------------
All tests passed!
----------------------------------------------------------------------
Total time for all tests: 106.7 seconds
    cpu time: 111.3 seconds
    cumulative wall time: 104.6 seconds

[tscrim]

Quick question: do you want to change the test in calculus.py? In particular, do you still want to test that Sage handles the case when something in sympy gives something Sage doesn't know how to handle?

[mforets]

Alright. I found this one:

sage: laplace(cos(-1/t), t, s, algorithm='sympy')
...
AttributeError: Unable to convert SymPy result (=meijerg(((), ()), ((-1/2, 0, 1/2), (0,)), s**2/16)/4) into Sage

(side note: i didn't read the syntax of sympy's meijerg function, but here is a quick comparison to ​WA's answer)

[rws]

Replying to mforets:

sage: laplace(cos(-1/t), t, s, algorithm='sympy')
...
AttributeError: Unable to convert SymPy result (=meijerg(((), ()), ((-1/2, 0, 1/2), (0,)), s**2/16)/4) into Sage

For that, there is #17970.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​f49ce18

add doctest with attributeerror

[rws]

When we have #17970 we'll likely convert the G expression to the Fresnel expression in Sage.