Computing Gauss Legendre quadrature for large $N$

There are asymptotic methods that essentially give you $N$ nodes and weights in $O(N)$ time if the precision is assumed to be fixed (e.g. at double precision).

See Nicholas Hale and Alex Townsend, "Fast and Accurate Computation of Gauss-Legendre and Gauss-Jacobi Quadrature Nodes and Weights", SIAM J. Sci. Comput., 35(2) (a PDF is available at http://eprints.maths.ox.ac.uk/1629/1/finalOR79.pdf).

They claim that their algorithm achieves double precision accuracy for $N \ge 100$.

For $N < 100$, you may as well precompute a big table with perfect accuracy using a computer algebra system or arbitrary precision library of your choice (or look up tables that others have published).

As to computing Legendre polynomials in a numerically stable way, use the three-term recurrence $(n+1) P_{n+1}(x) = (2n+1) x P_n(x) - n P_{n-1}(x)$ to evaluate $P(x)$ directly instead of computing the coefficients of the polynomial and using Horner's rule (similarly for $P'(x)$).

Update (2019): the issue of efficient computation for large N with arbitrary precision (and also with rigorous error bounds) is addressed in new work by myself and Marc Mezzarobba: SIAM Journal on Scientific Computing, 2018, Vol. 40, No. 6 : pp. C726-C747 https://doi.org/10.1137/18M1170133 (https://arxiv.org/abs/1802.03948).


To add to Fredrik Johansson's answer: A nice history of algorithms for computing Gauss quadrature rules can be found in this SIAM News article by Alex Townsend. Therein, it is stated that the "final chapter" was written by Ignace Bogaert in this SISC paper, which gives an algorithm that is even faster and more accurate than the algorithm of Hale & Townsend.

There is a free open-source implementation at http://sourceforge.net/projects/fastgausslegendrequadrature/.