Hermite polynomials in MATLAB

Have spent the afternoon ascertaining that

1. MATLAB has no built in function to generate Hermite polynomials, and
2. Various programmes on the InterWeb that purport to generate them do NOTHING of the sort, yet
3. This does not stop Randoms posting them, so
4. I have been forced to write my own.
   

If anyone finds it useful, please use the code (at your own risk etc). To normalise on the interval [+infinity..-infinity] with respect to the weighing function exp(-alpha x^2) just multiply by\sqrt{\frac{alpha}{\sqrt{\pi}2^n n!}}

The parameter k determines the largest hermite polynomial calculated (k-1)
The coefficients in x^0 x^1 x^3 are along the columns of the output matrix, and the
nth row is the (n-1)th polynomial. Try it out!

Code (Tested MATLAB 6.1)

% Function herm
% Generates the coefficients for a hemite polynomial of order 'k' from the recursion relation
% H_(n+1) = 2xH_n - 2nH_(n-1)

% Function output is a matrix
function A = herm(k)
% Generate a k+1 by k+1 matrix
A=zeros(k+1);
if k == 0; %only generate H_0
A(1,1)=1 ; % Defines H_0=1

elseif k == 1; % only generate H_1
A(1,1)=1 ; % Defines H_0=1
A(2,1)=0; % Defines H_1=2x
A(2,2)=2; % ditto


else

A(1,1)=1 ; % Defines H_0=1
A(2,1)=0; % Defines H_1=2x
A(2,2)=2; % ditto

for i=2:k
A(i+1,1)=-2*(i-1)*A(i-1,1) % -2nH_(n-1)

for j=2:k+1
A(i+1,j)=2*A(i,j-1) -2*(i-1)*A(i-1,j) % 2xH_n - 2nH_(n-1)
end
end

end