The relation between the exponential integrals E1 and Ei is
(x real and positive):
Ei(x) = -E1(-x) - pi*1i
and so the following results are correct:
E1( 3.0) = 0.01304838...
E1(-3.0) = -9.933833... - 3.141593i
Ei( 3.0) = 9.933833...
Ei(-3.0) = -0.013048... - 3.141593i
In your function I cannot calculate eint(big(-3.0)), it returns
"nan with 256 bits of precision" (never seen `nan` with so many digits).
So I assume your function computes Ei, and not the E1 exponential integral
while, e.g., The MATLAB function expint computes E1.
The documentation of MPFR for function mpfr_eint is not quite clear.
I also got a bit worried by "... a future version might use another
definition".
On Tuesday, June 3, 2014 8:31:38 AM UTC+2, anonymousnoobie wrote:
>
>
>
> On Monday, June 2, 2014 12:19:00 PM UTC-6, anonymousnoobie wrote:
>>
>> Does julia have a package or library to compute exponential integral
>>
>> http://en.wikipedia.org/wiki/Exponential_integral
>> <http://en.wikipedia.org/wiki/Exponential_integral>
>>
>> for real positive argument arbitrary precision ?? I am new to julia so
>>
>> explain appropriately. thanks.
>>
>
>
> function eint(x::BigFloat)
>
> z = BigFloat()
>
> ccall((:mpfr_eint, :libmpfr), Int32, (Ptr{BigFloat}, Ptr{BigFloat},
> Int32), &z, &x, Base.MPFR.ROUNDING_MODE[end])
>
> return z
>
> end
>
> eint(big(3.)),exp(-3.)*log(1.+1./3)#upper bound of E1(3.)
>
> #lower bound is (exp(-3.)/2)*log(1 + 2./3)
>
> Out[17]:
> (9.933832570625416558008336019216765262990653022987582119693034065242496228332754e+00
>
> with 256 bits of precision,0.014322847009365602)
>
>
> =================================================================================
>
> I implemented the function THANKs
> But I don't think the numerical answer is correct. E1(3) is about 9.9
> from the program but E1(x) is bracket-
> ed between (exp(-x)/2)*log(1+2/x) and exp(-x)*log(1+1/x).
> Also the online website
> http://www.miniwebtool.com/exponential-integral-calculator/?x=-3
>
> calculates Ei(-3). as *-0.013048381094197. This is really confusing but
> E1(x)=-Ei(-x)*
> *So I believe the numerical answer to the program eint(3) is not E1(3)*
>
>