In https://www.johndcook.com/blog/2012/02/21/care-and-treatment-of-singularities, the author explains the subtraction method to get rid of singularities when performing numerical integration.
The example he gives is the integral $\int_{0}^{1}\frac{1}{\sqrt{\sin x}}dx$.
To accurately compute this integral numerically, one first finds the approximation of integrand near $x=0$, which is $\frac{1}{\sqrt{x}}$ for this case.
Then, the original integral can be written like $\int_{0}^{1}(\frac{1}{\sqrt{\sin x}}-\frac{1}{\sqrt{x}})dx+\int_{0}^{1}\frac{1}{\sqrt{x}}dx$.
The latter term can be computed analytically by hand, so the author mentions that the work only left is computing the former one with numerical integration method such as trapezoidal, Simpson, etc.
I understand that the singularity is eliminated from the former integrand but I could not understand how one can actually compute it with numerical method.
For example, let's say I want to compute it with Simpson's method:
$\int_{0}^{2h}f(x)dx=\frac{h}{3}[f_{0}+4f_{1}+f_{2}]+\mathcal{O}(h^{5})$.
Then, I must evaluate $f(0)$ with the computer to fully compute the former integral. However, how can one compute $\frac{1}{\sqrt{sin(0)}}-\frac{1}{\sqrt{0}}$ with the computer? Shouldn't it give DivisionByZero error?
I want to understand how can this problem solved or am I misunderstanding something?
Best Answer
Obviously you may not attempt numerical evaluation of the function at the origin. But obviously also the value is $0$ so that you needn't even compute the term !