Two different answers.

The following appears to give two different answers.

Any suggestions?

 

> restart:
> Re(z);
> Re(a);
> Re(b);
> assume(a>0,b>0,z>0);
> H:=sqrt(b*(1+z)^3+a);
> bot:=int(1/H,z=0..infinity);
> evalf(subs(a=3/4,b=1/4,bot));
> H1:=sqrt((1/4)*(1+z)^3+3/4);
> int(1/H1,z=0..infinity);
> evalf(%);

Robert Israel's picture

Branch cuts

This seems to be another branch cut problem. 

I think it's a bit more convenient to take a = 1. 

> J := eval(Int(1/H, z=0..infinity), a=1);

= Int(1/((b+3*b*z+3*b*z^2+b*z^3+1)^(1/2)),z = 0 .. infinity)

> W := simplify(value(J));

= 1/3*1/b^(1/3)*3^(3/4)*EllipticF(2*3^(1/4)*(1+b^(1/3))^(1/2)/(1+b^(1/3)+3^(1/2)),1/4*2^(1/2)*3^(1/2)+1/4*2^(1/2))

 

Now EllipticF(z,k) has a branch point at z=1, corresponding to b = -10 + 6*sqrt(3) = .39230485..., and plotting J and W near there shows that J (the actual integral, as evaluated numerically) appears to be smooth there but W is not.

> plot([W,J], b= 0.35 .. 0.45, colour=[blue,red]);

My guess is that for b < -10 + 6*sqrt(3), the wrong branch of EllipticF is being used.

Comment viewing options

Select your preferred way to display the comments and click "Save settings" to activate your changes.
}