LabWindows/CVI

cancel
Showing results for 
Search instead for 
Did you mean: 

Gaussian hypergeometric function 2F1(a,b;c;z)

Solved!
Go to solution

Hi,

 

this is a somewhat specific question Smiley Happy

 

I was very happy when I realized that the Advanced Analysis Library of CVI provides the Gaussian or ordinary hypergeometric function 2F1(a,b;c;z), GaussHG.

 

Unfortunately I have not been able to get the expected results, so I have to assume that I am doing something wrong and I would be grateful for assistance.

 

Here is the simple code trying to calculate the function

 

    for ( index = 0;
          index <= 1000;
          index ++ )
    {
        x [ index ] = -10.0 + index * 0.0108;
        y_1 [ index ] = GaussHG ( x [ index ], 0.5, 2.0, 1.5 );
    }

 

The first 833 elements of the y_1 array yield NaNs, the following numbers are oscillating in the range from +0.57 ... +3.29

 

Here is what I would have expected, using the online function plotter from Wolfram Research:

 

23478i2A2D7E05F3A10E16

 

Help is very much appreciated!

Thanks,

Wolfgang

0 Kudos
Message 1 of 13
(5,750 Views)

Hello NI,

 

I have added some extra data to assist you in calrifying this issue, the data calculated by CVI are depeicted here:

23570iDDB41CEC9DB91005

 

As you can see, the function works reasonably only for not too negative x values. I would appreciate it a lot if this behaviour could be confirmed and fixed as soon as possible - I really need this function and do not want to implement it myself.

 

Thanks!

0 Kudos
Message 2 of 13
(5,723 Views)

Hi,

 

for easier reproducing the issue I am attaching a simple demo project - would be nice to receive some feedback Smiley Happy

0 Kudos
Message 3 of 13
(5,694 Views)

Hi Wolfgang,

 

I think this is a bug.  

 

Thanks a lot for reporting this to us!

 

Joyce

 

0 Kudos
Message 4 of 13
(5,670 Views)

Hi JoyceXe,

 

Thank you for addressing this topic. Sorry for being impatient this time but I need this function urgently. Is there a chance that you might provide a fix prior to CVI2010?? Otherwise I would need to implement it myself and obviously I'd like to avoid that... So please keep me updated about the progress with this issue.  Thanks again,

 

Wolfgang

0 Kudos
Message 5 of 13
(5,665 Views)
Solution
Accepted by topic author Wolfgang

Hi Wolfgang,

 

We have reported this issue as a bug. For a quick workaround, I suggest you take advantage of the following transform:

F(a,b,c,x) =  F(a, c - b, c, x/(x-1)) / (1-x)^a, which could map x from (-Inf, 0] to [0,1). 

 

You could revise the original code a bit if performance is not a concern:

 

for ( index = 0;index <= 1000;index ++ )
    {
        x [ index ] = -10.0 + index * 0.0108;
        if(x[ index ]>=0)
              y_1 [ index ] = GaussHG ( x [ index ], 0.5, 2.0, 1.5 );
        else //use the transform to map negative x value to [0,1)
              y_1[ index ] = GaussHG ( x [ index ] / (x [ index ] - 1.0), 0.5, -0.5, 1.5 ) / pow(1.0-x [ index ], 0.5);
     }

for ( index = 0;index <= 1000;index ++ )

    {

        x [ index ] = -10.0 + index * 0.0108;

        if(x[ index ]>=0)

              y_1 [ index ] = GaussHG ( x [ index ], 0.5, 2.0, 1.5 );

        else //use the transform to map negative x value to [0,1)

              y_1[ index ] = GaussHG ( x [ index ] / (x [ index ] - 1.0), 0.5, -0.5, 1.5 ) / pow(1.0-x [ index ], 0.5);

    }

 

Hope this workaround could help. Please feel free to let us know if you have any more questions, thanks!

 

Arthur

Message 6 of 13
(5,635 Views)

Hi Wolfgang,

 

We have reported this issue as a bug. For a quick workaround, I suggest you take advantage of the following transform:

F(a,b,c,x) =  F(a, c - b, c, x/(x-1)) / (1-x)^a, which could map x from (-Inf, 0] to [0,1). 

 

You could revise the original code a bit if performance is not a concern:

 

for ( index = 0;index <= 1000;index ++ )
    {
        x [ index ] = -10.0 + index * 0.0108;
        if(x[ index ]>=0)
              y_1 [ index ] = GaussHG ( x [ index ], 0.5, 2.0, 1.5 );
        else //use the transform to map negative x value to [0,1)
              y_1[ index ] = GaussHG ( x [ index ] / (x [ index ] - 1.0), 0.5, -0.5, 1.5 ) / pow(1.0-x [ index ], 0.5);
     }

for ( index = 0;index <= 1000;index ++ )

    {

        x [ index ] = -10.0 + index * 0.0108;

        if(x[ index ]>=0)

              y_1 [ index ] = GaussHG ( x [ index ], 0.5, 2.0, 1.5 );

        else //use the transform to map negative x value to [0,1)

              y_1[ index ] = GaussHG ( x [ index ] / (x [ index ] - 1.0), 0.5, -0.5, 1.5 ) / pow(1.0-x [ index ], 0.5);

    }

 

Hope this workaround could help. Please feel free to let us know if you have any more questions, thanks!

 

Arthur

0 Kudos
Message 7 of 13
(5,633 Views)

Hi Arthur,

 

Well, this workaround was a very nice surprise, thanks! Unfortunately, it is not yet fully working as expected, I got the following result:

 

23904i9C25042F065A5ACB

 

As you can see, for negative x values the function now works as expected, but for x values from 0... ~0.4 still NaN are the result...

 

Of course, I could first probe x and then choose one of the two functions, the original or the transform, but this would make things quite slow (still much better than not working).  

 

Thanks again,

 

Wolfgang

0 Kudos
Message 8 of 13
(5,622 Views)

Hi Wolfgang,

 

I guess you were using the transform when x>=0. However, we only use the transform when x is negative to avoid the bug(because x/(x-1) is positive). So the snippet I pasted actually probe x first and then choose between the original and the transform. I am sorry that will reduce the performance. But it is a quick workaround to get the correct result. We will try to investigate if there is other workaround to get both correctness and performance. Thanks!

 

Arthur

0 Kudos
Message 9 of 13
(5,601 Views)

Hi Arthur,

 

sorry for the late reply, I just returned from travel.

 

If I got you right you suggested using the original function for x=0 and positive x values. To me it seems that this approach does not give correct results.

 

For example, y=GaussHG (0.0, 0.5, 2.0, 1.5) returns y=8.948603854199588E-1, while I would have expected the result to be y=1.0, see here.

 

Accordingly, the function is not continuous at x=0, see the figure below:

 

hgf.jpg

 

Did I miss something? The above figure has been calculated according to your code snippet.

 

Thanks,

 

Wolfgang

0 Kudos
Message 10 of 13
(5,473 Views)