Dieter Kaiser wrote:
> Hello Ray,
>
> again thank you very much for your report.
>
> The Exponential Integrals were a first work of me. Thus, I have seen there could
> be done a lot much better. (Today I have implemented the noun, verb, alias and
> reversealias property and improved the checks for the special values 0 and 1).
>
> The value of expintegral_e(3,0.5) -> .2216043642751785 is definitly correct. I
> have checked it again with the value from functions.wolfram.com.
>
> When I transform the function to a representation in terms of expintegral_ei and
> do the numerically calculation I can verify the result:
>
> (%i36) ratsimp(expintegral_e(3,1/2));
> (%o36) (sqrt(%e)*(2*log(2)+log(-1/2)-2*expintegral_ei(-1/2)-log(-2))+4)
> /(16*sqrt(%e))
> (%i37) %,numer;
> (%o37) 0.22160436427518
>
> With CMUL and ECL we get an additional complex term. I think I have to look to
> the code to find all the places where we have explicit to declare the variables
> to be a float or complex value.
If you haven't figure it out already, the complex term comes from
expintegral_ei(-1/2). The other parts in %o36 actually come out to
1/4/sqrt(%e), both symbolically and numerically.
The issue is the expression
(- (* 0.5 (- (log z) (log (/ 1 z)))) (log (- z)))
in expintegral-ei.
If z is -0.5, this is 0. If z is #c(-0.5 0), it is #c(0 3.14159).
I'll have to think some more about which is correct....
Ray