--- /dev/null
+function P = legendre_p(n)
+ ## Return the nth legendre polynomial.
+ ##
+ ## INPUTS:
+ ##
+ ## * ``n`` - The index of the polynomial that we want.
+ ##
+ ## OUTPUTS:
+ ##
+ ## * ``P`` - A polynomial function of one argument.
+ ##
+ if (n < 0)
+ ## Can't do anything here. Return nothing.
+ P = NA;
+ elseif (n == 0)
+ ## One of our base cases.
+ P = @(x) 1
+ elseif (n == 1)
+ ## The second base case.
+ P = @(x) x
+ else
+ ## Compute recursively.
+ prev = legendre_p(n-1)
+ prev_prev = legendre_p(n-2)
+ P = @(x) (1/n)*( (2*n - 1)*prev(x) - (n-1)*prev_prev(x) )
+ end
+end