can not do integral of sqrt(sin(x))

[pallab]

It seems sympy can not integrate sqrt(sin(x)).

I did the following:

import sympy as sm
from sympy.abc import x,y,z

tointegrate=sm.sqrt(sm.sin(y))
sm.integrate(tointegrate)

output is : Integral(sqrt(sin(y)), y)

After a simple change of variable the integral is doable:

def ytoz(z):
    return(sm.asin(z))

def dydz(z):
    return(sm.diff(ytoz(z),z))
          
subdict={y:cvytoz(z)}
sm.integrate(dydz(z)*tointegrate.subs(subdict))

output is : z**(3/2)*gamma(3/4)*hyper((1/2, 3/4), (7/4,), z**2*exp_polar(2*I*pi))/(2*gamma(7/4))

best,

Pallab

[Aaron Meurer]

WolframAlpha gives an answer in terms of an elliptic integral:

http://www.wolframalpha.com/input/?i=integrate(sqrt(sin(x)),%20x)

So I imagine that if we added those that we would be able to compute this.

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To view this discussion on the web visit
> https://groups.google.com/d/msg/sympy/-/-aC1_0EHC8YJ.
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to
> sympy+un...@googlegroups.com.
> For more options, visit this group at
> http://groups.google.com/group/sympy?hl=en.

[rl]

The problem I see is *how* to add these to the integration routines.
There is no nice Meijer-G representation. And I suppose that Risch
can not handle these. At least not w/o major extensions.

[Aaron Meurer]

But if Pallab's answer is correct, then it can be, right, at least for
this specific one?

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To view this discussion on the web visit

> https://groups.google.com/d/msg/sympy/-/9ztiJTc6CxEJ.

[Aaron Meurer]

According to http://en.wikipedia.org/wiki/Elliptic_integral, every
elliptic integral can be reduced to an expression containing only the
three Legendre canonical forms. So I wonder if that is algorithmic.

Aaron Meurer

[pallab]

The answer seems to be wrong.
I do the following to numerically integral evaluation:

import numpy as np
import scipy as sp
import scipy.integrate

f=lambda x:np.sqrt(np.sin(x))

sp.integrate.quad(f,0.6,0.7)

get the following value:

0.07776347731181982

which matches with mathematica.

whereas

fI=sm.lambdify(z,sm.integrate(dydz(z)*tointegrate.subs(subdict)).subs(sm.exp_polar(2*sm.I*sm.pi),1.0))

and

fI(0.6)-fI(0.7)

gives:

mpf('-0.10638922100995812')

[pallab]

By answer wrong I mean the answer given by sympy,

[Joachim Durchholz]

Am 07.07.2012 01:18, schrieb Aaron Meurer:
> According to http://en.wikipedia.org/wiki/Elliptic_integral, every
> elliptic integral can be reduced to an expression containing only the
> three Legendre canonical forms. So I wonder if that is algorithmic.

The article says you need to find the "appropriate reduction formula".
Following that link game list of known structural patterns and a
reduction formula for each.

So it seems heuristic, with a set of well-known heuristics.

[Tom Bachmann]

Hi,

There are a couple of things to say here. First I think you have
confused when you give boundaries via x or via z. After substitution,
the integrand becomes sqrt(z)/sqrt(1-z**2), so let's discuss this. Then
we have:

In [26]: i = Integral(sqrt(z)/sqrt(1-z**2), (z, 0.3, 0.4))
In [29]: i.doit().subs(exp_polar, exp).n()
Out[29]: 0.0631741978430877
In [30]: i.n()
Out[30]: 0.0631741978430876

So this seems to work.

Some notes. First we need to replace exp_polar by exp here. This is
permitted, since |z| < 1. Arguably this is a shortcoming of evalf of
hyper. In fact the meijerg() evalf can evaluate for arbitrary polar
numbers, hijacking some obscure mpmath functionality. It is possible to
implement similar functionality for hyper(), but that's more work. At
the very least it should work for |z| < 1...

Since the integral can be done in terms of elliptic integrals this
hypergeometric function should have a closed form, but I cannot find it
documented anywhere.

Tom

> --
> You received this message because you are subscribed to the Google
> Groups "sympy" group.
> To view this discussion on the web visit

> https://groups.google.com/d/msg/sympy/-/867mUPbpSJAJ.

[rl]

A first start:

http://github.com/sympy/sympy/pull/1408

[pallab]

It does work. Sorry that I made a mistake.