This is the mail archive of the gsl-discuss@sources.redhat.com mailing list for the GSL project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]

Re: Complete Elliptic Integrals


trophime christophe wrote:
> 
> I have a few basic questions about complete elliptic integrals :
> 
> 1. E(1) should return 1, right?
> result = gsl_sf_ellint_Ecomp_impl(m(r1,z11), GSL_PREC_DOUBLE, &E)
> returns :
> E.val = 0 and  result == GSL_EDOM

Yes, I missed an endpoint condition. I fixed it.
In principle, the same sort of problem will occur
with the incomplete E, when phi is near pi/2, so
I fixed that as well (although you could argue
the fix is a hack). I don't think a similar
cancellation problem will occur in the way
that I used the Carlson representations for
the other functions (fingers crossed).


> 2.Looking in the source I have noticed on line 246
> (gsl_sf_ellint_RJ_impl) :
> if(x < 0.0 || y < 0.0 || y < 0.0) {
> Should I read z instead of the second y??

Yup. I fixed that too.


> 3.What are gsl_sf_ellint_Dcomp_impl? I can't find this notation anywhere

Unfortunately I don't have Carlson's paper handy,
so I can't quote any definition other than the
one in terms of R_D. I put the reference to
Carlson's paper in the header file and quoted
that definition; that's the best I can do right now.
My guess is that this is the only definition.


> Is the result as accurate as for any phi != Pi/2?

As long as there is no endpoint cancellation problem
(and I don't think there is for this case), then it
should be just as accurate. But if you find otherwise,
please let me know.


> why not providing a complete version of the third kind integrals?

Hmmm. I don't remember if there was a good reason
for that. If I look at Carlson's paper, I might
remember. Probably I just left it out.


I just committed the changes to CVS.

Thanks for the report.

-- 
G. Jungman

Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]