bugGNU Octave - Bugs: bug #43551, Matlab incompatibility for gamma ...

 
 

bug #43551: Matlab incompatibility for gamma function for negative integer arguments

Submitted by:  John W. Eaton <jwe>
Submitted on:  Thu 06 Nov 2014 08:07:42 PM UTC  
 
Category: Octave FunctionSeverity: 3 - Normal
Priority: 5 - NormalItem Group: Matlab Compatibility
Status: FixedAssigned to: None
Originator Name: jweOpen/Closed: Closed
Release: devOperating System: Any

Add a New Comment(Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

(Jump to the original submission Jump to the original submission)

Wed 12 Nov 2014 03:54:29 PM UTC, comment #9:

@jwe: I think you can go ahead and apply the patch and close down this report.

Rik <rik5>
Project Administrator
Sun 09 Nov 2014 01:26:16 AM UTC, comment #8:

Or, to quote Abramowitz and Stegun:

Gamma(z) is single valued and analytic over the entire
complex plane, save for the points z = -n (n+0, 1, 2,...)
where it possesses simple poles with residue (-1)^n/n!.

As said before, Inf will do.

Michael Godfrey <godfrey>
Project Member
Fri 07 Nov 2014 09:07:06 PM UTC, comment #7:

The value is undefined at negative integer values in the same way that 1/0 is undefined.

lim 1/x == -Inf
x->0-

lim
x->0+ 1/x == +Inf

The magnitude is infinite, but the phase can't be determined. We could return NaN in this case to be mathematically precise, but nobody bothers to. This suggests that we drop mathematical correctness and just go with the rest of the crowd.

Rik <rik5>
Project Administrator
Fri 07 Nov 2014 08:53:26 PM UTC, comment #6:

John:

The patch looks right to me, assuming that it
uses Inf or the negative poles. This is the usual
assumption in the math/computation (ACM, etc.)
literature and coded algorithms. Not just Matlab.

NaN is not obviously a good substitute for "not defined"
in this case.

Michael Godfrey <godfrey>
Project Member
Fri 07 Nov 2014 08:36:13 PM UTC, comment #5:

I propose the attached patch. Comments?

(file #32391)

John W. Eaton <jwe>
Project Administrator
Fri 07 Nov 2014 05:11:04 PM UTC, comment #4:

Isn't the value really undefined at exact negative integer values? So for the best math, we should have NaN. If Inf is more convenient (but mathematically incorrect regardless of sign) then we might as well just be compatible with Matlab and make it +Inf always.

John W. Eaton <jwe>
Project Administrator
Fri 07 Nov 2014 05:06:29 PM UTC, comment #3:

That's kind of what I expected. It does seem like it would be more useful to return Inf and have terms go to zero.

I was just picking up the patch that someone had supplied so I don't know much about the code.

The mod 2 thing seems to be trying to get the phase to be consistent for negative integers. The result of the gamma function is only infinite when x is a negative integer (or when X itself is +Inf). So the code, by checking xisinf (result), is effectively deciding that x itself was a negative integer. The next test will determine whether x was an even or odd integer and will flip the phase (+1,-1,+1,-1) to go with the absolute magnitude of infinity. This corresponds to consistently taking the limit from the right.

Rik <rik5>
Project Administrator
Fri 07 Nov 2014 04:37:26 PM UTC, comment #2:

The example that resulted in someone pointing out the change to me was one in which gamma(x) was in the denominator of an expression inside a summation. Previously, these terms were then zero and not causing trouble. With the new behavior, they are NaN and cause trouble because then the summation is NaN as well.

If we choose Inf as compatible with Matlab and perhaps more useful, then I think we can still treat -0 specially.

One more question: what is the purpose of the lines

now? Given the current checks at the top of the function, is this condition ever met? And why would we make the result negative only at odd values?

John W. Eaton <jwe>
Project Administrator
Fri 07 Nov 2014 04:25:51 PM UTC, comment #1:

I guess we should think about which form is more useful.

The Wikipedia entry for the gamma function has a nice plot of the function (http://en.wikipedia.org/wiki/Gamma_function). The magnitude of the function at a negative integer is infinite. But, if you choose a consistent direction, such as always approaching the pole from the left, then the sign of the result flips back and forth between +1 and -1. Wolfram Alpha, a.k.a Mathematica, returns complex infinity for gamma (negative integer). A complex infinity has infinite magnitude, but no determined phase. This is consistent with the plot which shows that you are either getting +Inf or -Inf in the limit.

So, is it more useful to simply pick a phase (+1) and return Inf as Matlab does? Or is it more useful to throw a NaN which will clog up downstream calculations and indicate that we have gone beyond our ability to truly calculate the result?

Rik <rik5>
Project Administrator
Thu 06 Nov 2014 08:07:42 PM UTC, original submission:

The change made because of the following patch submission results in Matlab-incompatible behavior.

https://savannah.gnu.org/patch/?func=detailitem&item_id=8101

Apparently Matlab's gamma function returns Inf for negative integer values and for 0 (Matlab has no -0).

Is that important for us?

If we decide to keep the current behavior, then I think we should at least mention that the difference with Matlab is intentional.

John W. Eaton <jwe>
Project Administrator

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach File(s):
   
   
Comment:
   

Attached Files
file #32391:  diffs.txt added by jwe (2KiB - text/plain)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -unavailable- added by godfrey (Posted a comment)
  • -unavailable- added by rik5 (Posted a comment)
  • -unavailable- added by jwe (Submitted the item)
  •  

    Do you think this task is very important?
    If so, you can click here to add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

    Please enter the title of George Orwell's famous dystopian book (it's a date):

     

     

    Follow 3 latest changes.

    Date Changed By Updated Field Previous Value => Replaced By
    Mon 01 Dec 2014 07:20:24 PM UTCrik5StatusNone=>Fixed
      Open/ClosedOpen=>Closed
    Fri 07 Nov 2014 08:36:13 PM UTCjweAttached File-=>Added diffs.txt, #32391

    Back to the top


    Powered by Savane 3.1-cleanup1