A Modern-Fortran module that calculates the Legendre-Gauss-Lobatto Quadrature as well as a useful set of functions for the Jacobi polynomials
Algorithms for computing these nodes and weights efficiently can be found in here.
The work from which this fortran implementation is inspired can be found in here.