I have an ugly numerical function $f(x)$ that diverges at the boundaries $x_a$ and $x_b$, but I want to calculate the integral, $$\int_{x_a}^{x_b} f(x) dx$$ For physical reasons I know that the integral needs to be convergent, but I don't know how to calculate it.
I have a series of values $g_1, g_2, ..., g_n$ corresponding to points $x_1, x_2, ..., x_n$. This is a smooth function which I can spline interpolate to be $g(x)$. This function is zero at the points $g(x_a) = g(x_b) = 0$ --- in particular it changes sign at these points, with $g(x_a < x < x_b) > 0$.
The function I need to integrate is $f(x) = 1/\sqrt{g(x)}$, which then diverges at $x_a$ and $x_b$.
edit:
My initial approach was to divide the region $[x_a,x_b]$ into three regions, $R_1:[x_a,x_a + \delta]$, $R_2:[x_a + \delta, x_b - \delta]$, $R_3:[x_b - \delta, x_b]$ and break each of these regions into $n$ intervals. Region $R_2$ is very well behaved, and I can easily calculate the integral using a Riemann sum with something like $n = 100$ segments.
To evaluate the accuracy of my integrals over the boundary regions, I compared left and right Riemann sums. As I decreased $\delta$, however, the difference between my Riemann sums diverged.