is this integration problem due to analytical continuation?

[Peter Pein]

Dear Maplers,

when trying to get a usable plot of

absint := unapply(abs(int(sqrt(sin(c*x)), x = 0 .. 2*Pi)), c);

I got convinced there's a bug in Maple. I'm using

> interface(version);
`Standard Worksheet Interface, Maple 15.01, Windows 7, June 2 2011 Build
ID 636299`


First, it is clear that for nonnegative integer parameters c the value
of the integral should be the same because absint(c) is or should be the
same as abs(Int(sqrt(sin(x))/c, x=0..2*Pi*c));

Numerical verification:

fai := c -> abs(int(sqrt(sin(c*x)), x = 0 .. 2*Pi, numeric = true));
`$`(fai(i), i = 1 .. 3)

gives three times the value evalf(4*GAMMA(3/4)^2/sqrt(Pi)) (modulo
floating point noise)

3.388852339, 3.388852339, 3.388852339

OK

But I am interested how the plot of absint looks like for real values of c.

If I try a plot of the function absint in the range c=1..2, I get a
weird result ( http://dl.dropbox.com/u/3030567/Maple/absint12.png );
increasing Digits to e.g. 30 does not help.

Using Mathematica, I got a result which looks promising to me (
http://dl.dropbox.com/u/3030567/Maple/mmaplot.png ). It has got values
around 3.88 at c=1,2 and is smooth in between.

Converting the Mathematica expression to Maple led to a surprise:

mmaint := convert("Abs[(-2*EllipticE[(Pi - 4*c*Pi)/4, 2] +
Sqrt[2/Pi]*Gamma[3/4]^2)/c]", FromMma);

(converted to mmaint := abs((-2*EllipticE(sin((1/4)*Pi-c*Pi),
sqrt(2))+sqrt(2/Pi)*GAMMA(3/4)^2)/c), because the parameters of the
elliptic functions are handled differently)

The values at integer parameters are:

evalf(`$`(mmaint, c = 1 .. 3));

2.396280469, 0.4e-9, .7987601563

not a too obvious relation to the correct values above :-(

and the plot does not show values near (or almost near) 3.8 for c=1 or
c=2 ( http://dl.dropbox.com/u/3030567/Maple/mmaint12.png )

Using limit on absint gives a useless result:

limit(absint(i), i = 1);

limit(2^(1/2)*abs((-GAMMA(3/4)^2-I*GAMMA(3/4)^2+(2*I)*Pi^(1/2)*hypergeom([-1/4,
1/2], [3/4], -(1/4)/(cos(Pi*i)^2*(-1+cos(Pi*i)^2)


I lost patience after a few minutes while trying to do the same plot for
the function "fai" above. But the attempt to evaluate it for c=17 led to:

foo := CodeTools:-Usage(fai(17), output = all);
Error, (in property/LinearProp/+) too many levels of recursion

and the value at c=1.7

bar := CodeTools:-Usage(fai(1.7), output = all);

bar := Record(realtime = 9.300, cputime = 9.267, bytesused = 534493216,
bytesalloc = 0, output = 3.417123397)

has nothing to do with

evalf(eval(mmaint, c = 17/10));

.5215125561

but is the same as calculated by Mathematica (here "maint" is the
Mma-expression converted above):

In[10]:= mmaint/.c->1.7
Out[10]= 3.41712


Do I have to reinstall Maple or is this 'normal' behaviour?

Thanks,
Peter

b.t.w.: Maxima 5.26 fails to get a numeric or symbolic result for

quad_qags(sqrt(sin(x)), x, 0, 2*%pi);
and
integrate(sqrt(sin(x)), x, 0, 2*%pi);
respectively

[G. A. Edgar]

In article <jihll2$j7c$1...@online.de>, Peter Pein <pet...@dordos.net>
wrote:

>
> I got convinced there's a bug in Maple. I'm using
>

Taking the bug to its essence, I came to an integral
that does not involve square-roots of negatives:

int(sqrt(sin(x/4)),x=0..2*Pi);

-4*2^(1/2)*EllipticK((1/2)*2^(1/2))+8*2^(1/2)*EllipticE((1/2)*2^(1/2))

evalf(%);
4.79256094

But of course that integrand is <= 1, so that integral cannot possibly
be larger than 2*Pi.

--
G. A. Edgar http://www.math.ohio-state.edu/~edgar/

[G. A. Edgar]

In article <290220120741092665%ed...@math.ohio-state.edu.invalid>, G.

A. Edgar <ed...@math.ohio-state.edu.invalid> wrote:

> In article <jihll2$j7c$1...@online.de>, Peter Pein <pet...@dordos.net>
> wrote:
>
> >
> > I got convinced there's a bug in Maple. I'm using
> >
>
> Taking the bug to its essence, I came to an integral
> that does not involve square-roots of negatives:
>
> int(sqrt(sin(x/4)),x=0..2*Pi);
>
> -4*2^(1/2)*EllipticK((1/2)*2^(1/2))+8*2^(1/2)*EllipticE((1/2)*2^(1/2))
>
> evalf(%);
> 4.79256094
>
> But of course that integrand is <= 1, so that integral cannot possibly
> be larger than 2*Pi.

Silly me, it *is* smaller than 2*Pi.

OK, backing up... Maple evaluates int(sqrt(sin(c*x)), x = 0 .. 2*Pi)

as

-2^(1/2)*(-GAMMA(3/4)^2-I*GAMMA(3/4)^2+(2*I)*Pi^(1/2)
*hypergeom([-1/4, 1/2], [3/4],-(1/4)/(cos(Pi*c)^2*(-1+cos(Pi*c)^2)))
*(sin(Pi*c)*cos(Pi*c))^(1/2))/(Pi^(1/2)*c)

and that is larger than 2*Pi for c close to zero. But maybe
this is still because of square-roots of negatives in there?
If I do instead
int(sqrt(sin(c*x)), x = 0 .. 2*Pi) assuming c>0, c<1/2
I get

-2^(1/2)*EllipticF(2*(sin(Pi*c)*cos(Pi*c)/(2*sin(Pi*c)
*cos(Pi*c)+1))^(1/2), (1/2)*2^(1/2))/c+2^(1/2)
*EllipticPi(2*(sin(Pi*c)*cos(Pi*c)/(2*sin(Pi*c)
*cos(Pi*c)+1))^(1/2), 1/2, (1/2)*2^(1/2))/c

which does stay below 2*Pi but for c=1/2 gives me the erroneous answer
0.

[Axel Vogt]

On 28.02.2012 05:36, Peter Pein wrote:
> Dear Maplers,
>
> when trying to get a usable plot of
>
> absint := unapply(abs(int(sqrt(sin(c*x)), x = 0 .. 2*Pi)), c);
>
> I got convinced there's a bug in Maple. I'm using
>

...

> Using Mathematica, I got a result which looks promising to me ( http://dl.dropbox.com/u/3030567/Maple/mmaplot.png ). It
> has got values around 3.88 at c=1,2 and is smooth in between.
>
> Converting the Mathematica expression to Maple led to a surprise:
>
> mmaint := convert("Abs[(-2*EllipticE[(Pi - 4*c*Pi)/4, 2] + Sqrt[2/Pi]*Gamma[3/4]^2)/c]", FromMma);

...

Despite the symbolic errors shown so far you can plot it in
a brute (Digits:=15), using evalc for the integrand and then

G:= c -> 1/2*abs( evalf(
Int( x -> abs(sin(c*x))^(1/2)*(1+signum(sin(c*x))), 0 .. 2*Pi,
method = _d01ajc, epsilon=eps) ) +
I* evalf(
Int( x -> abs(sin(c*x))^(1/2)*(1-signum(sin(c*x))), 0 .. 2*Pi,
method = _d01ajc, epsilon=eps) )
);

For eps:=1e-10 I get

[1,2,3]; map(G, %);

[3.38885233917592, 3.38885233917590, 3.38885233917591]

Plotting is *very* slow, but gives an impression, at least by

plot(G, 0 .. 2, numpoints=10);

[Axel Vogt]

On 29.02.2012 16:35, G. A. Edgar wrote:
> In article<290220120741092665%ed...@math.ohio-state.edu.invalid>, G.
> A. Edgar<ed...@math.ohio-state.edu.invalid> wrote:
>
>> In article<jihll2$j7c$1...@online.de>, Peter Pein<pet...@dordos.net>
>> wrote:

...

How about writing that as 1/c*Int(sin(xi)^(1/2),xi = 0 .. 2*Pi*c), c*x=xi ?

Then - up to damping by 1/c - this seems to be periodic via xi and it reduces
to the fractional part of c. And thus gives the sign change in c= 1/2.

And would reduce c=1.7 to c_new = 0.7 = 1/2 + 0.2 (and the correct damping).

Now for c=1/2 a direct call for 'value' seems to give a correct result:

Int(sin(xi)^(1/2),xi = 0 .. 2*Pi* 1/2);
evalf(%), evalf(value(%));

2.39628046947118, 2.39628046947114

While

Int(sin(xi)^(1/2),xi = 0 .. 2*Pi* t);
value(%) assuming 0<t, t<1/2; combine(%);

returns a result in elliptic functions, which equals zero for c=1/2.

But perhaps that is, what G.A. Edgar already said.

[Axel Vogt]

On 28.02.2012 05:36, Peter Pein wrote:

> Dear Maplers,
>
> when trying to get a usable plot of
>
> absint := unapply(abs(int(sqrt(sin(c*x)), x = 0 .. 2*Pi)), c);

First: for c=0.9 I get 2.68... using the expression 'mmaint', while a direct
numerical integration for Int( sqrt(sin(c*x)), x=0 .. 2*Pi) gives me the value
2.66253385496798+2.29880150861876*I, so abs(ofThat) = 3.51..., Digits:=15.

The integral is the same as 1/c*Int(sqrt(sin(x)),x = 0 .. 2*c*Pi)

I used Peter's idea to ask MMA:

It gives "-2*EllipticE[(Pi - 2*x)/4, 2]" as an anti-derivative for sqrt(sin(x)),
with the surprise, that "EllipticE[z, a]" is EllipticE(sin(z),sqrt(a)) in Maple.

This is continous in x and at least for 0 < sin(x) I convinced myself by diff
and some manual work, that it _is_ an anti-derivative. From that one has the
integral for 0 < c < 1/2.

However for x=3/2*Pi it has a 'constant' jump, in x=Pi*5/3 one directly can see,
that it is not an anti-drivative. But after manual work I got an continous anti-
derivate

A := 2*PIECEWISE(
[-EllipticE(cos(1/4*Pi+1/2*x),2^(1/2)),
x <= 3/2*Pi],
[EllipticE(cos(1/4*Pi+1/2*x),2^(1/2))+(2+2*I)*EllipticE(1/2*2^(1/2),2^(1/2)),
3/2*Pi < x])

A proof is missing for that. But from there one finds

J := c -> piecewise(
c <= 3/4,
2*EllipticE(1/2*2^(1/2),2^(1/2))-2*EllipticE(cos(1/4*Pi+c*Pi),2^(1/2)),
3/4 < c,
2*EllipticE(cos(1/4*Pi+c*Pi),2^(1/2))+(6+4*I)*EllipticE(1/2*2^(1/2),2^(1/2)));

with J(c) = Int(sqrt(sin(x)),x = 0 .. 2*c*Pi) for 0 <= c <= 1.

If c is an integer then (by periodics of sin(x)) one just has integer multiples
of J1:= Int( sqrt(sin(1*x)), x=0 .. 2*Pi) and Maple find that as the value of
J1 = (4+4*I)*(-1/2*EllipticK(1/2*2^(1/2))+EllipticE(1/2*2^(1/2)))*2^(1/2).

From that Int(sin(x)^(1/2),x = 0 .. 2*c*Pi) = floor(c)*J1 +
Int(sin(x)^(1/2),x = 0 .. 2*frac(c)*Pi), just using 0 <= c=integer+fraction

and finally

Int( sqrt(sin(c*x)), x=0 .. 2*Pi) = '(floor(c)*J1 + J(frac(c)))/c';

For numerical tests against direct integration one should use evalf[Digits-2],
otherwise already for moderate c the numerical intgration returns unevaluated,
indicating errors to high.

Plotting now is reasonable fast:

'abs((floor(c)*J1 + J(frac(c)))/c)';
plot(%, c=0 .. 7);

[Peter Pein]

Wow,

thank you and G. A. Edgar very much :)

Peter