Take a Fredholm integral equation $$ u(x) - \lambda \int_{-1}^{1} K(x,y)u(y) \, \mathrm{d} y = f(x) $$ and discretize it via (say) Gaussian quadrature with nodes $\{x_j\}$ and weights $\{w_j\}$ to get $$ u(x) - \lambda \sum_{j} w_j K(x,x_j)u(x_j) = f(x) $$ It feels natural to use the same set of quadrature nodes over $f$ to get a linear system, so $$ u(x_k) - \lambda \sum_{j} w_j K(x_k,x_j)u(x_j) = f(x_k) $$ But if $K$ is singular, we have a big problem on out hands.
What are the normal methods of dealing with this?
Obviously we could use Gaussian quadrature nodes over $y$ and Chebyshev nodes over $x$, but that feels a little awkward; are there better ways?