Gamma Function

The factorial $n!=1\cdot2\cdots n$ answers a counting question — how many ways to order $n$ objects — and so it is defined only at the whole numbers. Yet the same product appears everywhere a physicist does not want to be confined to integers: in the moments of a Gaussian, in the volume of a sphere in $d$ dimensions, in the normalization of a Bessel function, in the loop integrals of quantum field theory. What is $\left(\tfrac12\right)!$ ? What is the "volume of a $4.7$ -dimensional ball"? To ask these questions is to demand a smooth curve that passes through the factorial points and still obeys the one rule the factorial lives by — climb one step by multiplying by the next number. The gamma function $\Gamma(z)$ is that curve. Euler found it around 1730, and the Bohr–Mollerup theorem later proved it is the unique such interpolation once you add log-convexity: the only function with $\Gamma(1)=1$ , $\Gamma(z+1)=z\,\Gamma(z)$ , and a gently curving logarithm. Euler's integral $\Gamma(z)=\int_0^\infty t^{z-1}e^{-t}\,\mathrm{d}t$ gives a concrete formula in the right half-plane; the recurrence marches that definition leftward until $\Gamma$ becomes a meromorphic function on the whole complex plane, with simple poles at $0,-1,-2,\dots$ and no zeros at all. The reflection formula ties $\Gamma$ to the sine and hands us $\Gamma(1/2)=\sqrt{\pi}$ — the number behind every Gaussian normalization — while Stirling's asymptotic series governs its explosive growth and underlies entropy, ensembles, and the saddle-point method. The gamma function is, in short, the connective tissue of special-function theory.

Key takeaways

The gamma function $\Gamma(z)$ is the smooth interpolation of the factorial: $\Gamma(n)=(n-1)!$ for positive integers, but it is defined for every complex $z$ except $0,-1,-2,\dots$
Euler's integral $\Gamma(z)=\int_0^\infty t^{z-1}e^{-t}\,\mathrm{d}t$ converges and is holomorphic for $\mathrm{Re}(z)>0$ .
The functional equation $\Gamma(z+1)=z\,\Gamma(z)$ is the factorial recurrence; iterating it backward continues $\Gamma$ to the whole plane with simple poles of residue $(-1)^n/n!$ at $z=-n$ .
The reflection formula $\Gamma(z)\Gamma(1-z)=\pi/\sin(\pi z)$ yields $\Gamma(1/2)=\sqrt{\pi}\approx1.7725$ , the source of Gaussian normalization.
Stirling's series $\ln\Gamma(z)\sim(z-\tfrac12)\ln z-z+\tfrac12\ln 2\pi$ governs large arguments; its leading term already gives $10!$ to within $0.83\%$ .
The gamma function normalizes the beta function, Bessel and Legendre functions, the volume $V_n=\pi^{n/2}/\Gamma(n/2+1)$ of the $n$ -ball, and the dimensionally-regularized integrals of QFT.
Key idea: the gamma function is the unique smooth interpolation of the factorial obeying $\Gamma(z+1)=z\,\Gamma(z)$ ; that single recurrence reproduces the factorial, continues $\Gamma$ to a meromorphic function with simple poles at $0,-1,-2,\dots$ , and fixes its residues.

The core insight

The gamma function is the unique smooth interpolation of the factorial obeying $\Gamma(z+1)=z\,\Gamma(z)$ ; this single recurrence simultaneously reproduces the factorial, continues $\Gamma$ to a meromorphic function with simple poles at $0,-1,-2,\dots$ , and generates its residues.

Where this lives today

The gamma function is not a museum piece — it is live infrastructure. In quantum field theory, dimensional regularization tames ultraviolet divergences by continuing loop integrals to $d=4-\epsilon$ dimensions, where the divergences reappear as poles of $\Gamma$ near its non-positive integers; the modern computation of the electron's anomalous magnetic moment and of collider cross-sections runs through gamma functions of $\epsilon$ . In statistical mechanics, Stirling's series makes the thermodynamic limit tractable, turning $\ln N!$ into $N\ln N-N$ and delivering the Sackur–Tetrode entropy. Across probability and machine learning, the gamma, beta, Dirichlet, and Student- $t$ distributions are all built from $\Gamma$ . It even seeded string theory: Veneziano's 1968 scattering amplitude was literally a beta function $B(-\alpha(s),-\alpha(t))$ , and reading its gamma-pole structure as an infinite tower of resonances launched the dual-resonance model. Numerically, the Lanczos approximation computes $\Gamma$ to full double precision in a few terms, so every scientific library carries an implementation.

The physical picture

Picture the factorial as a staircase — isolated points $0!,1!,2!,\dots$ with nothing between them. The gamma function draws the smooth banister threading exactly through those points (shifted by one: $\Gamma(n)=(n-1)!$ ) while never breaking the rule that connects each step to the next, $\Gamma(z+1)=z\,\Gamma(z)$ . That single rule is all the gamma function remembers about the factorial, and it forces everything else — the dip between $1$ and $2$ , the runaway growth to the right, and the tower of poles marching off to the left. The first figure is that curve.

Figure

The gamma function Γ(x) on the real axis. It threads the factorial points Γ(1)=0!=1, Γ(2)=1!=1, Γ(3)=2!=2, Γ(4)=3!=6, dips to a single minimum Γ(1.4616)=0.8856 between 1 and 2, and blows up to simple poles at x = 0, −1, −2, −3 whose lobes alternate in sign

The poles are not a defect — they are forced. Read the recurrence backward, $\Gamma(z)=\Gamma(z+1)/z$ : each leftward step divides by a number that shrinks toward zero, so at $z=0,-1,-2,\dots$ you are dividing a finite value by zero and a simple pole is born. Their residues, $(-1)^n/n!$ , alternate in sign, which is exactly why the lobes between the poles flip up, down, up in the figure. To the right the curve does the opposite of pathology: it grows faster than any exponential, and how fast is the content of Stirling's formula. The second figure makes that asymptotics quantitative — it shows how quickly the simplest Stirling estimate closes in on the true factorial as the argument grows.

Figure

How good is Stirling? The relative error of the leading term √(2πn)(n/e)ⁿ against the exact n! falls as ≈ 1/(12n): 7.79% at n=1, 4.05% at n=2, 1.65% at n=5, and just 0.83% at n=10 — the blue measured error tracking the amber 1/(12n) prediction

The two figures tell the whole story of the function's shape: a gentle interpolation with a single minimum, framed on the left by an infinite ladder of poles and stretching to the right with a growth law so regular that one crude formula pins it to a fraction of a percent.

The core equations

The starting point is a definite integral that agrees with the factorial and makes sense for continuous arguments. Euler's integral of the second kind converges wherever the integrand is both integrable at the origin ( $t^{z-1}$ must not diverge too fast) and decaying at infinity (the $e^{-t}$ guarantees it), which pins the domain to the right half-plane.

\Gamma(z)=\int_0^\infty t^{z-1}e^{-t}\,\mathrm{d}t,\qquad \mathrm{Re}(z)>0 \qquad (1)

where: $z$ — complex argument, must have $\mathrm{Re}(z)>0$ ; $t$ — real integration variable on $(0,\infty)$ . Reading it: $\Gamma$ is the "continuous sum" of $t^{z-1}$ weighted by a decaying exponential; differentiating under the integral shows it is holomorphic throughout $\mathrm{Re}(z)>0$ . Everything structural flows from one integration by parts. Move the derivative from $e^{-t}$ onto $t^{z}$ and the boundary term vanishes, leaving the recurrence that encodes the factorial's "multiply by the next integer" rule.

\Gamma(z+1)=z\,\Gamma(z) \qquad (2)

where: $z$ — any argument where both sides are defined. Reading it: with $\Gamma(1)=1$ this forces $\Gamma(2)=1,\ \Gamma(3)=2,\ \Gamma(4)=6,\dots$ , i.e. $\Gamma(n)=(n-1)!$ — and, run backward, it is the engine of analytic continuation. The next identity is the deepest elementary fact about $\Gamma$ . It multiplies the value at $z$ by the value at the reflected point $1-z$ and returns something whose only singularities are the integer poles of $\Gamma$ — precisely $\pi/\sin(\pi z)$ .

\Gamma(z)\,\Gamma(1-z)=\frac{\pi}{\sin(\pi z)} \qquad (3)

where: $z$ — any non-integer complex number. Reading it: at the self-dual point $z=\tfrac12$ the two factors coincide, giving $\Gamma(\tfrac12)^2=\pi/\sin(\pi/2)=\pi$ , hence $\Gamma(\tfrac12)=\sqrt{\pi}$ ; more generally it links the "left" and "right" halves of the plane. Two product/limit forms make the analytic structure manifest. Gauss's limit rebuilds $\Gamma$ from finite factorials, and Weierstrass's canonical product exhibits $1/\Gamma$ as an entire function whose zeros are exactly the poles of $\Gamma$ .

\Gamma(z)=\lim_{n\to\infty}\frac{n!\,n^{z}}{z(z+1)\cdots(z+n)},\qquad \frac{1}{\Gamma(z)}=z\,e^{\gamma z}\prod_{k=1}^{\infty}\left(1+\frac{z}{k}\right)e^{-z/k} \qquad (4)

where: $\gamma\approx0.5772156649$ — the Euler–Mascheroni constant; $k,n$ — positive integers. Reading it: because $1/\Gamma$ is entire with simple zeros at $0,-1,-2,\dots$ , the gamma function itself has no zeros anywhere and only simple poles there. Half-integer arguments proliferate in physics (angular momentum, $d$ -dimensional integrals), and Legendre's duplication formula is the tool that tames them, relating $\Gamma(2z)$ to the pair $\Gamma(z),\Gamma(z+\tfrac12)$ .

\Gamma(2z)=\frac{2^{2z-1}}{\sqrt{\pi}}\,\Gamma(z)\,\Gamma\!\left(z+\tfrac12\right) \qquad (5)

where: $z$ — any argument avoiding the poles of the three factors. Reading it: set $z=\tfrac12$ and it collapses to $\Gamma(1)=1$ , a quick consistency check; it is the $N=2$ case of the Gauss multiplication theorem. The gamma function's closest relative is the beta function, the Euler integral of the first kind. The two are joined by an identity that turns a product of gammas into a single normalized integral — the workhorse of Feynman-parameter and probability calculations.

B(x,y)=\int_0^1 t^{x-1}(1-t)^{y-1}\,\mathrm{d}t=\frac{\Gamma(x)\,\Gamma(y)}{\Gamma(x+y)} \qquad (6)

where: $x,y$ — complex with positive real parts; $B$ — the beta function. Reading it: $B$ is symmetric, $B(x,y)=B(y,x)$ , and $B(1,1)=1=\Gamma(1)^2/\Gamma(2)$ — the identity is derived below. Finally, for large arguments the growth is captured by Stirling's asymptotic series — a divergent expansion whose first few terms nonetheless give extraordinary accuracy, and the single most-used approximation in statistical physics.

\ln\Gamma(z)\sim\left(z-\tfrac12\right)\ln z-z+\tfrac12\ln(2\pi)+\frac{1}{12z}-\frac{1}{360z^{3}}+\cdots \qquad (7)

where: $z$ — large argument with $|\arg z|<\pi$ . Reading it: exponentiating the leading terms gives $\Gamma(z)\sim\sqrt{2\pi}\,z^{z-1/2}e^{-z}$ , and for the factorial $\ln N!\approx N\ln N-N$ — the form that appears in entropy.

Derivation AAnalytic continuation and the tower of poles

The goal of this derivation is to turn the recurrence $\Gamma(z+1)=z\,\Gamma(z)$ into a machine for extending $\Gamma$ from its native strip $\mathrm{Re}(z)>0$ to the entire complex plane, and to read off — with no extra input — where the singularities sit and how strong each one is. We will iterate the functional equation backward, watch the domain grow one unit to the left at every step, and then compute the residue at a generic pole $z=-n$ by isolating the offending factor. The payoff is the complete analytic character of $\Gamma$ : meromorphic, simple poles at the non-positive integers, residues $(-1)^n/n!$ , and never a zero.

Step 1. Rearrange the recurrence as $\Gamma(z)=\Gamma(z+1)/z$ . The right-hand side is defined wherever $\Gamma(z+1)$ is — that is, for $\mathrm{Re}(z+1)>0$ , i.e. $\mathrm{Re}(z)>-1$ — except at $z=0$ where we divide by zero. Why: this single rearrangement already extends $\Gamma$ leftward by one unit and exposes the first pole, at the origin, as nothing more than a division by $z$ .

Step 2. Iterate $n+1$ times to reach $\Gamma(z)=\dfrac{\Gamma(z+n+1)}{z(z+1)(z+2)\cdots(z+n)}$ . The numerator $\Gamma(z+n+1)$ is holomorphic for $\mathrm{Re}(z)>-n-1$ , and the denominator is a polynomial vanishing simply at $z=0,-1,-2,\dots,-n$ . Why: because $n$ is arbitrary, this formula defines $\Gamma$ on the whole plane, and it displays every singularity as a simple zero of an explicit polynomial — hence every pole is simple.

Step 3. Compute the residue at $z=-n$ by multiplying through by $(z+n)$ and taking the limit. Only the factor $(z+n)$ in the denominator is cancelled; the surviving factors are evaluated at $z=-n$ , and the numerator tends to $\Gamma(-n+n+1)=\Gamma(1)=1$ . This gives $\mathrm{Res}_{z=-n}\Gamma=\dfrac{1}{(-n)(-n+1)\cdots(-1)}$ . Why: the residue is precisely the coefficient that controls the strength of the pole, and it drops out cleanly because all other factors stay finite there.

Step 4. Simplify the denominator: the product $(-n)(-n+1)\cdots(-1)$ has $n$ factors, each negative, with magnitudes $1,2,\dots,n$ , so it equals $(-1)^n\,n!$ . Therefore $\mathrm{Res}_{z=-n}\Gamma=\dfrac{(-1)^n}{n!}$ . Why: the alternating sign is exactly what makes successive lobes of $\Gamma(x)$ in the first figure flip from positive to negative, and it is checkable at $n=0$ , where $\mathrm{Res}_{z=0}=1$ matches $\Gamma(z)\approx1/z$ near the origin. *As a quick check: at $n=0$ the formula gives $\Gamma(z)\approx 1/z$ near the origin, and at $n=1$ it gives residue $-1$ at $z=-1$ , both of which are easy to confirm directly from the continued expression.*

Go deeper · Derivation B · 1 minThe value at one half and the Gaussian integral

Derivation BThe value at one half and the Gaussian integral

The single most consequential number in the theory is $\Gamma(1/2)$ , because it is the bridge between the gamma function and the Gaussian that pervades statistics and quantum mechanics. The goal here is to evaluate $\Gamma(1/2)$ directly from Euler's integral, by the substitution that turns the awkward $t^{-1/2}$ weight into a clean Gaussian, and then to cross-check the answer against the reflection formula. The result, $\Gamma(1/2)=\sqrt{\pi}\approx1.7725$ , is what forces the $\sqrt{\pi}$ into every normal-distribution normalization and every half-integer gamma value.

Step 1. Set $z=\tfrac12$ in Euler's integral: $\Gamma(\tfrac12)=\int_0^\infty t^{-1/2}e^{-t}\,\mathrm{d}t$ . The integrand behaves like $t^{-1/2}$ near $0$ , which is integrable since $\tfrac12<1$ , so the object is well defined. Why: half-integer arguments are exactly where the plain factorial fails, so this integral is the first genuinely new value the gamma function supplies.

Step 2. Substitute $t=u^2$ , so $\mathrm{d}t=2u\,\mathrm{d}u$ and $t^{-1/2}=u^{-1}$ ; as $t$ runs over $(0,\infty)$ so does $u$ . The integral becomes $\int_0^\infty u^{-1}e^{-u^2}\,(2u\,\mathrm{d}u)=2\int_0^\infty e^{-u^2}\,\mathrm{d}u$ . Why: the substitution is engineered so the $u^{-1}$ from the weight and the $u$ from $\mathrm{d}t$ cancel, leaving a pure Gaussian with no leftover powers.

Step 3. Evaluate with the standard Gaussian integral $\int_{-\infty}^{\infty}e^{-u^2}\,\mathrm{d}u=\sqrt{\pi}$ , so $\int_0^\infty e^{-u^2}\,\mathrm{d}u=\tfrac12\sqrt{\pi}$ and $\Gamma(\tfrac12)=2\cdot\tfrac12\sqrt{\pi}=\sqrt{\pi}$ . Why: this is the moment the transcendental constant $\pi$ enters the factorial's world — through the area under a bell curve.

Step 4. Confirm independently with the reflection formula at $z=\tfrac12$ : $\Gamma(\tfrac12)\Gamma(\tfrac12)=\pi/\sin(\pi/2)=\pi$ , so $\Gamma(\tfrac12)=\sqrt{\pi}$ (positive root, since the integrand is everywhere positive). Why: two independent derivations agreeing is the working physicist's proof, and it also fixes the sign, which the squaring alone would leave ambiguous. *Numerically $\Gamma(0.5)=1.7724538509=\sqrt{\pi}$ ; the recurrence then gives $\Gamma(3/2)=\tfrac12\sqrt{\pi}=0.8862$ and $\Gamma(5/2)=\tfrac34\sqrt{\pi}=1.3293$ .*

Go deeper · Derivation C · 1 minThe beta function factorizes into gammas

Derivation CThe beta function factorizes into gammas

Many physics integrals are not gamma integrals but beta integrals — anything of the form $\int_0^1 t^{a-1}(1-t)^{b-1}\,\mathrm{d}t$ , which shows up in Feynman parameterization, in the Dirichlet and Student distributions, and historically in the Veneziano amplitude. The goal of this derivation is to prove the beta–gamma relation $B(x,y)=\Gamma(x)\Gamma(y)/\Gamma(x+y)$ , which lets any such integral be read straight off a table of gamma values. The strategy is to multiply two gamma integrals, interpret the result over the first quadrant, and change to a radial-plus-fraction coordinate that separates the two effects.

Step 1. Write the product as a double integral: $\Gamma(x)\Gamma(y)=\int_0^\infty\!\int_0^\infty s^{x-1}t^{y-1}e^{-(s+t)}\,\mathrm{d}s\,\mathrm{d}t$ , over the first quadrant $s,t>0$ . Why: turning a product of one-dimensional integrals into a single two-dimensional one is what allows a change of variables to mix them — the same trick that evaluates the Gaussian integral by going to polar coordinates.

Step 2. Introduce $r=s+t$ and $u=s/(s+t)$ , i.e. $s=ru$ and $t=r(1-u)$ , with $r\in(0,\infty)$ and $u\in(0,1)$ . The Jacobian is $\partial(s,t)/\partial(r,u)=r$ , so $\mathrm{d}s\,\mathrm{d}t=r\,\mathrm{d}r\,\mathrm{d}u$ . Why: $r$ captures the total "size" that the exponential $e^{-r}$ decays with, while $u$ captures the proportion between $s$ and $t$ — the two roles the beta and gamma factors will inherit.

Step 3. Substitute: the integrand $s^{x-1}t^{y-1}=r^{x+y-2}u^{x-1}(1-u)^{y-1}$ , and multiplying by the Jacobian $r$ raises the power of $r$ to $x+y-1$ . The double integral separates as $\left(\int_0^\infty r^{x+y-1}e^{-r}\,\mathrm{d}r\right)\left(\int_0^1 u^{x-1}(1-u)^{y-1}\,\mathrm{d}u\right)$ . Why: the exponential depends only on $r$ and the polynomial weight factors into an $r$ -part and a $u$ -part, so the integral is a clean product — the whole point of choosing these coordinates.

Step 4. Identify the pieces: the $r$ -integral is $\Gamma(x+y)$ by definition, and the $u$ -integral is $B(x,y)$ . Hence $\Gamma(x)\Gamma(y)=\Gamma(x+y)\,B(x,y)$ , which rearranges to the boxed relation. Why: with $\Gamma$ tabulated, every beta integral is now a division — for instance $B(\tfrac12,\tfrac12)=\Gamma(\tfrac12)^2/\Gamma(1)=\pi$ , the area formula behind the arcsine law. *As consistency checks: the symmetry $B(x,y)=B(y,x)$ holds by inspection; $B(1,1)=\Gamma(1)^2/\Gamma(2)=1$ ; and $B(\tfrac12,\tfrac12)=\pi$ .*

Limiting cases

Special values and asymptotic regimes worth committing to memory:

Case	Result	Meaning
Integer argument	$\Gamma(n)=(n-1)!\quad(n=1,2,3,\dots)$	Reproduces the factorial with the index shifted down by one — the perennial off-by-one.
Value at one half	$\Gamma\!\left(\tfrac12\right)=\sqrt{\pi}\approx1.7725$	Source of Gaussian normalization; seeds all half-integers via $\Gamma(z+1)=z\Gamma(z)$ .
Near a pole	$\Gamma(z)\sim\dfrac{(-1)^n}{n!}\,\dfrac{1}{z+n}\quad(z\to -n)$	Simple pole at every non-positive integer, residue $(-1)^n/n!$ alternating in sign.
No zeros	$1/\Gamma(z)$ entire; $\Gamma(z)\neq0$ for all $z$	Weierstrass product puts the zeros on $1/\Gamma$ , so $\Gamma$ itself never vanishes.
Large argument (Stirling)	$\Gamma(z)\sim\sqrt{2\pi}\,z^{z-\frac12}e^{-z}\quad(\|z\|\to\infty,\ \|\arg z\|<\pi)$	Leading term of a divergent asymptotic series; the backbone of saddle-point methods.

Worked examples

1 · Half-integer values and the volume of a hypersphere (fully worked)

Physics constantly needs $\Gamma$ at half-integers — for the surface area of spheres in $d$ dimensions and for Gaussian moments. Compute $\Gamma(7/2)$ from scratch, then use $V_n=\pi^{n/2}/\Gamma(n/2+1)$ to find the volume of the unit ball in $n=3$ and $n=4$ dimensions.

Step 1. Build $\Gamma(7/2)$ up from $\Gamma(1/2)=\sqrt{\pi}$ using the recurrence $\Gamma(z+1)=z\,\Gamma(z)$ three times: $\Gamma(3/2)=\tfrac12\Gamma(1/2)=\tfrac12\sqrt{\pi}$ , then $\Gamma(5/2)=\tfrac32\Gamma(3/2)=\tfrac34\sqrt{\pi}$ , then $\Gamma(7/2)=\tfrac52\Gamma(5/2)=\tfrac{15}{8}\sqrt{\pi}$ . Numerically $\tfrac{15}{8}\sqrt{\pi}=1.875\times1.77245=\mathbf{3.3234}$ .

Step 2. For $n=3$ the formula needs $\Gamma(3/2+1)=\Gamma(5/2)=\tfrac34\sqrt{\pi}$ . Then $V_3=\dfrac{\pi^{3/2}}{\tfrac34\sqrt{\pi}}=\dfrac{4}{3}\cdot\dfrac{\pi^{3/2}}{\sqrt{\pi}}=\dfrac{4}{3}\pi$ . This is the familiar $\tfrac43\pi r^3$ at $r=1$ : $V_3=4\pi/3\approx4.1888$ — the gamma formula reproduces the schoolbook sphere.

Step 3. For $n=4$ the formula needs $\Gamma(4/2+1)=\Gamma(3)=2!=2$ . Then $V_4=\dfrac{\pi^{2}}{\Gamma(3)}=\dfrac{\pi^2}{2}=\mathbf{4.9348}$ . Notice the pattern: the $4$ -ball has larger volume than the $3$ -ball, but ball volumes eventually peak (near $n=5$ ) and then fall toward zero — a fact only visible because $\Gamma$ lets $n$ vary continuously. The half-integer machinery ( $\Gamma(7/2)=\tfrac{15}{8}\sqrt{\pi}$ ) and the integer case ( $\Gamma(3)=2$ ) are the same function evaluated at different points.

2 · Stirling's approximation for a physicist's factorial

In statistical mechanics one constantly needs $\ln N!$ for large $N$ . Estimate $10!$ with the leading Stirling term, quantify the error, and show the first correction nearly eliminates it.

Step 1. Leading Stirling for the factorial is $N!\approx\sqrt{2\pi N}\,(N/e)^N$ . With $N=10$ : $\sqrt{2\pi\cdot10}=\sqrt{62.832}=7.9266$ and $(10/e)^{10}=10^{10}/e^{10}=10^{10}/22026.47=4.5400\times10^{5}$ . Multiplying, $10!\approx7.9266\times4.5400\times10^{5}=\mathbf{3.5987\times10^{6}}$ .

Step 2. The exact value is $10!=3\,628\,800$ . The relative error is $(3\,628\,800-3\,598\,696)/3\,628\,800=30\,104/3\,628\,800=\mathbf{0.83\%}$ , matching the rule of thumb $1/(12N)=1/120=0.833\%$ read off the second figure.

Step 3. Apply the first correction factor $1+1/(12N)$ : $3\,598\,696\times(1+1/120)=3\,598\,696\times1.008333=\mathbf{3\,628\,685}$ , now within $115$ of the exact $3\,628\,800$ — a relative error of $3\times10^{-5}$ . One crude formula plus one correction term nails a seven-digit factorial to five significant figures, which is exactly why $\ln N!\approx N\ln N-N$ is trusted throughout thermodynamics.

3 · A reflection-formula check

Verify the reflection formula numerically at $z=1/4$ .

Step 1. The formula predicts $\Gamma(1/4)\,\Gamma(3/4)=\pi/\sin(\pi/4)=\pi/(\tfrac{\sqrt2}{2})=\pi\sqrt2=\mathbf{4.4429}$ .

Step 2. The tabulated values are $\Gamma(1/4)=3.62561$ and $\Gamma(3/4)=1.22542$ ; their product is $3.62561\times1.22542=\mathbf{4.4429}$ . The two agree to all quoted digits, confirming that the reflection formula ties together values a full half-unit apart. $\Gamma(1/4)\Gamma(3/4)=\pi\sqrt2$ , and neither factor alone has a closed form in elementary constants — the product does.

Physics of the portrait

The gamma function is Euler's, born from a question Daniel Bernoulli and Goldbach were tossing around in the 1720s: what is the factorial of a non-integer? Euler answered it twice — with an infinite product and with the integral we still use — and in doing so opened the whole subject of special functions.

Origins · 1 minLeonhard Euler and how this was discovered

Leonhard Euler

1707–1783 · the interpolation of the factorial and the birth of special-function theory

In a letter to Christian Goldbach dated 13 October 1729, the twenty-two-year-old Leonhard Euler solved a problem circulating in the Basel–St. Petersburg correspondence: how to interpolate $1,2,6,24,120,\dots$ to non-integer indices. His first answer was an infinite product; within a year he recast it as the integral $\int_0^1(-\ln s)^{z}\,\mathrm{d}s$ , which under $s=e^{-t}$ becomes the modern $\int_0^\infty t^{z}e^{-t}\,\mathrm{d}t$ . Euler was not chasing a definition for its own sake — he needed fractional factorials to sum series and evaluate integrals in mechanics and number theory. He computed $\Gamma(1/2)=\sqrt{\pi}$ , found the reflection formula, and connected the object to the beta function. The name and the symbol $\Gamma$ came later, from Adrien-Marie Legendre in 1811 (Gauss preferred a shifted $\Pi(z)=\Gamma(z+1)$ ). But the conceptual leap — that the discrete factorial is the restriction of a single smooth, complex-analytic function — is Euler's, and it set the template for every special function that followed: define by an integral, extend by a functional equation, understand through its singularities. Karl Weierstrass completed the analytic picture in the 1850s with the entire product for $1/\Gamma$ , and the Bohr–Mollerup theorem (1922) proved Euler's function was the only reasonable answer all along.

History

1729

— Euler interpolates the factorial to non-integer arguments in a letter to Goldbach

1730

— Stirling publishes the asymptotic growth law later exponentiated into Γ

1811

— Legendre introduces the symbol Γ and proves the duplication formula

1812

— Gauss studies the shifted Π-function for general complex argument

1922

— Bohr and Mollerup characterize Γ uniquely by log-convexity

Going further

Texts. Arfken, Weber & Harris, Mathematical Methods for Physicists, 7th ed., Ch. 13 (Gamma Function) · Hassani, Mathematical Physics, Ch. on the gamma and beta functions · Stone & Goldbart, Mathematics for Physics, special-functions chapter · Whittaker & Watson, A Course of Modern Analysis, Ch. XII · Abramowitz & Stegun, Handbook of Mathematical Functions, Ch. 6 · Byron & Fuller, Mathematics of Classical and Quantum Physics, §6. Key numerical anchors used throughout this article: the integration-by-parts proof establishes $\Gamma(z+1)=z\Gamma(z)$ ; the residues are $(-1)^n/n!$ ; the $t=u^2$ substitution gives $\Gamma(1/2)=\sqrt{\pi}=1.77245$ ; the Jacobian- $r$ change of variables yields $B(x,y)=\Gamma(x)\Gamma(y)/\Gamma(x+y)$ ; the half-integer chain gives $\Gamma(7/2)=\tfrac{15}{8}\sqrt{\pi}=3.3234$ ; the ball volumes are $V_3=4\pi/3=4.1888$ and $V_4=\pi^2/2=4.9348$ ; the Stirling leading-term error at $N=10$ is $0.83\%$ ; and the Euler–Mascheroni constant is $\gamma\approx0.5772156649$ .

Prerequisites. Improper Integrals and Convergence · Analytic Continuation of the Gamma Function · Integration by Parts. Leads to. Beta Function · Analytic Continuation and Poles of the Gamma Function · Bessel Functions · Legendre Polynomials · Recurrence Relations for Special Functions.

Ask the tutor: Ask me to continue any half-integer chain, expand $\Gamma$ near a pole for a dimensional-regularization calculation, or push Stirling's series to the accuracy your problem needs.

Why is Γ(n) equal to (n−1)! and not n!?

It is a bookkeeping consequence of where Euler anchored the integral. The Euler integral gives $\Gamma(1)=\int_0^\infty e^{-t}\,\mathrm{d}t=1$ , and the recurrence $\Gamma(z+1)=z\Gamma(z)$ then climbs as $\Gamma(2)=1\cdot\Gamma(1)=1$ , $\Gamma(3)=2\cdot\Gamma(2)=2$ , $\Gamma(4)=3\cdot\Gamma(3)=6$ , so $\Gamma(n)=(n-1)!$ . If you want a function $\Pi$ with $\Pi(n)=n!$ exactly, that is Gauss's shifted convention $\Pi(z)=\Gamma(z+1)$ . Physicists mostly use $\Gamma$ and simply remember the one-unit shift; the safest habit is to write $\Gamma(z+1)=z!$ only as informal shorthand and to carry the shift explicitly in any real calculation.

If the gamma function blows up at 0, −1, −2, …, is it useless there?

Not at all — the poles are informative, not fatal. Near $z=-n$ the function behaves like $\Gamma(z)\approx\dfrac{(-1)^n}{n!}\,\dfrac{1}{z+n}$ , a clean simple pole with a known residue, so any expression containing $\Gamma$ can be Laurent-expanded and its singular part handled analytically. This is precisely how dimensional regularization works: a divergent loop integral is written as a $\Gamma(\epsilon)$ -type object, its $1/\epsilon$ pole is isolated using exactly these residues, and the pole is absorbed into a renormalization. Far from a defect, the pole structure is what makes $\Gamma$ the natural language for divergences.

Does the gamma function ever equal zero?

Never — a useful property and a common point of confusion. The Weierstrass product writes $1/\Gamma(z)=z\,e^{\gamma z}\prod_{k\ge1}(1+z/k)e^{-z/k}$ , an entire function (no poles) whose zeros are exactly $z=0,-1,-2,\dots$ Since $\Gamma=1/(1/\Gamma)$ , those zeros become the poles of $\Gamma$ , and $\Gamma$ has no zeros of its own anywhere. Practically, you may freely divide by $\Gamma$ — the $1/\Gamma(z)$ ratios in Bessel-function series are perfectly regular, and at the negative integers $1/\Gamma=0$ simply kills those terms.

Why does √π appear in a function built from factorials?

Because the gamma function at a half-integer is secretly a Gaussian integral. Setting $z=1/2$ in Euler's integral and substituting $t=u^2$ turns $\Gamma(1/2)$ into $2\int_0^\infty e^{-u^2}\,\mathrm{d}u=\sqrt{\pi}$ — the area under a bell curve. Once $\Gamma(1/2)=\sqrt{\pi}$ is in hand, the recurrence propagates the $\sqrt{\pi}$ to every half-integer: $\Gamma(3/2)=\tfrac12\sqrt{\pi}$ , $\Gamma(5/2)=\tfrac34\sqrt{\pi}$ , and so on. This is why $\sqrt{\pi}$ (or $\sqrt{2\pi}$ ) shows up in the normalization of the normal distribution, in Gaussian moment integrals $\int_0^\infty x^n e^{-ax^2}\mathrm{d}x=\tfrac12 a^{-(n+1)/2}\Gamma(\tfrac{n+1}{2})$ , and in the volumes of even-dimensional balls.

Stirling's series is divergent — how can it be so accurate?

Stirling's expansion is an asymptotic series, not a convergent one: for fixed $z$ its terms eventually grow without bound, so adding all of them is meaningless. What makes it useful is that for large $z$ the terms first shrink dramatically before they start growing, and truncating just before the smallest term gives a superb approximation. For $N=10$ the leading term alone is already good to $0.83\%$ , and one correction term reaches five significant figures. The optimal truncation error is roughly the size of the first omitted term. Physicists exploit exactly this in the saddle-point (steepest-descent) method, where Stirling-type asymptotics control the leading exponential and the Gaussian fluctuations around it.

What is the relationship between the gamma and beta functions?

The beta function $B(x,y)=\int_0^1 t^{x-1}(1-t)^{y-1}\,\mathrm{d}t$ is the Euler integral of the first kind, and it factorizes completely into gammas: $B(x,y)=\Gamma(x)\Gamma(y)/\Gamma(x+y)$ . This makes beta a "normalized" partner to gamma — symmetric, finite on the unit interval, and appearing wherever a probability weight $t^{a-1}(1-t)^{b-1}$ lives (the beta and Dirichlet distributions), in Feynman-parameter integrals, and in the Veneziano amplitude. Whenever you meet an integral over $[0,1]$ with two power-law endpoints, it is a beta integral and the gamma table finishes the job in one line.