Gamma function: Difference between revisions
Content added Content deleted
m (Reduced number of digits in coefficients to reflect float64 math.) |
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
||
Line 21: | Line 21: | ||
{{trans|Python}} |
{{trans|Python}} |
||
< |
<syntaxhighlight lang="11l">V _a = [ 1.00000000000000000000, 0.57721566490153286061, -0.65587807152025388108, |
||
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
||
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
||
Line 40: | Line 40: | ||
L(i) 1..10 |
L(i) 1..10 |
||
print(‘#.14’.format(gamma(i / 3.0)))</ |
print(‘#.14’.format(gamma(i / 3.0)))</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 58: | Line 58: | ||
=={{header|360 Assembly}}== |
=={{header|360 Assembly}}== |
||
For maximum compatibility, this program uses only the basic instruction set. |
For maximum compatibility, this program uses only the basic instruction set. |
||
< |
<syntaxhighlight lang="360asm">GAMMAT CSECT |
||
USING GAMMAT,R13 |
USING GAMMAT,R13 |
||
SAVEAR B STM-SAVEAR(R15) |
SAVEAR B STM-SAVEAR(R15) |
||
Line 216: | Line 216: | ||
EQUREGS |
EQUREGS |
||
EQUREGS REGS=FPR |
EQUREGS REGS=FPR |
||
END GAMMAT</ |
END GAMMAT</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre style="height:20ex"> |
<pre style="height:20ex"> |
||
Line 256: | Line 256: | ||
The coefficients are taken from ''Mathematical functions and their approximations'' |
The coefficients are taken from ''Mathematical functions and their approximations'' |
||
by [[wp:Yudell Luke|Yudell L. Luke]]. |
by [[wp:Yudell Luke|Yudell L. Luke]]. |
||
< |
<syntaxhighlight lang="ada">function Gamma (X : Long_Float) return Long_Float is |
||
A : constant array (0..29) of Long_Float := |
A : constant array (0..29) of Long_Float := |
||
( 1.00000_00000_00000_00000, |
( 1.00000_00000_00000_00000, |
||
Line 296: | Line 296: | ||
end loop; |
end loop; |
||
return 1.0 / Sum; |
return 1.0 / Sum; |
||
end Gamma;</ |
end Gamma;</syntaxhighlight> |
||
Test program: |
Test program: |
||
< |
<syntaxhighlight lang="ada">with Ada.Text_IO; use Ada.Text_IO; |
||
with Gamma; |
with Gamma; |
||
Line 306: | Line 306: | ||
Put_Line (Long_Float'Image (Gamma (Long_Float (I) / 3.0))); |
Put_Line (Long_Float'Image (Gamma (Long_Float (I) / 3.0))); |
||
end loop; |
end loop; |
||
end Test_Gamma;</ |
end Test_Gamma;</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre> |
<pre> |
||
Line 329: | Line 329: | ||
{{works with|ALGOL 68G|Any - tested with release mk15-0.8b.fc9.i386}} |
{{works with|ALGOL 68G|Any - tested with release mk15-0.8b.fc9.i386}} |
||
<!-- {{does not work with|ELLA ALGOL 68|Any (with appropriate job cards) - tested with release 1.8.8d.fc9.i386 - formatted transput is missing, gamma is an extension, and various LONG vs INT operators missing}} --> |
<!-- {{does not work with|ELLA ALGOL 68|Any (with appropriate job cards) - tested with release 1.8.8d.fc9.i386 - formatted transput is missing, gamma is an extension, and various LONG vs INT operators missing}} --> |
||
< |
<syntaxhighlight lang="algol68"># Coefficients used by the GNU Scientific Library # |
||
[]LONG REAL p = ( LONG 0.99999 99999 99809 93, |
[]LONG REAL p = ( LONG 0.99999 99999 99809 93, |
||
LONG 676.52036 81218 851, |
LONG 676.52036 81218 851, |
||
Line 459: | Line 459: | ||
FOR x FROM 10 BY 10 TO 70 DO |
FOR x FROM 10 BY 10 TO 70 DO |
||
printf((sample exp fmt, x, sample(x))) |
printf((sample exp fmt, x, sample(x))) |
||
OD</ |
OD</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 496: | Line 496: | ||
{{trans|BBC Basic}} - Lanczos method. |
{{trans|BBC Basic}} - Lanczos method. |
||
< |
<syntaxhighlight lang="ansi standard basic">100 DECLARE EXTERNAL FUNCTION FNlngamma |
||
110 |
110 |
||
120 DEF FNgamma(z) = EXP(FNlngamma(z)) |
120 DEF FNgamma(z) = EXP(FNlngamma(z)) |
||
Line 521: | Line 521: | ||
330 NEXT i |
330 NEXT i |
||
340 LET FNlngamma = (LOG(SQR(2*PI)) + LOG(a) - b) + LOG(b) * (z+0.5) |
340 LET FNlngamma = (LOG(SQR(2*PI)) + LOG(a) - b) + LOG(b) * (z+0.5) |
||
350 END FUNCTION</ |
350 END FUNCTION</syntaxhighlight> |
||
=={{header|Arturo}}== |
=={{header|Arturo}}== |
||
< |
<syntaxhighlight lang="rebol">A: @[ |
||
1.00000000000000000000 0.57721566490153286061 neg 0.65587807152025388108 |
1.00000000000000000000 0.57721566490153286061 neg 0.65587807152025388108 |
||
neg 0.04200263503409523553 0.16653861138229148950 neg 0.04219773455554433675 |
neg 0.04200263503409523553 0.16653861138229148950 neg 0.04219773455554433675 |
||
Line 556: | Line 556: | ||
pad (to :string v1-v2) 30 |
pad (to :string v1-v2) 30 |
||
] |
] |
||
]</ |
]</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 574: | Line 574: | ||
{{AutoHotkey case}} |
{{AutoHotkey case}} |
||
Source: [http://www.autohotkey.com/forum/topic44657.html AutoHotkey forum] by Laszlo |
Source: [http://www.autohotkey.com/forum/topic44657.html AutoHotkey forum] by Laszlo |
||
< |
<syntaxhighlight lang="autohotkey">/* |
||
Here is the upper incomplete Gamma function. Omitting or setting |
Here is the upper incomplete Gamma function. Omitting or setting |
||
the second parameter to 0 we get the (complete) Gamma function. |
the second parameter to 0 we get the (complete) Gamma function. |
||
Line 655: | Line 655: | ||
1.000000000e+000 1.190639348e+000 1.504575489e+000 2.000000000e+000 2.778158479e+000 |
1.000000000e+000 1.190639348e+000 1.504575489e+000 2.000000000e+000 2.778158479e+000 |
||
1.386831185e+080 1.711224524e+098 8.946182131e+116 1.650795516e+136 9.332621544e+155 |
1.386831185e+080 1.711224524e+098 8.946182131e+116 1.650795516e+136 9.332621544e+155 |
||
*/</ |
*/</syntaxhighlight> |
||
=={{header|AWK}}== |
=={{header|AWK}}== |
||
<syntaxhighlight lang="awk"> |
|||
<lang AWK> |
|||
# syntax: GAWK -f GAMMA_FUNCTION.AWK |
# syntax: GAWK -f GAMMA_FUNCTION.AWK |
||
BEGIN { |
BEGIN { |
||
Line 677: | Line 677: | ||
return exp(b*log(a)) |
return exp(b*log(a)) |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 707: | Line 707: | ||
{{trans|Phix}} - Lanczos method. |
{{trans|Phix}} - Lanczos method. |
||
< |
<syntaxhighlight lang="freebasic">print " x Stirling Lanczos" |
||
print |
print |
||
for i = 1 to 20 |
for i = 1 to 20 |
||
Line 735: | Line 735: | ||
next i |
next i |
||
return sqr(2.0 * pi) * (t ^ (x + 0.5)) * exp(-t) * a |
return sqr(2.0 * pi) * (t ^ (x + 0.5)) * exp(-t) * a |
||
end function</ |
end function</syntaxhighlight> |
||
=={{header|BBC BASIC}}== |
=={{header|BBC BASIC}}== |
||
{{works with|BBC BASIC for Windows}} |
{{works with|BBC BASIC for Windows}} |
||
Uses the Lanczos approximation. |
Uses the Lanczos approximation. |
||
< |
<syntaxhighlight lang="bbcbasic"> *FLOAT64 |
||
INSTALL @lib$+"FNUSING" |
INSTALL @lib$+"FNUSING" |
||
Line 762: | Line 762: | ||
a += lz(i%) / (z + i%) |
a += lz(i%) / (z + i%) |
||
NEXT |
NEXT |
||
= (LNSQR(2*PI) + LN(a) - b) + LN(b) * (z+0.5)</ |
= (LNSQR(2*PI) + LN(a) - b) + LN(b) * (z+0.5)</syntaxhighlight> |
||
'''Output:''' |
'''Output:''' |
||
<pre> |
<pre> |
||
Line 792: | Line 792: | ||
This implements [[wp:Stirling's approximation|Stirling's approximation]] and [[wp:Spouge's approximation|Spouge's approximation]]. |
This implements [[wp:Stirling's approximation|Stirling's approximation]] and [[wp:Spouge's approximation|Spouge's approximation]]. |
||
< |
<syntaxhighlight lang="c">#include <stdio.h> |
||
#include <stdlib.h> |
#include <stdlib.h> |
||
#include <math.h> |
#include <math.h> |
||
Line 843: | Line 843: | ||
} |
} |
||
return 0; |
return 0; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> Stirling Spouge GSL libm |
<pre> Stirling Spouge GSL libm |
||
Line 860: | Line 860: | ||
=={{header|C sharp}}== |
=={{header|C sharp}}== |
||
This is just rewritten from the Wikipedia Lanczos article. Works with complex numbers as well as reals. |
This is just rewritten from the Wikipedia Lanczos article. Works with complex numbers as well as reals. |
||
< |
<syntaxhighlight lang="csharp">using System; |
||
using System.Numerics; |
using System.Numerics; |
||
Line 887: | Line 887: | ||
} |
} |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|C++}}== |
=={{header|C++}}== |
||
< |
<syntaxhighlight lang="cpp">#include <math.h> |
||
#include <numbers> |
#include <numbers> |
||
#include <stdio.h> |
#include <stdio.h> |
||
Line 948: | Line 948: | ||
} |
} |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 974: | Line 974: | ||
=={{header|Clojure}}== |
=={{header|Clojure}}== |
||
< |
<syntaxhighlight lang="clojure">(defn gamma |
||
"Returns Gamma(z + 1 = number) using Lanczos approximation." |
"Returns Gamma(z + 1 = number) using Lanczos approximation." |
||
[number] |
[number] |
||
Line 988: | Line 988: | ||
(Math/exp (- (+ n 7 0.5))) |
(Math/exp (- (+ n 7 0.5))) |
||
(+ (first c) |
(+ (first c) |
||
(apply + (map-indexed #(/ %2 (+ n %1 1)) (next c))))))))</ |
(apply + (map-indexed #(/ %2 (+ n %1 1)) (next c))))))))</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
< |
<syntaxhighlight lang="clojure">(map #(printf "%.1f %.4f\n" % (gamma %)) (map #(float (/ % 10)) (range 1 31)))</syntaxhighlight> |
||
<pre> |
<pre> |
||
0.1 9.5135 |
0.1 9.5135 |
||
Line 1,025: | Line 1,025: | ||
=={{header|Common Lisp}}== |
=={{header|Common Lisp}}== |
||
< |
<syntaxhighlight lang="lisp">; Taylor series coefficients |
||
(defconstant tcoeff |
(defconstant tcoeff |
||
'( 1.00000000000000000000 0.57721566490153286061 -0.65587807152025388108 |
'( 1.00000000000000000000 0.57721566490153286061 -0.65587807152025388108 |
||
Line 1,050: | Line 1,050: | ||
(loop for i from 1 to 10 |
(loop for i from 1 to 10 |
||
do ( |
do ( |
||
format t "~12,10f~%" (gamma (/ i 3.0))))</ |
format t "~12,10f~%" (gamma (/ i 3.0))))</syntaxhighlight> |
||
{{out|Produces}} |
{{out|Produces}} |
||
<pre> |
<pre> |
||
Line 1,068: | Line 1,068: | ||
====Taylor Series | Lanczos Method | Builtin Function==== |
====Taylor Series | Lanczos Method | Builtin Function==== |
||
{{trans|Taylor Series from Ruby; Lanczos Method from C#}} |
{{trans|Taylor Series from Ruby; Lanczos Method from C#}} |
||
< |
<syntaxhighlight lang="ruby"># Taylor Series |
||
def a |
def a |
||
[ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
[ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
||
Line 1,110: | Line 1,110: | ||
puts " Taylor Series Lanczos Method Builtin Function" |
puts " Taylor Series Lanczos Method Builtin Function" |
||
(1..27).each { |i| n = i/3.0; puts "gamma(%.2f) = %.14e %.14e %.14e" % [n, taylor_gamma(n), lanczos_gamma(n), Math.gamma(n)] } |
(1..27).each { |i| n = i/3.0; puts "gamma(%.2f) = %.14e %.14e %.14e" % [n, taylor_gamma(n), lanczos_gamma(n), Math.gamma(n)] } |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,144: | Line 1,144: | ||
=={{header|D}}== |
=={{header|D}}== |
||
< |
<syntaxhighlight lang="d">import std.stdio, std.math, std.mathspecial; |
||
real taylorGamma(in real x) pure nothrow @safe @nogc { |
real taylorGamma(in real x) pure nothrow @safe @nogc { |
||
Line 1,205: | Line 1,205: | ||
x.taylorGamma, x.lanczosGamma, x.gamma); |
x.taylorGamma, x.lanczosGamma, x.gamma); |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>0.333333: 2.6789385347077476335e+00 2.6789385347077470551e+00 2.6789385347077476339e+00 |
<pre>0.333333: 2.6789385347077476335e+00 2.6789385347077470551e+00 2.6789385347077476339e+00 |
||
Line 1,221: | Line 1,221: | ||
{{libheader| System.Math}} |
{{libheader| System.Math}} |
||
{{Trans|Go}} |
{{Trans|Go}} |
||
<syntaxhighlight lang="delphi"> |
|||
<lang Delphi> |
|||
program Gamma_function; |
program Gamma_function; |
||
Line 1,247: | Line 1,247: | ||
writeln(format('%5.1f %24.16g', [x, lanczos7(x)])); |
writeln(format('%5.1f %24.16g', [x, lanczos7(x)])); |
||
readln; |
readln; |
||
end.</ |
end.</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> x Lanczos7 |
<pre> x Lanczos7 |
||
Line 1,263: | Line 1,263: | ||
=={{header|Elixir}}== |
=={{header|Elixir}}== |
||
{{trans|Ruby}} |
{{trans|Ruby}} |
||
< |
<syntaxhighlight lang="elixir">defmodule Gamma do |
||
@a [ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
@a [ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
||
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
||
Line 1,283: | Line 1,283: | ||
Enum.each(Enum.map(1..10, &(&1/3)), fn x -> |
Enum.each(Enum.map(1..10, &(&1/3)), fn x -> |
||
:io.format "~f ~18.15f~n", [x, Gamma.taylor(x)] |
:io.format "~f ~18.15f~n", [x, Gamma.taylor(x)] |
||
end)</ |
end)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 1,303: | Line 1,303: | ||
Solved using the Lanczos Coefficients described in Numerical Recipes (Press et al) |
Solved using the Lanczos Coefficients described in Numerical Recipes (Press et al) |
||
<syntaxhighlight lang="f sharp"> |
|||
<lang F Sharp> |
|||
open System |
open System |
||
Line 1,320: | Line 1,320: | ||
seq { for i in 1 .. 10 do yield ((double)i*10.0) } |> Seq.iter ( fun v -> System.Console.WriteLine("{0} : {1}", v, gamma v ) ) |
seq { for i in 1 .. 10 do yield ((double)i*10.0) } |> Seq.iter ( fun v -> System.Console.WriteLine("{0} : {1}", v, gamma v ) ) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
Line 1,357: | Line 1,357: | ||
=={{header|Factor}}== |
=={{header|Factor}}== |
||
< |
<syntaxhighlight lang="factor">! built in |
||
USING: picomath prettyprint ; |
USING: picomath prettyprint ; |
||
0.1 gamma . ! 9.513507698668723 |
0.1 gamma . ! 9.513507698668723 |
||
2.0 gamma . ! 1.0 |
2.0 gamma . ! 1.0 |
||
10. gamma . ! 362880.0</ |
10. gamma . ! 362880.0</syntaxhighlight> |
||
=={{header|Forth}}== |
=={{header|Forth}}== |
||
Cristinel Mortici describes this method in Applied Mathematics Letters. "A substantial improvement of the Stirling formula". This algorithm is said to give about 7 good digits, but becomes more inaccurate close to zero. Therefore, a "shift" is performed to move the value returned into the more accurate domain. |
Cristinel Mortici describes this method in Applied Mathematics Letters. "A substantial improvement of the Stirling formula". This algorithm is said to give about 7 good digits, but becomes more inaccurate close to zero. Therefore, a "shift" is performed to move the value returned into the more accurate domain. |
||
< |
<syntaxhighlight lang="forth">8 constant (gamma-shift) |
||
: (mortici) ( f1 -- f2) |
: (mortici) ( f1 -- f2) |
||
Line 1,379: | Line 1,379: | ||
fdup (gamma-shift) s>f f+ (mortici) fswap |
fdup (gamma-shift) s>f f+ (mortici) fswap |
||
1 s>f (gamma-shift) 0 do fover i s>f f+ f* loop fswap fdrop f/ |
1 s>f (gamma-shift) 0 do fover i s>f f+ f* loop fswap fdrop f/ |
||
;</ |
;</syntaxhighlight> |
||
<pre> |
<pre> |
||
0.1e gamma f. 9.51348888533932 ok |
0.1e gamma f. 9.51348888533932 ok |
||
Line 1,387: | Line 1,387: | ||
</pre> |
</pre> |
||
This is a word, based on a formula of Ramanujan's famous "lost notebook", which was rediscovered in 1976. His formula contained a constant, which had a value between 1/100 and 1/30. In 2008, E.A. Karatsuba described the function, which determines the value of this constant. Since it contains the gamma function itself, it can't be used in a word calculating the gamma function, so here it is emulated by two symmetrical sigmoidals. |
This is a word, based on a formula of Ramanujan's famous "lost notebook", which was rediscovered in 1976. His formula contained a constant, which had a value between 1/100 and 1/30. In 2008, E.A. Karatsuba described the function, which determines the value of this constant. Since it contains the gamma function itself, it can't be used in a word calculating the gamma function, so here it is emulated by two symmetrical sigmoidals. |
||
< |
<syntaxhighlight lang="forth">2 constant (gamma-shift) \ don't change this |
||
\ an approximation of the d(x) function |
\ an approximation of the d(x) function |
||
: ~d(x) ( f1 -- f2) |
: ~d(x) ( f1 -- f2) |
||
Line 1,411: | Line 1,411: | ||
fdup (gamma-shift) 1- s>f f+ (ramanujan) fswap |
fdup (gamma-shift) 1- s>f f+ (ramanujan) fswap |
||
1 s>f (gamma-shift) 0 do fover i s>f f+ f* loop fswap fdrop f/ |
1 s>f (gamma-shift) 0 do fover i s>f f+ f* loop fswap fdrop f/ |
||
;</ |
;</syntaxhighlight> |
||
<pre> |
<pre> |
||
0.1e gamma f. 9.51351721918848 ok |
0.1e gamma f. 9.51351721918848 ok |
||
Line 1,423: | Line 1,423: | ||
{{works with|Fortran|2008}} |
{{works with|Fortran|2008}} |
||
{{works with|Fortran|95 with extensions}} |
{{works with|Fortran|95 with extensions}} |
||
< |
<syntaxhighlight lang="fortran">program ComputeGammaInt |
||
implicit none |
implicit none |
||
Line 1,520: | Line 1,520: | ||
end function lacz_gamma |
end function lacz_gamma |
||
end program ComputeGammaInt</ |
end program ComputeGammaInt</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> Simpson Lanczos Builtin |
<pre> Simpson Lanczos Builtin |
||
Line 1,536: | Line 1,536: | ||
=={{header|FreeBASIC}}== |
=={{header|FreeBASIC}}== |
||
{{trans|Java}} |
{{trans|Java}} |
||
< |
<syntaxhighlight lang="freebasic">' FB 1.05.0 Win64 |
||
Const pi = 3.1415926535897932 |
Const pi = 3.1415926535897932 |
||
Line 1,582: | Line 1,582: | ||
Print |
Print |
||
Print "Press any key to quit" |
Print "Press any key to quit" |
||
Sleep</ |
Sleep</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 1,611: | Line 1,611: | ||
=={{header|Go}}== |
=={{header|Go}}== |
||
< |
<syntaxhighlight lang="go">package main |
||
import ( |
import ( |
||
Line 1,637: | Line 1,637: | ||
1.5056327351493116e-7/(z+7) |
1.5056327351493116e-7/(z+7) |
||
return math.Sqrt2 * math.SqrtPi * math.Pow(t, z-.5) * math.Exp(-t) * x |
return math.Sqrt2 * math.SqrtPi * math.Pow(t, z-.5) * math.Exp(-t) * x |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,655: | Line 1,655: | ||
=={{header|Groovy}}== |
=={{header|Groovy}}== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="groovy">a = [ 1.00000000000000000000, 0.57721566490153286061, -0.65587807152025388108, |
||
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
||
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
||
Line 1,669: | Line 1,669: | ||
(1..10).each{ printf("% 1.9e\n", gamma(it / 3.0)) } |
(1..10).each{ printf("% 1.9e\n", gamma(it / 3.0)) } |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> 2.678938535e+00 |
<pre> 2.678938535e+00 |
||
Line 1,685: | Line 1,685: | ||
Based on [http://www.haskell.org/haskellwiki/?title=Gamma_and_Beta_function&oldid=25546 HaskellWiki] ([http://www.haskell.org/haskellwiki/HaskellWiki:Copyrights compatible license]): |
Based on [http://www.haskell.org/haskellwiki/?title=Gamma_and_Beta_function&oldid=25546 HaskellWiki] ([http://www.haskell.org/haskellwiki/HaskellWiki:Copyrights compatible license]): |
||
:The Gamma and Beta function as described in 'Numerical Recipes in C++', the approximation is taken from [Lanczos, C. 1964 SIAM Journal on Numerical Analysis, ser. B, vol. 1, pp. 86-96] |
:The Gamma and Beta function as described in 'Numerical Recipes in C++', the approximation is taken from [Lanczos, C. 1964 SIAM Journal on Numerical Analysis, ser. B, vol. 1, pp. 86-96] |
||
< |
<syntaxhighlight lang="haskell">cof :: [Double] |
||
cof = |
cof = |
||
[ 76.18009172947146 |
[ 76.18009172947146 |
||
Line 1,705: | Line 1,705: | ||
main :: IO () |
main :: IO () |
||
main = mapM_ print $ gammaln <$> [0.1,0.2 .. 1.0]</ |
main = mapM_ print $ gammaln <$> [0.1,0.2 .. 1.0]</syntaxhighlight> |
||
Or equivalently, as a point-free applicative expression: |
Or equivalently, as a point-free applicative expression: |
||
< |
<syntaxhighlight lang="haskell">import Control.Applicative |
||
cof :: [Double] |
cof :: [Double] |
||
Line 1,729: | Line 1,729: | ||
main :: IO () |
main :: IO () |
||
main = mapM_ print $ gammaln <$> [0.1,0.2 .. 1.0]</ |
main = mapM_ print $ gammaln <$> [0.1,0.2 .. 1.0]</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre>2.252712651734255 |
<pre>2.252712651734255 |
||
Line 1,746: | Line 1,746: | ||
This works in Unicon. Changing the <tt>!10</tt> into <tt>(1 to 10)</tt> would enable it |
This works in Unicon. Changing the <tt>!10</tt> into <tt>(1 to 10)</tt> would enable it |
||
to work in Icon. |
to work in Icon. |
||
< |
<syntaxhighlight lang="unicon">procedure main() |
||
every write(left(i := !10/10.0,5),gamma(.i)) |
every write(left(i := !10/10.0,5),gamma(.i)) |
||
end |
end |
||
Line 1,752: | Line 1,752: | ||
procedure gamma(z) # Stirling's approximation |
procedure gamma(z) # Stirling's approximation |
||
return (2*&pi/z)^0.5 * (z/&e)^z |
return (2*&pi/z)^0.5 * (z/&e)^z |
||
end</ |
end</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
Line 1,772: | Line 1,772: | ||
=={{header|J}}== |
=={{header|J}}== |
||
This code shows the built-in method, which works for any value (positive, negative and complex numbers -- but note that negative integer arguments give infinite results). |
This code shows the built-in method, which works for any value (positive, negative and complex numbers -- but note that negative integer arguments give infinite results). |
||
< |
<syntaxhighlight lang="j">gamma=: !@<:</syntaxhighlight> |
||
<code><:</code> subtracts one from a number. It's sort of like <code>--lvalue</code> in C, except it always accepts an "rvalue" as an argument (which means it does not modify that argument). And <code>!value</code> finds the factorial of value if value is a positive integer. This illustrates the close relationship between the factorial and gamma functions. |
<code><:</code> subtracts one from a number. It's sort of like <code>--lvalue</code> in C, except it always accepts an "rvalue" as an argument (which means it does not modify that argument). And <code>!value</code> finds the factorial of value if value is a positive integer. This illustrates the close relationship between the factorial and gamma functions. |
||
The following direct coding of the task comes from the [[J:Essays/Stirling's%20Approximation|Stirling's approximation essay]] on the J wiki: |
The following direct coding of the task comes from the [[J:Essays/Stirling's%20Approximation|Stirling's approximation essay]] on the J wiki: |
||
< |
<syntaxhighlight lang="j">sbase =: %:@(2p1&%) * %&1x1 ^ ] |
||
scorr =: 1 1r12 1r288 _139r51840 _571r2488320&p.@% |
scorr =: 1 1r12 1r288 _139r51840 _571r2488320&p.@% |
||
stirlg=: sbase * scorr</ |
stirlg=: sbase * scorr</syntaxhighlight> |
||
Checking against <code>!@<:</code> we can see that this approximation loses accuracy for small arguments |
Checking against <code>!@<:</code> we can see that this approximation loses accuracy for small arguments |
||
< |
<syntaxhighlight lang="j"> (,. stirlg ,. gamma) 10 1p1 1x1 1.5 1 |
||
10 362880 362880 |
10 362880 362880 |
||
3.14159 2.28803 2.28804 |
3.14159 2.28803 2.28804 |
||
2.71828 1.56746 1.56747 |
2.71828 1.56746 1.56747 |
||
1.5 0.886155 0.886227 |
1.5 0.886155 0.886227 |
||
1 0.999499 1</ |
1 0.999499 1</syntaxhighlight> |
||
(Column 1 is the argument, column 2 is the stirling approximation and column 3 uses the builtin support for gamma.) |
(Column 1 is the argument, column 2 is the stirling approximation and column 3 uses the builtin support for gamma.) |
||
=={{header|Java}}== |
=={{header|Java}}== |
||
Implementation of Stirling's approximation and Lanczos approximation. |
Implementation of Stirling's approximation and Lanczos approximation. |
||
< |
<syntaxhighlight lang="java">public class GammaFunction { |
||
public double st_gamma(double x){ |
public double st_gamma(double x){ |
||
Line 1,820: | Line 1,820: | ||
} |
} |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,848: | Line 1,848: | ||
=={{header|JavaScript}}== |
=={{header|JavaScript}}== |
||
Implementation of Lanczos approximation. |
Implementation of Lanczos approximation. |
||
< |
<syntaxhighlight lang="javascript">function gamma(x) { |
||
var p = [0.99999999999980993, 676.5203681218851, -1259.1392167224028, |
var p = [0.99999999999980993, 676.5203681218851, -1259.1392167224028, |
||
771.32342877765313, -176.61502916214059, 12.507343278686905, |
771.32342877765313, -176.61502916214059, 12.507343278686905, |
||
Line 1,867: | Line 1,867: | ||
return Math.sqrt(2 * Math.PI) * Math.pow(t, x + 0.5) * Math.exp(-t) * a; |
return Math.sqrt(2 * Math.PI) * Math.pow(t, x + 0.5) * Math.exp(-t) * a; |
||
}</ |
}</syntaxhighlight> |
||
=={{header|jq}}== |
=={{header|jq}}== |
||
Line 1,873: | Line 1,873: | ||
====Taylor Series==== |
====Taylor Series==== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="jq">def gamma: |
||
[ |
[ |
||
1.00000000000000000000, 0.57721566490153286061, -0.65587807152025388108, -0.04200263503409523553, |
1.00000000000000000000, 0.57721566490153286061, -0.65587807152025388108, -0.04200263503409523553, |
||
Line 1,888: | Line 1,888: | ||
| reduce range(2; 1+$n) as $an |
| reduce range(2; 1+$n) as $an |
||
($a[$n-1]; (. * $y) + $a[$n - $an]) |
($a[$n-1]; (. * $y) + $a[$n - $an]) |
||
| 1 / . ;</ |
| 1 / . ;</syntaxhighlight> |
||
====Lanczos Approximation==== |
====Lanczos Approximation==== |
||
< |
<syntaxhighlight lang="jq"># for reals, but easily extended to complex values |
||
def gamma_by_lanczos: |
def gamma_by_lanczos: |
||
def pow(x): if x == 0 then 1 elif x == 1 then . else x * log | exp end; |
def pow(x): if x == 0 then 1 elif x == 1 then . else x * log | exp end; |
||
Line 1,905: | Line 1,905: | ||
($p[0]; . + ($p[$i]/($x + $i) )) |
($p[0]; . + ($p[$i]/($x + $i) )) |
||
* ((2*$pi) | sqrt) * ($t | pow($x+0.5)) * ((-$t)|exp) |
* ((2*$pi) | sqrt) * ($t | pow($x+0.5)) * ((-$t)|exp) |
||
end;</ |
end;</syntaxhighlight> |
||
====Stirling's Approximation==== |
====Stirling's Approximation==== |
||
< |
<syntaxhighlight lang="jq">def gamma_by_stirling: |
||
def pow(x): if x == 0 then 1 elif x == 1 then . else x * log | exp end; |
def pow(x): if x == 0 then 1 elif x == 1 then . else x * log | exp end; |
||
((1|atan) * 8) as $twopi |
((1|atan) * 8) as $twopi |
||
| . as $x |
| . as $x |
||
| (($twopi/$x) | sqrt) * ( ($x / (1|exp)) | pow($x));</ |
| (($twopi/$x) | sqrt) * ( ($x / (1|exp)) | pow($x));</syntaxhighlight> |
||
====Examples==== |
====Examples==== |
||
Stirling's method produces poor results, so to save space, the examples |
Stirling's method produces poor results, so to save space, the examples |
||
contrast the Taylor series and Lanczos methods with built-in tgamma: |
contrast the Taylor series and Lanczos methods with built-in tgamma: |
||
< |
<syntaxhighlight lang="jq">def pad(n): tostring | . + (n - length) * " "; |
||
" i: gamma lanczos tgamma", |
" i: gamma lanczos tgamma", |
||
(range(1;11) |
(range(1;11) |
||
| . / 3.0 |
| . / 3.0 |
||
| "\(pad(18)): \(gamma|pad(18)) : \(gamma_by_lanczos|pad(18)) : \(tgamma)")</ |
| "\(pad(18)): \(gamma|pad(18)) : \(gamma_by_lanczos|pad(18)) : \(tgamma)")</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
< |
<syntaxhighlight lang="sh">$ jq -M -r -n -f Gamma_function_Stirling.jq |
||
i: gamma lanczos tgamma |
i: gamma lanczos tgamma |
||
0.3333333333333333: 2.6789385347077483 : 2.6789385347077483 : 2.678938534707748 |
0.3333333333333333: 2.6789385347077483 : 2.6789385347077483 : 2.678938534707748 |
||
Line 1,933: | Line 1,933: | ||
2.6666666666666665: 1.5045754882515399 : 1.5045754882515576 : 1.5045754882515558 |
2.6666666666666665: 1.5045754882515399 : 1.5045754882515576 : 1.5045754882515558 |
||
3 : 1.9999999999939684 : 2.0000000000000013 : 2 |
3 : 1.9999999999939684 : 2.0000000000000013 : 2 |
||
3.3333333333333335: 2.778158479338573 : 2.778158480437665 : 2.7781584804376647</ |
3.3333333333333335: 2.778158479338573 : 2.778158480437665 : 2.7781584804376647</syntaxhighlight> |
||
=={{header|Jsish}}== |
=={{header|Jsish}}== |
||
{{trans|Javascript}} |
{{trans|Javascript}} |
||
< |
<syntaxhighlight lang="javascript">#!/usr/bin/env jsish |
||
/* Gamma function, in Jsish, using the Lanczos approximation */ |
/* Gamma function, in Jsish, using the Lanczos approximation */ |
||
function gamma(x) { |
function gamma(x) { |
||
Line 1,992: | Line 1,992: | ||
5.5 +5.234278e+01 |
5.5 +5.234278e+01 |
||
=!EXPECTEND!= |
=!EXPECTEND!= |
||
*/</ |
*/</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,027: | Line 2,027: | ||
'''Built-in function''': |
'''Built-in function''': |
||
<lang |
<syntaxhighlight lang="julia">@show gamma(1)</syntaxhighlight> |
||
'''By adaptive Gauss-Kronrod integration''': |
'''By adaptive Gauss-Kronrod integration''': |
||
< |
<syntaxhighlight lang="julia">using QuadGK |
||
gammaquad(t::Float64) = first(quadgk(x -> x ^ (t - 1) * exp(-x), zero(t), Inf, reltol = 100eps(t))) |
gammaquad(t::Float64) = first(quadgk(x -> x ^ (t - 1) * exp(-x), zero(t), Inf, reltol = 100eps(t))) |
||
@show gammaquad(1.0)</ |
@show gammaquad(1.0)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,040: | Line 2,040: | ||
{{works with|Julia|1.0}} |
{{works with|Julia|1.0}} |
||
'''Library function''': |
'''Library function''': |
||
< |
<syntaxhighlight lang="julia">using SpecialFunctions |
||
gamma(1/2) - sqrt(pi)</ |
gamma(1/2) - sqrt(pi)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,047: | Line 2,047: | ||
=={{header|Kotlin}}== |
=={{header|Kotlin}}== |
||
< |
<syntaxhighlight lang="scala">// version 1.0.6 |
||
fun gammaStirling(x: Double): Double = Math.sqrt(2.0 * Math.PI / x) * Math.pow(x / Math.E, x) |
fun gammaStirling(x: Double): Double = Math.sqrt(2.0 * Math.PI / x) * Math.pow(x / Math.E, x) |
||
Line 2,081: | Line 2,081: | ||
println("%17.15f".format(gammaLanczos(d))) |
println("%17.15f".format(gammaLanczos(d))) |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,114: | Line 2,114: | ||
A fairly straightforward port of the Go code. (It could almost have been done with sed). A few small differences are in the use of a tuple as a return value for the builtin gamma function, and we import a few functions from the math library so that we don't have to qualify them. |
A fairly straightforward port of the Go code. (It could almost have been done with sed). A few small differences are in the use of a tuple as a return value for the builtin gamma function, and we import a few functions from the math library so that we don't have to qualify them. |
||
< |
<syntaxhighlight lang="limbo">implement Lanczos7; |
||
include "sys.m"; sys: Sys; |
include "sys.m"; sys: Sys; |
||
Line 2,158: | Line 2,158: | ||
return sqrt(2.0) * sqrt(Math->Pi) * pow(t, z - 0.5) * exp(-t) * x; |
return sqrt(2.0) * sqrt(Math->Pi) * pow(t, z - 0.5) * exp(-t) * x; |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{output}} |
{{output}} |
||
Line 2,176: | Line 2,176: | ||
=={{header|Lua}}== |
=={{header|Lua}}== |
||
Uses the [[wp:Reciprocal gamma function]] to calculate small values. |
Uses the [[wp:Reciprocal gamma function]] to calculate small values. |
||
< |
<syntaxhighlight lang="lua">gamma, coeff, quad, qui, set = 0.577215664901, -0.65587807152056, -0.042002635033944, 0.16653861138228, -0.042197734555571 |
||
function recigamma(z) |
function recigamma(z) |
||
return z + gamma * z^2 + coeff * z^3 + quad * z^4 + qui * z^5 + set * z^6 |
return z + gamma * z^2 + coeff * z^3 + quad * z^4 + qui * z^5 + set * z^6 |
||
Line 2,186: | Line 2,186: | ||
else return (z - 1) * gammafunc(z-1) |
else return (z - 1) * gammafunc(z-1) |
||
end |
end |
||
end</ |
end</syntaxhighlight> |
||
=={{header|M2000 Interpreter}}== |
=={{header|M2000 Interpreter}}== |
||
<syntaxhighlight lang="m2000 interpreter"> |
|||
<lang M2000 Interpreter> |
|||
Module PrepareLambdaFunctions { |
Module PrepareLambdaFunctions { |
||
Const e = 2.7182818284590452@ |
Const e = 2.7182818284590452@ |
||
Line 2,229: | Line 2,229: | ||
Print $("") |
Print $("") |
||
clipboard doc$ |
clipboard doc$ |
||
</syntaxhighlight> |
|||
</lang> |
|||
χ Stirling Lanczos |
χ Stirling Lanczos |
||
0.10 5.697187148977170 9.5135076986687024462927178610 |
0.10 5.697187148977170 9.5135076986687024462927178610 |
||
Line 2,254: | Line 2,254: | ||
=={{header|Maple}}== |
=={{header|Maple}}== |
||
Built-in method that accepts any value. |
Built-in method that accepts any value. |
||
< |
<syntaxhighlight lang="maple">GAMMA(17/2); |
||
GAMMA(7*I); |
GAMMA(7*I); |
||
M := Matrix(2, 3, 'fill' = -3.6); |
M := Matrix(2, 3, 'fill' = -3.6); |
||
MTM:-gamma(M);</ |
MTM:-gamma(M);</syntaxhighlight> |
||
{{Out|Output}} |
{{Out|Output}} |
||
<pre>2027025*sqrt(Pi)*(1/256) |
<pre>2027025*sqrt(Pi)*(1/256) |
||
Line 2,265: | Line 2,265: | ||
=={{header|Mathematica}}/{{header|Wolfram Language}}== |
=={{header|Mathematica}}/{{header|Wolfram Language}}== |
||
This code shows the built-in method, which works for any value (positive, negative and complex numbers). |
This code shows the built-in method, which works for any value (positive, negative and complex numbers). |
||
<lang |
<syntaxhighlight lang="mathematica">Gamma[x]</syntaxhighlight> |
||
Output integers and half-integers (a space is multiplication in Mathematica): |
Output integers and half-integers (a space is multiplication in Mathematica): |
||
<pre> |
<pre> |
||
Line 2,312: | Line 2,312: | ||
=={{header|Maxima}}== |
=={{header|Maxima}}== |
||
< |
<syntaxhighlight lang="maxima">fpprec: 30$ |
||
gamma_coeff(n) := block([a: makelist(1, n)], |
gamma_coeff(n) := block([a: makelist(1, n)], |
||
Line 2,332: | Line 2,332: | ||
gamma_approx(12.3b0) - gamma(12.3b0); |
gamma_approx(12.3b0) - gamma(12.3b0); |
||
/* -9.25224705314470500985141176997b-15 */</ |
/* -9.25224705314470500985141176997b-15 */</syntaxhighlight> |
||
=={{header|МК-61/52}}== |
=={{header|МК-61/52}}== |
||
Line 2,345: | Line 2,345: | ||
=={{header|Modula-3}}== |
=={{header|Modula-3}}== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="modula3">MODULE Gamma EXPORTS Main; |
||
FROM IO IMPORT Put; |
FROM IO IMPORT Put; |
||
Line 2,381: | Line 2,381: | ||
Put(Extended(Taylor(FLOAT(i, EXTENDED) / 3.0X0), style := Style.Sci) & "\n"); |
Put(Extended(Taylor(FLOAT(i, EXTENDED) / 3.0X0), style := Style.Sci) & "\n"); |
||
END; |
END; |
||
END Gamma.</ |
END Gamma.</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,399: | Line 2,399: | ||
{{trans|Ada}} |
{{trans|Ada}} |
||
The algorithm is a translation of that from the Ada solution. We have added a comparison with the gamma function provided by the “math” module from Nim standard library (which is, in fact, the C “tgamma” function). |
The algorithm is a translation of that from the Ada solution. We have added a comparison with the gamma function provided by the “math” module from Nim standard library (which is, in fact, the C “tgamma” function). |
||
< |
<syntaxhighlight lang="nim">import math, strformat |
||
const A = [ |
const A = [ |
||
Line 2,425: | Line 2,425: | ||
let val1 = gamma(i.toFloat / 3) |
let val1 = gamma(i.toFloat / 3) |
||
let val2 = math.gamma(i.toFloat / 3) |
let val2 = math.gamma(i.toFloat / 3) |
||
echo &"{val1:18.16f} {val2:18.16f} {val1 - val2:11.4e}"</ |
echo &"{val1:18.16f} {val2:18.16f} {val1 - val2:11.4e}"</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
Line 2,442: | Line 2,442: | ||
=={{header|OCaml}}== |
=={{header|OCaml}}== |
||
< |
<syntaxhighlight lang="ocaml">let e = exp 1. |
||
let pi = 4. *. atan 1. |
let pi = 4. *. atan 1. |
||
let sqrttwopi = sqrt (2. *. pi) |
let sqrttwopi = sqrt (2. *. pi) |
||
Line 2,504: | Line 2,504: | ||
(Stirling.gamma z) |
(Stirling.gamma z) |
||
(Stirling2.gamma z) |
(Stirling2.gamma z) |
||
done</ |
done</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,538: | Line 2,538: | ||
=={{header|Octave}}== |
=={{header|Octave}}== |
||
< |
<syntaxhighlight lang="octave">function g = lacz_gamma(a, cg=7) |
||
p = [ 0.99999999999980993, 676.5203681218851, -1259.1392167224028, \ |
p = [ 0.99999999999980993, 676.5203681218851, -1259.1392167224028, \ |
||
771.32342877765313, -176.61502916214059, 12.507343278686905, \ |
771.32342877765313, -176.61502916214059, 12.507343278686905, \ |
||
Line 2,559: | Line 2,559: | ||
for i = 1:10 |
for i = 1:10 |
||
printf("%f %f\n", gamma(i/3.0), lacz_gamma(i/3.0)); |
printf("%f %f\n", gamma(i/3.0), lacz_gamma(i/3.0)); |
||
endfor</ |
endfor</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>2.678939 2.678939 |
<pre>2.678939 2.678939 |
||
Line 2,575: | Line 2,575: | ||
=={{header|Oforth}}== |
=={{header|Oforth}}== |
||
< |
<syntaxhighlight lang="oforth">import: math |
||
[ |
[ |
||
Line 2,591: | Line 2,591: | ||
2 Pi * sqrt * |
2 Pi * sqrt * |
||
t z 0.5 + powf * |
t z 0.5 + powf * |
||
t neg exp * ;</ |
t neg exp * ;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,618: | Line 2,618: | ||
=={{header|PARI/GP}}== |
=={{header|PARI/GP}}== |
||
===Built-in=== |
===Built-in=== |
||
<lang |
<syntaxhighlight lang="parigp">gamma(x)</syntaxhighlight> |
||
===Double-exponential integration=== |
===Double-exponential integration=== |
||
<code>[[+oo],k]</code> means that the function approaches <math>+\infty</math> as <math>\exp(-kx).</math> |
<code>[[+oo],k]</code> means that the function approaches <math>+\infty</math> as <math>\exp(-kx).</math> |
||
< |
<syntaxhighlight lang="parigp">Gamma(x)=intnum(t=0,[+oo,1],t^(x-1)/exp(t))</syntaxhighlight> |
||
===Romberg integration=== |
===Romberg integration=== |
||
< |
<syntaxhighlight lang="parigp">Gamma(x)=intnumromb(t=0,9,t^(x-1)/exp(t),0)+intnumromb(t=9,max(x,99)^9,t^(x-1)/exp(t),2)</syntaxhighlight> |
||
===Stirling approximation=== |
===Stirling approximation=== |
||
< |
<syntaxhighlight lang="parigp">Stirling(x)=x--;sqrt(2*Pi*x)*(x/exp(1))^x</syntaxhighlight> |
||
=={{header|Pascal}}== |
=={{header|Pascal}}== |
||
Line 2,636: | Line 2,636: | ||
Gamma(x)*Gamma(1 - x) = pi/sin(pi*x). |
Gamma(x)*Gamma(1 - x) = pi/sin(pi*x). |
||
< |
<syntaxhighlight lang="pascal"> |
||
program GammaTest; |
program GammaTest; |
||
{$mode objfpc}{$H+} |
{$mode objfpc}{$H+} |
||
Line 2,713: | Line 2,713: | ||
WriteLn( SysUtils.Format( '-(32/945)Sqrt(pi) = %g', [-k*Sqrt(PI)])); |
WriteLn( SysUtils.Format( '-(32/945)Sqrt(pi) = %g', [-k*Sqrt(PI)])); |
||
end. |
end. |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,736: | Line 2,736: | ||
=={{header|Perl}}== |
=={{header|Perl}}== |
||
< |
<syntaxhighlight lang="perl">use strict; |
||
use warnings; |
use warnings; |
||
use constant pi => 4*atan2(1, 1); |
use constant pi => 4*atan2(1, 1); |
||
Line 2,805: | Line 2,805: | ||
print join(' ', map { sprintf "%.12f", Gamma($_/3, $method) } 1 .. 10); |
print join(' ', map { sprintf "%.12f", Gamma($_/3, $method) } 1 .. 10); |
||
print "\n"; |
print "\n"; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> MPFR: 2.678938534708 1.354117939426 1.000000000000 0.892979511569 0.902745292951 1.000000000000 1.190639348759 1.504575488252 2.000000000000 2.778158480438 |
<pre> MPFR: 2.678938534708 1.354117939426 1.000000000000 0.892979511569 0.902745292951 1.000000000000 1.190639348759 1.504575488252 2.000000000000 2.778158480438 |
||
Line 2,814: | Line 2,814: | ||
=={{header|Phix}}== |
=={{header|Phix}}== |
||
{{trans|C}} |
{{trans|C}} |
||
<!--< |
<!--<syntaxhighlight lang="phix">(phixonline)--> |
||
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">c</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">12</span><span style="color: #0000FF;">)</span> |
<span style="color: #004080;">sequence</span> <span style="color: #000000;">c</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">12</span><span style="color: #0000FF;">)</span> |
||
Line 2,855: | Line 2,855: | ||
<span style="color: #000000;">sq</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma</span><span style="color: #0000FF;">(</span><span style="color: #000000;">7</span><span style="color: #0000FF;">/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">),</span><span style="color: #000000;">8</span><span style="color: #0000FF;">/</span><span style="color: #000000;">15</span><span style="color: #0000FF;">)</span> |
<span style="color: #000000;">sq</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma</span><span style="color: #0000FF;">(</span><span style="color: #000000;">7</span><span style="color: #0000FF;">/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">),</span><span style="color: #000000;">8</span><span style="color: #0000FF;">/</span><span style="color: #000000;">15</span><span style="color: #0000FF;">)</span> |
||
<span style="color: #000000;">si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma</span><span style="color: #0000FF;">(</span><span style="color: #000000;">4</span><span style="color: #0000FF;">))</span> |
<span style="color: #000000;">si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma</span><span style="color: #0000FF;">(</span><span style="color: #000000;">4</span><span style="color: #0000FF;">))</span> |
||
<!--</ |
<!--</syntaxhighlight>--> |
||
{{out}} |
{{out}} |
||
<small>(slightly different results under pwa/p2js...)</small> |
<small>(slightly different results under pwa/p2js...)</small> |
||
Line 2,873: | Line 2,873: | ||
Above translated to mpfr, with higher accuracy and more iterations as per REXX, and compared against the builtin. |
Above translated to mpfr, with higher accuracy and more iterations as per REXX, and compared against the builtin. |
||
{{libheader|Phix/mpfr}} |
{{libheader|Phix/mpfr}} |
||
<!--< |
<!--<syntaxhighlight lang="phix">--> |
||
<span style="color: #008080;">without</span> <span style="color: #008080;">javascript_semantics</span> <span style="color: #000080;font-style:italic;">-- (no mpfr_exp(), mpfr_gamma() in pwa/p2js)</span> |
<span style="color: #008080;">without</span> <span style="color: #008080;">javascript_semantics</span> <span style="color: #000080;font-style:italic;">-- (no mpfr_exp(), mpfr_gamma() in pwa/p2js)</span> |
||
<span style="color: #7060A8;">requires</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"1.0.0"</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (mpfr_set_default_prec[ision] has been renamed)</span> |
<span style="color: #7060A8;">requires</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"1.0.0"</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (mpfr_set_default_prec[ision] has been renamed)</span> |
||
Line 2,969: | Line 2,969: | ||
<span style="color: #000000;">sq</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma2</span><span style="color: #0000FF;">(</span><span style="color: #000000;">7</span><span style="color: #0000FF;">/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">),</span><span style="color: #000000;">8</span><span style="color: #0000FF;">,</span><span style="color: #000000;">15</span><span style="color: #0000FF;">)</span> |
<span style="color: #000000;">sq</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma2</span><span style="color: #0000FF;">(</span><span style="color: #000000;">7</span><span style="color: #0000FF;">/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">),</span><span style="color: #000000;">8</span><span style="color: #0000FF;">,</span><span style="color: #000000;">15</span><span style="color: #0000FF;">)</span> |
||
<span style="color: #000000;">si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma2</span><span style="color: #0000FF;">(</span><span style="color: #000000;">4</span><span style="color: #0000FF;">))</span> |
<span style="color: #000000;">si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gamma2</span><span style="color: #0000FF;">(</span><span style="color: #000000;">4</span><span style="color: #0000FF;">))</span> |
||
<!--</ |
<!--</syntaxhighlight>--> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,996: | Line 2,996: | ||
=={{header|Phixmonti}}== |
=={{header|Phixmonti}}== |
||
< |
<syntaxhighlight lang="phixmonti">0.577215664901 var gamma |
||
-0.65587807152056 var coeff |
-0.65587807152056 var coeff |
||
-0.042002635033944 var quad |
-0.042002635033944 var quad |
||
Line 3,036: | Line 3,036: | ||
for |
for |
||
dup print " = " print gammafunc print nl |
dup print " = " print gammafunc print nl |
||
endfor</ |
endfor</syntaxhighlight> |
||
=={{header|PicoLisp}}== |
=={{header|PicoLisp}}== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="picolisp">(scl 28) |
||
(de *A |
(de *A |
||
Line 3,059: | Line 3,059: | ||
(for A (cdr *A) |
(for A (cdr *A) |
||
(setq Sum (+ A (*/ Sum Y 1.0))) ) |
(setq Sum (+ A (*/ Sum Y 1.0))) ) |
||
(*/ 1.0 1.0 Sum) ) )</ |
(*/ 1.0 1.0 Sum) ) )</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>: (for I (range 1 10) |
<pre>: (for I (range 1 10) |
||
Line 3,075: | Line 3,075: | ||
=={{header|PL/I}}== |
=={{header|PL/I}}== |
||
< |
<syntaxhighlight lang="pl/i">/* From Rosetta Fortran */ |
||
test: procedure options (main); |
test: procedure options (main); |
||
Line 3,119: | Line 3,119: | ||
end lanczos_gamma; |
end lanczos_gamma; |
||
end test;</ |
end test;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,137: | Line 3,137: | ||
=={{header|PowerShell}}== |
=={{header|PowerShell}}== |
||
I would download the Math.NET Numerics dll(s). Documentation and download at: http://cyber-defense.sans.org/blog/2015/06/27/powershell-for-math-net-numerics/comment-page-1/ |
I would download the Math.NET Numerics dll(s). Documentation and download at: http://cyber-defense.sans.org/blog/2015/06/27/powershell-for-math-net-numerics/comment-page-1/ |
||
<syntaxhighlight lang="powershell"> |
|||
<lang PowerShell> |
|||
Add-Type -Path "C:\Program Files (x86)\Math\MathNet.Numerics.3.12.0\lib\net40\MathNet.Numerics.dll" |
Add-Type -Path "C:\Program Files (x86)\Math\MathNet.Numerics.3.12.0\lib\net40\MathNet.Numerics.dll" |
||
1..20 | ForEach-Object {[MathNet.Numerics.SpecialFunctions]::Gamma($_ / 10)} |
1..20 | ForEach-Object {[MathNet.Numerics.SpecialFunctions]::Gamma($_ / 10)} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{Out}} |
{{Out}} |
||
Line 3,169: | Line 3,169: | ||
=={{header|Prolog}}== |
=={{header|Prolog}}== |
||
This version matches Wolfram Alpha to within a few digits at the end, so the last few digits are a bit off. There's an early check to stop evaluating coefficients once the desired accuracy is reached. Removing this check does not improve accuracy vs. Wolfram Alpha. |
This version matches Wolfram Alpha to within a few digits at the end, so the last few digits are a bit off. There's an early check to stop evaluating coefficients once the desired accuracy is reached. Removing this check does not improve accuracy vs. Wolfram Alpha. |
||
<syntaxhighlight lang="prolog"> |
|||
<lang Prolog> |
|||
gamma_coefficients( |
gamma_coefficients( |
||
[ 1.00000000000000000000000, 0.57721566490153286060651, -0.65587807152025388107701, |
[ 1.00000000000000000000000, 0.57721566490153286060651, -0.65587807152025388107701, |
||
Line 3,206: | Line 3,206: | ||
X2 is X1 + C*Power, |
X2 is X1 + C*Power, |
||
recip_gamma(X, Power, Cs, X1, X2, Y). |
recip_gamma(X, Power, Cs, X1, X2, Y). |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{Out}} |
{{Out}} |
||
<pre> |
<pre> |
||
Line 3,253: | Line 3,253: | ||
* Natural Logarithm of the Complete Gamma function |
* Natural Logarithm of the Complete Gamma function |
||
* Factorial function |
* Factorial function |
||
< |
<syntaxhighlight lang="purebasic">Procedure.d Gamma(x.d) ; AKJ 01-May-10 |
||
; Complete Gamma function for x>0 and x<2 (approx) |
; Complete Gamma function for x>0 and x<2 (approx) |
||
; Extended outside this range via: Gamma(x+1) = x*Gamma(x) |
; Extended outside this range via: Gamma(x+1) = x*Gamma(x) |
||
Line 3,334: | Line 3,334: | ||
Procedure Factorial(x) ; AKJ 01-May-10 |
Procedure Factorial(x) ; AKJ 01-May-10 |
||
ProcedureReturn Gamma(x+1) |
ProcedureReturn Gamma(x+1) |
||
EndProcedure</ |
EndProcedure</syntaxhighlight> |
||
;Examples |
;Examples |
||
< |
<syntaxhighlight lang="purebasic">Debug "Gamma()" |
||
For i = 12 To 15 |
For i = 12 To 15 |
||
Debug StrD(i/3.0, 3)+" "+StrD(Gamma(i/3.0)) |
Debug StrD(i/3.0, 3)+" "+StrD(Gamma(i/3.0)) |
||
Line 3,343: | Line 3,343: | ||
Debug "Ln(Gamma(5.0)) = "+StrD(GammLn(5.0), 16) ; Ln(24) |
Debug "Ln(Gamma(5.0)) = "+StrD(GammLn(5.0), 16) ; Ln(24) |
||
Debug "" |
Debug "" |
||
Debug "Factorial 6 = "+StrD(Factorial(6), 0) ; 72</ |
Debug "Factorial 6 = "+StrD(Factorial(6), 0) ; 72</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>[Debug] Gamma(): |
<pre>[Debug] Gamma(): |
||
Line 3,358: | Line 3,358: | ||
===Procedural=== |
===Procedural=== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="python">_a = ( 1.00000000000000000000, 0.57721566490153286061, -0.65587807152025388108, |
||
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
-0.04200263503409523553, 0.16653861138229148950, -0.04219773455554433675, |
||
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
-0.00962197152787697356, 0.00721894324666309954, -0.00116516759185906511, |
||
Line 3,380: | Line 3,380: | ||
for i in range(1,11): |
for i in range(1,11): |
||
print " %20.14e" % gamma(i/3.0) |
print " %20.14e" % gamma(i/3.0) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> 2.67893853470775e+00 |
<pre> 2.67893853470775e+00 |
||
Line 3,396: | Line 3,396: | ||
In terms of fold/reduce: |
In terms of fold/reduce: |
||
{{Works with|Python|3.7}} |
{{Works with|Python|3.7}} |
||
< |
<syntaxhighlight lang="python">'''Gamma function''' |
||
from functools import reduce |
from functools import reduce |
||
Line 3,477: | Line 3,477: | ||
# MAIN --- |
# MAIN --- |
||
if __name__ == '__main__': |
if __name__ == '__main__': |
||
main()</ |
main()</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre> i -> gamma(i/3): |
<pre> i -> gamma(i/3): |
||
Line 3,495: | Line 3,495: | ||
Lanczos' approximation is loosely converted from the Octave code. |
Lanczos' approximation is loosely converted from the Octave code. |
||
{{trans|Octave}} |
{{trans|Octave}} |
||
< |
<syntaxhighlight lang="r">stirling <- function(z) sqrt(2*pi/z) * (exp(-1)*z)^z |
||
nemes <- function(z) sqrt(2*pi/z) * (exp(-1)*(z + (12*z - (10*z)^-1)^-1))^z |
nemes <- function(z) sqrt(2*pi/z) * (exp(-1)*(z + (12*z - (10*z)^-1)^-1))^z |
||
Line 3,544: | Line 3,544: | ||
all.equal(as.complex(gamma(z)), lanczos(z)) # TRUE |
all.equal(as.complex(gamma(z)), lanczos(z)) # TRUE |
||
all.equal(gamma(z), spouge(z)) # TRUE |
all.equal(gamma(z), spouge(z)) # TRUE |
||
data.frame(z=z, stirling=stirling(z), nemes=nemes(z), lanczos=lanczos(z), spouge=spouge(z), builtin=gamma(z))</ |
data.frame(z=z, stirling=stirling(z), nemes=nemes(z), lanczos=lanczos(z), spouge=spouge(z), builtin=gamma(z))</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
z stirling nemes lanczos spouge builtin |
z stirling nemes lanczos spouge builtin |
||
Line 3,559: | Line 3,559: | ||
=={{header|Racket}}== |
=={{header|Racket}}== |
||
< |
<syntaxhighlight lang="racket">#lang racket |
||
(define (gamma number) |
(define (gamma number) |
||
(if (> 1/2 number) |
(if (> 1/2 number) |
||
Line 3,587: | Line 3,587: | ||
; 1.1642297137253037 |
; 1.1642297137253037 |
||
; 1.068628702119319 |
; 1.068628702119319 |
||
; 1.0)</ |
; 1.0)</syntaxhighlight> |
||
=={{header|Raku}}== |
=={{header|Raku}}== |
||
(formerly Perl 6) |
(formerly Perl 6) |
||
<lang |
<syntaxhighlight lang="raku" line>sub Γ(\z) { |
||
constant g = 9; |
constant g = 9; |
||
z < .5 ?? pi/ sin(pi * z) / Γ(1 - z) !! |
z < .5 ?? pi/ sin(pi * z) / Γ(1 - z) !! |
||
Line 3,612: | Line 3,612: | ||
} |
} |
||
say Γ($_) for 1/3, 2/3 ... 10/3;</ |
say Γ($_) for 1/3, 2/3 ... 10/3;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>2.67893853470775 |
<pre>2.67893853470775 |
||
Line 3,632: | Line 3,632: | ||
Note: The Taylor series isn't much good above values of <big>'''6½'''</big>. |
Note: The Taylor series isn't much good above values of <big>'''6½'''</big>. |
||
< |
<syntaxhighlight lang="rexx">/*REXX program calculates GAMMA using the Taylor series coefficients; ≈80 decimal digits*/ |
||
/*The GAMMA function symbol is the Greek capital letter: Γ */ |
/*The GAMMA function symbol is the Greek capital letter: Γ */ |
||
numeric digits 90 /*be able to handle extended precision.*/ |
numeric digits 90 /*be able to handle extended precision.*/ |
||
Line 3,699: | Line 3,699: | ||
sum=sum*xm + #.k |
sum=sum*xm + #.k |
||
end /*k*/ |
end /*k*/ |
||
return 1/sum</ |
return 1/sum</syntaxhighlight> |
||
{{out|output|text= when using the input of: <tt> 0.5 </tt>}} |
{{out|output|text= when using the input of: <tt> 0.5 </tt>}} |
||
<pre> |
<pre> |
||
Line 3,723: | Line 3,723: | ||
Many of the "normal" high-level mathematical functions aren't available in REXX, so some of them (RYO) are included here. |
Many of the "normal" high-level mathematical functions aren't available in REXX, so some of them (RYO) are included here. |
||
< |
<syntaxhighlight lang="rexx">/*REXX program calculates the gamma function using Spouge's approximation with 87 digits*/ |
||
e=2.71828182845904523536028747135266249775724709369995957496696762772407663035354759457138 |
e=2.71828182845904523536028747135266249775724709369995957496696762772407663035354759457138 |
||
numeric digits length(e) - length(.) /*use the number of decimal digits in E*/ |
numeric digits length(e) - length(.) /*use the number of decimal digits in E*/ |
||
Line 3,791: | Line 3,791: | ||
do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/ |
do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/ |
||
do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/ |
do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/ |
||
numeric digits d; return g/1</ |
numeric digits d; return g/1</syntaxhighlight> |
||
{{out|output|text= when using the default input:}} |
{{out|output|text= when using the default input:}} |
||
<pre> |
<pre> |
||
Line 3,807: | Line 3,807: | ||
=={{header|Ring}}== |
=={{header|Ring}}== |
||
< |
<syntaxhighlight lang="ring"> |
||
decimals(3) |
decimals(3) |
||
gamma = 0.577 |
gamma = 0.577 |
||
Line 3,826: | Line 3,826: | ||
but fabs(z) <= 0.5 return 1 / recigamma(z) |
but fabs(z) <= 0.5 return 1 / recigamma(z) |
||
else return (z - 1) * gammafunc(z-1) ok |
else return (z - 1) * gammafunc(z-1) ok |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|RLaB}}== |
=={{header|RLaB}}== |
||
Line 3,841: | Line 3,841: | ||
====Taylor series==== |
====Taylor series==== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="ruby">$a = [ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
||
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
||
-0.00962_19715_27876_97356, 0.00721_89432_46663_09954, -0.00116_51675_91859_06511, |
-0.00962_19715_27876_97356, 0.00721_89432_46663_09954, -0.00116_51675_91859_06511, |
||
Line 3,857: | Line 3,857: | ||
end |
end |
||
(1..10).each {|i| puts format("%.14e", gamma(i/3.0))}</ |
(1..10).each {|i| puts format("%.14e", gamma(i/3.0))}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>2.67893853470775e+00 |
<pre>2.67893853470775e+00 |
||
Line 3,870: | Line 3,870: | ||
2.77815847933857e+00</pre> |
2.77815847933857e+00</pre> |
||
====Built in==== |
====Built in==== |
||
< |
<syntaxhighlight lang="ruby">(1..10).each{|i| puts Math.gamma(i/3.0)}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>2.678938534707748 |
<pre>2.678938534707748 |
||
Line 3,885: | Line 3,885: | ||
=={{header|Scala}}== |
=={{header|Scala}}== |
||
< |
<syntaxhighlight lang="scala">import java.util.Locale._ |
||
object Gamma { |
object Gamma { |
||
Line 3,910: | Line 3,910: | ||
println("%.1f -> %.16f %.16f".formatLocal(ENGLISH, x, stGamma(x), laGamma(x))) |
println("%.1f -> %.16f %.16f".formatLocal(ENGLISH, x, stGamma(x), laGamma(x))) |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>Gamma Stirling Lanczos |
<pre>Gamma Stirling Lanczos |
||
Line 3,939: | Line 3,939: | ||
{{trans|Ruby}} for Taylor |
{{trans|Ruby}} for Taylor |
||
< |
<syntaxhighlight lang="scheme"> |
||
(import (scheme base) |
(import (scheme base) |
||
(scheme inexact) |
(scheme inexact) |
||
Line 4,002: | Line 4,002: | ||
(number->string (gamma-taylor i)) |
(number->string (gamma-taylor i)) |
||
"\n"))) |
"\n"))) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
Line 4,089: | Line 4,089: | ||
=={{header|Scilab}}== |
=={{header|Scilab}}== |
||
<lang>function x=gammal(z) // Lanczos' |
<syntaxhighlight lang="text">function x=gammal(z) // Lanczos' |
||
lz=[ 1.000000000190015 .. |
lz=[ 1.000000000190015 .. |
||
76.18009172947146 -86.50532032941677 24.01409824083091 .. |
76.18009172947146 -86.50532032941677 24.01409824083091 .. |
||
Line 4,110: | Line 4,110: | ||
x=i/10 |
x=i/10 |
||
printf("%4.1f %9f %9f\n",x,gamma(x),gammal(x)) |
printf("%4.1f %9f %9f\n",x,gamma(x),gammal(x)) |
||
end</ |
end</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre style="height:20ex"> |
<pre style="height:20ex"> |
||
Line 4,148: | Line 4,148: | ||
=={{header|Seed7}}== |
=={{header|Seed7}}== |
||
{{trans|Ada}} |
{{trans|Ada}} |
||
< |
<syntaxhighlight lang="seed7">$ include "seed7_05.s7i"; |
||
include "float.s7i"; |
include "float.s7i"; |
||
Line 4,189: | Line 4,189: | ||
writeln((gamma(flt(I) / 3.0)) digits 15); |
writeln((gamma(flt(I) / 3.0)) digits 15); |
||
end for; |
end for; |
||
end func;</ |
end func;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>2.678937911987305 |
<pre>2.678937911987305 |
||
Line 4,204: | Line 4,204: | ||
=={{header|Sidef}}== |
=={{header|Sidef}}== |
||
{{trans|Ruby}} |
{{trans|Ruby}} |
||
< |
<syntaxhighlight lang="ruby">var a = [ 1.00000_00000_00000_00000, 0.57721_56649_01532_86061, -0.65587_80715_20253_88108, |
||
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
-0.04200_26350_34095_23553, 0.16653_86113_82291_48950, -0.04219_77345_55544_33675, |
||
-0.00962_19715_27876_97356, 0.00721_89432_46663_09954, -0.00116_51675_91859_06511, |
-0.00962_19715_27876_97356, 0.00721_89432_46663_09954, -0.00116_51675_91859_06511, |
||
Line 4,222: | Line 4,222: | ||
for i in 1..10 { |
for i in 1..10 { |
||
say ("%.14e" % gamma(i/3)) |
say ("%.14e" % gamma(i/3)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 4,238: | Line 4,238: | ||
Lanczos approximation: |
Lanczos approximation: |
||
< |
<syntaxhighlight lang="ruby">func gamma(z) { |
||
var epsilon = 0.0000001 |
var epsilon = 0.0000001 |
||
func withinepsilon(x) { |
func withinepsilon(x) { |
||
Line 4,274: | Line 4,274: | ||
for i in 1..10 { |
for i in 1..10 { |
||
say ("%.14e" % gamma(i/3)) |
say ("%.14e" % gamma(i/3)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 4,291: | Line 4,291: | ||
A simpler implementation: |
A simpler implementation: |
||
< |
<syntaxhighlight lang="ruby">define ℯ = Num.e |
||
define τ = Num.tau |
define τ = Num.tau |
||
Line 4,301: | Line 4,301: | ||
for i in (1..10) { |
for i in (1..10) { |
||
say ("%.14e" % Γ(i/3)) |
say ("%.14e" % Γ(i/3)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 4,321: | Line 4,321: | ||
The results are compared to Mata's '''[https://www.stata.com/help.cgi?mf_gamma gamma]''' function for each real between 1/100 and 100, by steps of 1/100. |
The results are compared to Mata's '''[https://www.stata.com/help.cgi?mf_gamma gamma]''' function for each real between 1/100 and 100, by steps of 1/100. |
||
< |
<syntaxhighlight lang="stata">mata |
||
_gamma_coef = 1.0, |
_gamma_coef = 1.0, |
||
5.772156649015328606065121e-1, |
5.772156649015328606065121e-1, |
||
Line 4,381: | Line 4,381: | ||
v=map(&gamma_(),x) |
v=map(&gamma_(),x) |
||
max(abs((v-u):/u)) |
max(abs((v-u):/u)) |
||
end</ |
end</syntaxhighlight> |
||
'''Output''' |
'''Output''' |
||
Line 4,389: | Line 4,389: | ||
Here follows the Python program to compute coefficients. |
Here follows the Python program to compute coefficients. |
||
< |
<syntaxhighlight lang="python">from mpmath import mp |
||
mp.dps = 50 |
mp.dps = 50 |
||
Line 4,419: | Line 4,419: | ||
for x in gc: |
for x in gc: |
||
print(mp.nstr(x, 25))</ |
print(mp.nstr(x, 25))</syntaxhighlight> |
||
=={{header|Tcl}}== |
=={{header|Tcl}}== |
||
Line 4,425: | Line 4,425: | ||
{{tcllib|math}} |
{{tcllib|math}} |
||
{{tcllib|math::calculus}} |
{{tcllib|math::calculus}} |
||
< |
<syntaxhighlight lang="tcl">package require math |
||
package require math::calculus |
package require math::calculus |
||
Line 4,459: | Line 4,459: | ||
set f2 [expr $f2 * $x] |
set f2 [expr $f2 * $x] |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> 1.0 1.000000 1.000000 1.000000 |
<pre> 1.0 1.000000 1.000000 1.000000 |
||
Line 4,486: | Line 4,486: | ||
As far as Gamma(n)=(n-1)! , we have everything needed. |
As far as Gamma(n)=(n-1)! , we have everything needed. |
||
===Stirling's approximation=== |
===Stirling's approximation=== |
||
< |
<syntaxhighlight lang="ti83b">for(I,1,10) |
||
I/2→X |
I/2→X |
||
X^(X-1/2)e^(-X)√(2π)→Y |
X^(X-1/2)e^(-X)√(2π)→Y |
||
Disp X,(X-1)!,Y |
Disp X,(X-1)!,Y |
||
Pause |
Pause |
||
End</ |
End</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
The output display for x=0.5 to 5 by 0.5 : x, (x-1)!, Y(x) . |
The output display for x=0.5 to 5 by 0.5 : x, (x-1)!, Y(x) . |
||
Line 4,529: | Line 4,529: | ||
===Lanczos' approximation=== |
===Lanczos' approximation=== |
||
< |
<syntaxhighlight lang="ti83b">for(I,1,10) |
||
I/2→X |
I/2→X |
||
e^(ln((1.0 |
e^(ln((1.0 |
||
Line 4,542: | Line 4,542: | ||
Disp X,(X-1)!,Y |
Disp X,(X-1)!,Y |
||
Pause |
Pause |
||
End</ |
End</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
The output display for x=0.5 to 5 by 0.5 : x, (x-1)!, Y(x) . |
The output display for x=0.5 to 5 by 0.5 : x, (x-1)!, Y(x) . |
||
Line 4,581: | Line 4,581: | ||
=={{header|Visual FoxPro}}== |
=={{header|Visual FoxPro}}== |
||
Translation of BBC Basic but with OOP extensions. Also some ideas from Numerical Methods (Press ''et al''). |
Translation of BBC Basic but with OOP extensions. Also some ideas from Numerical Methods (Press ''et al''). |
||
< |
<syntaxhighlight lang="vfp"> |
||
LOCAL i As Integer, x As Double, o As lanczos |
LOCAL i As Integer, x As Double, o As lanczos |
||
CLOSE DATABASES ALL |
CLOSE DATABASES ALL |
||
Line 4,642: | Line 4,642: | ||
ENDDEFINE |
ENDDEFINE |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 4,670: | Line 4,670: | ||
=={{header|Vlang}}== |
=={{header|Vlang}}== |
||
{{trans|go}} |
{{trans|go}} |
||
< |
<syntaxhighlight lang="vlang">import math |
||
fn main() { |
fn main() { |
||
println(" x math.Gamma Lanczos7") |
println(" x math.Gamma Lanczos7") |
||
Line 4,690: | Line 4,690: | ||
1.5056327351493116e-7/(z+7) |
1.5056327351493116e-7/(z+7) |
||
return math.sqrt2 * math.sqrt_pi * math.pow(t, z-.5) * math.exp(-t) * x |
return math.sqrt2 * math.sqrt_pi * math.pow(t, z-.5) * math.exp(-t) * x |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> x math.Gamma Lanczos7 |
<pre> x math.Gamma Lanczos7 |
||
Line 4,709: | Line 4,709: | ||
{{libheader|Wren-math}} |
{{libheader|Wren-math}} |
||
The ''gamma'' method in the Math class is based on the Lanczos approximation. |
The ''gamma'' method in the Math class is based on the Lanczos approximation. |
||
< |
<syntaxhighlight lang="ecmascript">import "/fmt" for Fmt |
||
import "/math" for Math |
import "/math" for Math |
||
Line 4,720: | Line 4,720: | ||
System.write("%(Fmt.f(16, stirling.call(d), 14))\t") |
System.write("%(Fmt.f(16, stirling.call(d), 14))\t") |
||
System.print("%(Fmt.f(16, Math.gamma(d), 14))") |
System.print("%(Fmt.f(16, Math.gamma(d), 14))") |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 4,750: | Line 4,750: | ||
=={{header|Yabasic}}== |
=={{header|Yabasic}}== |
||
{{trans|Phix}} |
{{trans|Phix}} |
||
< |
<syntaxhighlight lang="yabasic">dim c(12) |
||
sub gamma(z) |
sub gamma(z) |
||
Line 4,779: | Line 4,779: | ||
for i = 0.1 to 2.1 step .1 |
for i = 0.1 to 2.1 step .1 |
||
print i, " = "; : si(gamma(i)) |
print i, " = "; : si(gamma(i)) |
||
next</ |
next</syntaxhighlight> |
||
=={{header|zkl}}== |
=={{header|zkl}}== |
||
{{trans|D}} but without a built in gamma function. |
{{trans|D}} but without a built in gamma function. |
||
< |
<syntaxhighlight lang="zkl">fcn taylorGamma(x){ |
||
var table = T( |
var table = T( |
||
0x1p+0, 0x1.2788cfc6fb618f4cp-1, |
0x1p+0, 0x1.2788cfc6fb618f4cp-1, |
||
Line 4,805: | Line 4,805: | ||
return(1.0 / sm); |
return(1.0 / sm); |
||
}</ |
}</syntaxhighlight> |
||
< |
<syntaxhighlight lang="zkl">fcn lanczosGamma(z) { z = z.toFloat(); |
||
// Coefficients used by the GNU Scientific Library. |
// Coefficients used by the GNU Scientific Library. |
||
// http://en.wikipedia.org/wiki/Lanczos_approximation |
// http://en.wikipedia.org/wiki/Lanczos_approximation |
||
Line 4,833: | Line 4,833: | ||
return((2.0 * PI).sqrt() * t.pow(z + 0.5) * exp(-t) * x); |
return((2.0 * PI).sqrt() * t.pow(z + 0.5) * exp(-t) * x); |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |