Home › Forums › Get Help › How do I ...? with Maple

Hypergeometric function with small arguments

Posted on 2008-02-20 00:42 By vishwa_r (9)
Categories:

I am trying to evaluate the hypergeometric function : hypergeom([1],[1.5],-z) where z is in the order of 1e-6. Its taking very long time where as if I use maple through matlab, then its quite fast !!. Evaluation of Gauss-hypergeometric functions for the similar arguments in maple is quite fast.  Is there any trick I can use for the fast computation of the confluent hypergeom in the case above?.

Thanks for any help

‹ Generate bitstrings of length n that do not contain two consecutives ones Interpret 2 symbolic solve solutions when only 1 numerical solution exists? ›
Robert Israel's picture

Long time?

Comment on 2008-02-20 01:33 from Robert Israel (Maple Rating 4 1555)

I don't find it slow at all.  For example:

Q:= z -> hypergeom([1],[1.5],-z);

L:= [seq(RandomTools[Generate](float(range=0.5e-6 .. 5e-6)), 
    count=1..10000)]:

ti:= time(): for t in L do Q(t) end do:  time()-ti;

         0.875

So on my computer it took 0.875 seconds of CPU time to do 10000 of these hypergeom evaluations.  I wouldn't consider that slow at all, as Maple functions go.  Compare that to the equivalent "closed form" version:

Q2:= z -> evalf(-1/2*Pi^(1/2)*(-1+erfc((-z)^(1/2)))/exp(z)
     /(-z)^(1/2));

ti:= time(): for t in L do Q2(t) end do:  time()-ti;

      12.453

I wonder if there's really something else going on.  Perhaps if you posted your actual code we could have something to suggest.

Hypergeometric function with small arguments : Long time?

Comment on 2008-02-20 02:40 from vishwa_r (9)

Hi Robert,

Sorry, I did a mistake while posting. Small number I meant a negetive number of order 1e6. My apologies.

My code is a loop with a part of it is just evaluation of confluent hypergeom of type : hypergeom([1],[1.5],-1e6) . Maple just seems like hanged with this input whereas matlab seems to be quite quick. 

In matlab   tic;hypergeom([1],[1.5],-1e6);toc
Elapsed time is 0.075050 seconds.

what am I missing here?

Thanks a lot for the help

Robert Israel's picture

Large negative vs small

Comment on 2008-02-20 03:30 from Robert Israel (Maple Rating 4 1555)

Ah... we generally call that "a large negative number", not "a small number".

It looks like Maple is basically adding up the hypergeometric series term by term. 
Of course with z so large (whether positive or negative) it's going to take an awful lot of terms to do that.  I have no idea what matlab does with this.  But here in Maple it's much better to use the closed-form expression

-1/2*Pi^(1/2)*(-1+erfc((-z)^(1/2)))/exp(z)/(-z)^(1/2)

or, even better, the asymptotic series of this as z -> infinity:

1/(2*z)+1/(4*z^2)+3/8/z^3+15/16/z^4+O((1/z)^(9/2))
JacquesC's picture

Surprised

Comment on 2008-02-20 11:18 from JacquesC (Maple Rating 4 1699)

I am really surprised by how badly Maple is performing here.  Fast and accurate numerical evaluation of even very complex mathematical functions is one of the unsung strengths of Maple.  [Mathematica's basic arithmetic is way faster, but Maple's implementations for special functions are much more accurate and yet still fast].  Especially since I know that a huge amount of effort has been spent on the hypergeometric family.  One of the standard tricks is to use different algorithms over different regions (for both accuracy and speed), so that for this function I fully expected the asymptotic expansion to be used.  I guess that 1F1 is still on the 'to do' pile!

Slow computation of 1F1(a,b,c)

Comment on 2008-11-07 16:38 from fernando (2)

To date the computation of regular confluent hypergeometric function is slow for large |az| values. I found that Matlab does not perform well when a or z is large. The only program I know of, on Mathworks was a normal series expansion of the 1F1, which I think is not actually a good method, as it requires alot of terms to be added up. So I  cant actually see that Matlab gives more accurate or even faster results than Maple or even Mathematica, anyway

What I do want to know is where does one get the closed form solution Robert gave as

-1/2*Pi^(1/2)*(-1+erfc((-z)^(1/2)))/exp(z)/(-z)^(1/2)

I am still a masters student and still learning and I can't seem to see or find an article where that closed form solution are derived from. Probably missing something since I dont see any the values of a or b in it, maybe already incorporated?

Can somebody help.

And dont worry next year a new method will be suggested for solving the 1F1(a,b,c) for all parameters a, b, and c!  :)

Cheers.

 

already specialized

Comment on 2008-11-07 17:25 from jakubi (Maple Rating 4 943)

It is a=1 and b=3/2:

convert(hypergeom([1],[3/2],-z),erfc);
                          1/2                1/2
                        Pi    (-1 + erfc((-z)   ))
                   -1/2 --------------------------
                                  1/2
                              (-z)    exp(z)

<p>Thanks...,</p> <p>Okay so

Comment on 2008-11-07 18:23 from fernando (2)

<p>Thanks...,</p>
<p>Okay so it only counts for certain instances. When a and b are complex as in for example the problem for pricing Asian options where a and b are complex&nbsp; with z large negative, we again have the problem of slow computational time since |az| increases. </p>
<p>Can&nbsp; you advise me to a paper, if any,&nbsp; where one can compute 1F1 faster for the given a,b and z.</p>
<p>Thank you.</p>

Axel Vogt's picture

rationalize inputs?

Comment on 2008-11-07 18:48 from Axel Vogt (Maple Rating 4 994) 
If you convert to rationals before, you may go the following way:

  KummerM(mu, nu, z); convert(%,hypergeom); convert(%,KummerM);

                          KummerM(mu, nu, z)
                       hypergeom([mu], [nu], z)
                          KummerM(mu, nu, z)

This says: you can use 'KummerM' instead.

And rationals will give you what you want:

For the (once) posted example I have no problems with an immediate
answer (while for floats it is a pain):

  KummerM(a,b,z); subs(a=1,b=1.5,z= -1/10^6,%); evalf(%);
                             0.9999993333

  hypergeom([1],[1.5],-1/10^6);
                             0.9999993333

for large abs(z)

Comment on 2008-11-07 19:10 from jakubi (Maple Rating 4 943)

ie

KummerM(a,b,z); subs(a=1,b=3/2,z= -10^6,%); evalf(%);

 
it is also slow.

Mathematica computes faster

Comment on 2008-11-07 20:01 from jakubi (Maple Rating 4 943)

In Mathematica 5:

N[Hypergeometric1F1[1+I,2*I,-10^6]]

6.739100053151782*10^(-9) + 2.0757700378175314*10^(-7)*I

is very fast.

Comment viewing options

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