Numerical integration of cauchy principal value

I have an integral of the form:

I = \int\limits_{-1}^{1} f(x)*w(x) dx,

The function f(x) has a singularity of the form 1/(x-a) inside the interval [-1, 1], so that the integral I exists only as a cauchy principal value. The weight function w(x) also has a singularity of the form log|x-b| at some point b which may be different from the point a.

Could anyone advise a good quadrature rule to deal with the cauchy principal value integral?
Is such a quadrature rule sensitive to the presence of other integrable singularities?
Maybe this rule is implemented in some Julia package?

Thank you in advance.

You can use SingularIntegralEquations.jl for Hilbert and log integrals:

julia> using ApproxFun, SingularIntegralEquations

julia> x = Fun()
Fun(Chebyshev(),[0.0, 1.0])

julia> hilbert(exp(x), 0.1) # 1/π * ∫ exp(t)/(t-x) dt
0.6363129376196711

julia> logkernel(exp(x), 0.1) # 1/π * ∫ exp(t)*log(abs(t-x)) dt
-0.7386375437288631

julia> hilbert(exp(x)/sqrt(1-x^2), 0.1) # works with Chebyshev singularities
1.1404096104609651
3 Likes

As explained here, you can do Cauchy principal-value integrals easily using any integration package (e.g. QuadGK.jl) using a simple transformation that eliminates the singularity: Cauchy principal value integrals by fmeirinhos · Pull Request #44 · JuliaMath/QuadGK.jl · GitHub

1 Like

Can SingularIntegralEquations.jl deal with the situation when there are both 1/x and log|x| singularities present? Also, is it possible to get the points and weights of the quadrature separately (I know f(x) only numerically, although I have an analytical estimate for the singularity)?

It’s not really a “quadrature rule” as it first expands in a basis then evaluates the relevant singular integrals exactly. If you have data you can use a least squares fit: Frequently Asked Questions · ApproxFun.jl

It might be possible to get a quadrature rule out of this a la Clenshaw–Curtis


Thank you.

In my case, I would prefer not to evaluate x*f(x) at singularity.
I suppose, I should separate a symmetric interval centered on 1/x singularity and construct a symmetric quadrature on this interval. Then, the odd part of the integrand won’t contribute to the integral.

However, there is the problem stemming from the fact that the integrand looks like f(x)*w(x), and the weight function has its own singularity. In normal situation, I can account for the singularity of w(x) by constructing the appropriate gauss quadrature with this weight function.
So, I either suit-tailor the quadrature to the w(x), but it is non-symmetric with respect to 1/x singularity, or I construct symmetric quadrature for unity weight (f(x)*w(x) = (f(x)*w(x))*1), but it won’t be very accurate in the vicinity of the singularity of w(x).

In the situation when the singularities of f(x) and w(x) are well separated, I could treat them separately. What can be done in the situation when they are close to each other?

In any case where you have multiple singularities, you must typically break the integral into pieces at each singularity. (Assuming you don’t want to build a special quadrature rule that takes all of your singularities into account. This is possible, but is probably more trouble than it is worth unless you are evaluating a lot of integrals for exactly the same singularities.)

(With QuadGK.jl, you can pass multiple integration segments to quadgk in a single call, which is useful because that means that adaptivity is based on a global estimate for the error of the whole integral, not on the errors of individual segments treated independently.)

5 Likes

2 posts were split to a new topic: Relative tolerances for near-zero integrals