This .zip file contains 3 mfiles for computing the nodes and weights for Legendre Laguerre and Hermite - Gauss Quadrature of any order n.
Contrary to most of the files in use, the included files are all based on a symmetrical companion matrix, the eigenvalues of which are always real and correspond to the roots of the respective polynomials.
Hence, there is a 100% certainty of avoiding complex roots at high orders. It is known that the latter issue often causes huge numerical troubles.
Can I ask a question? Can I use this for multivariate normal case?
Very clever and precise.
Is it me or these routines are not working correctly in Matlab 2013b?
PS. I have been using them since 2010.
Simple and fast.
I added this on the end to refine the eigenvalues via newtons method
success = false;
for its=1:50 %maxit=50, usually will take 2
p1 = pi^(-1/4);
for j=1:n %make hermite we need
pp = sqrt(2*n).*p2;
if all(abs(z-z1)< 20*eps)
success = true;
warning('failed to converge to desired accuracy')
w(x>=0) = 2./pp.^2;
w(x<=0) = flipud(w(x>=0));
x(x>=0) = z;
x(x<=0) = -flipud(z);
Oh and I forgot to add this required you insist x = 0 is the middle value for odd via
Very nice submission, simple and fast
One addition I did is this
% Linear map from[-1,1] to [a,b]
a and b are user provided, and I think now that you can shift from [1 1] to [a b]
Very useful program.
How can I calculate the GaussLegendre in the [0 1] space?
Some minor changes to the code's comments.
No updates were made. Only the description was changed: it should read "the eigenvalues of the companion matrix" instead of "the roots of the companion matrix". Obviously, a matrix has no roots.