Created function returning nodes and weights for 1D Gaussian quadrature.
I think the GL_nodes_and_weights works now
Added handling of n==1-case.
1D quad now takes nodes and weights as arguments.
Solver now gets nodes and weights from the get-go and uses these throughout the routine.