gonum.org/v1/gonum@v0.14.0/mathext/internal/gonum/beta.go (about)

     1  // Copyright ©2016 The Gonum Authors. All rights reserved.
     2  // Use of this source code is governed by a BSD-style
     3  // license that can be found in the LICENSE file.
     4  
     5  package gonum
     6  
     7  import (
     8  	"math"
     9  )
    10  
    11  // Beta returns the value of the complete beta function B(a, b). It is defined as
    12  //
    13  //	Γ(a)Γ(b) / Γ(a+b)
    14  //
    15  // Special cases are:
    16  //
    17  //	B(a,b) returns NaN if a or b is Inf
    18  //	B(a,b) returns NaN if a and b are 0
    19  //	B(a,b) returns NaN if a or b is NaN
    20  //	B(a,b) returns NaN if a or b is < 0
    21  //	B(a,b) returns +Inf if a xor b is 0.
    22  //
    23  // See http://mathworld.wolfram.com/BetaFunction.html for more detailed information.
    24  func Beta(a, b float64) float64 {
    25  	return math.Exp(Lbeta(a, b))
    26  }
    27  
    28  // Lbeta returns the natural logarithm of the complete beta function B(a,b).
    29  // Lbeta is defined as:
    30  //
    31  //	Ln(Γ(a)Γ(b)/Γ(a+b))
    32  //
    33  // Special cases are:
    34  //
    35  //	Lbeta(a,b) returns NaN if a or b is Inf
    36  //	Lbeta(a,b) returns NaN if a and b are 0
    37  //	Lbeta(a,b) returns NaN if a or b is NaN
    38  //	Lbeta(a,b) returns NaN if a or b is < 0
    39  //	Lbeta(a,b) returns +Inf if a xor b is 0.
    40  func Lbeta(a, b float64) float64 {
    41  	switch {
    42  	case math.IsInf(a, +1) || math.IsInf(b, +1):
    43  		return math.NaN()
    44  	case a == 0 && b == 0:
    45  		return math.NaN()
    46  	case a < 0 || b < 0:
    47  		return math.NaN()
    48  	case math.IsNaN(a) || math.IsNaN(b):
    49  		return math.NaN()
    50  	case a == 0 || b == 0:
    51  		return math.Inf(+1)
    52  	}
    53  
    54  	la, _ := math.Lgamma(a)
    55  	lb, _ := math.Lgamma(b)
    56  	lab, _ := math.Lgamma(a + b)
    57  	return la + lb - lab
    58  }