Cauchy singular integral equation with variable coefficients

Singular integral equations with variable coefficients and Cauchy type kernel known as Carleman type singular integral equations (CTSIE). These equations arise in the integral equation formulation of a large class of mixed BVPs in science and engineering, especially when two-dimensional problems are considered. These equations also appear in contact problems in fracture mechanics, mainly crack problems in elasticity (Muskhelishvili, 1953; Gakhov, 1966; Erdogan et ah, 1973; Duduchava, 1979; Estrada and Kanwal, 1987; Peters and Helsing, 1998). Here a(x), ш(х), and f(x) are known input functions defined on (0,1) and a is unknown to be obtained. The integral is in the sense of CPV. The numerical solution of CTSIE was investigated using a number of methods such as dominant equation method by using orthogonal polynomial approximation (Dow and Elliott, 1979), Galerkin method (Ioakimidis, 1981), piecewise-polynomial method (Gerasoulis, 1986), collocation method (Miel, 1986; Junghanns and Kaiser, 2013), method based on Riemann-Hilbert problem of complex variable theory (Estrada and Kanwal, 1987), Jacobi-polynomial method (Karczmarek et ah, 2006), Taylor series expansion and Legendre polynomial method (Arzhang, 2010), so on.

In the previous section a method based on LMW has been developed to approximate solution of SIEs with Cauchy type kernel and constant coefficients. Illustrative examples considered there have solutions bounded within the domain of definition. However, in many physical problems, viz., the contact problems in fracture mechanics, integral equations with Cauchy singular kernel arise whose solution may be unbounded at the endpoints. It is difficult to approximate solution to such problems by applying the method developed in the previous section in a straightforward way.

In this section, a scheme based on LMWs proposed in Paul et ah (Paul et ah, 2018) will be discussed to get approximation of solution of the IE ( even if their solutions are unbounded at the end points of their domain. In the process of development of the scheme evaluat ion of integrals for obtaining the MSR of integral operator associated with CTSIE will be discussed. Subsequently, MSR of integral operator involving the product of functions in Holder class and the Cauclqy singular kernel have been derived. These results have been used in the next subsection to transform the IE ( into systems of linear algebraic equations to be solved numerically. An estimate of the Holder exponent of the solution at the end points is proposed. Error estimation of the proposed method is also presented. The numerical scheme has been tested through its application to a number of examples arising in mathematical analysis of contact problems in fracture mechanics.

Evaluation of integrals involving function, Cauchy singular kernel and elements in LMW basis

We write Carleman type singular integral equation in the operator form as where

For the evaluation of the MSR of the integral operator wK-c ■ we have to evaluate the CPV integrals involving product of elements of basis and their images under loJCc- Here the known function ui(x)

in Holder class may be regular or singular, e.g., to(x) = wa>p(x) = xa(1 — x)&, (a,/3 > —1), known as Jacobi weight. The method to evaluate such integrals is now discussed in some details. As in previous section singular integrals, whenever they appear, are in the sense of CPV. We first consider the case: b(x) = uj(x),g(t) = 1.

Scale functions

V 'e use the notation

where фг(х) is a polynomial of degree i, which is defined in ( The explicit expression of J) t-x ^ can eas'br calculated by the property of CPV integral. If we denote

then its explicit form for К = 5 can be found as

Note that each of the functions Cx), 1 = 0, ...,4 has logarithmic singularity at, x = 0 and x = 1. Now by using notations ( in (, we get

The values of PuKc(h,h), h,h = 0,1 ,..,K — 1 for given uj(x) can be obtained either integrating analytically or by using a suitable quardrature rule. Explicit forms of these integrals for Jacobi weight Lu'n f-i(x) = xa(1 — x)P,a, /3 > —1 with К = 5 can be obtained by introducing the expression

into the formula

where Hn denotes the nth harmonic number. Г is the gamma function (Olver et ah, 2010) and is the coeffcient of xm in ф'л(х). For other forms of ш{х) (nonsingular positive semidefinite), values of ptilKc(li,h) can be obtained by using suitable quadrature rule.

Scale functions and wavelets

Here we use the notation,

If we denote

then by using an appropriate transformation, we obtain

where v = 2:lx — k. The integral is defined in the sense of CPV whenever v € (0,1).

Now, using the relation ( in ( for j = 0 and after some algebraic simplifications we


By using ( in (, we get the formula for a(li;l,j,k) as

Two integrals in ( for ш(х) in polynomial form can be evaluated explicitly interms of j. k. l,l for all admissible values of к for given j. For other forms of w(a;), a suitable quadrature rule may be used to evaluate the integral in (

Wavelets and scale functions

Here we use the notations

Now by using the relation ( in (, we get

The integral in ( when ш(х) in polynomial form can be evaluated explicitly in terms of j, к, l, /2 for all admissible values of к for given j. For other forms of ш(х), we use appropriate library functions available in MATHEMATICA or in other software to evaluate the integral in (


W 7e use the notation

which involves integration of the product of wavelets in LMW basis. Rescale of variables x. t followed by translation leads to

In ф>2 (2J2-Jl t — A?2 — 2-b-A/c1), if we write then, ( can be recast into the form

Now, we use the transformation t = r + r in the above integral and get

Due to the finite support of ipjl2_j1 Г2 (r), the above formula can be recast into the form

Given hi, the values of rappearing in the RHS are unique and can be obtained by using the formula ( For rq = 0, the value of the integral can be obtained by using suitable library function available in MATHEMATICA, or any other softwares like MATLAB. MAPLES, etc. while for rq = 0, the integral can be written as

By using ( and (, we get

For 72 < Л, where,

For jq = 0, the value of the integral can be obtained by an appropriate library function mentioned above, while for jq = 0, the integral can be written as

By using ( and (, we get

Thus, all integrals p„Kc{h,h) otsjKc(hh,j,k),ш)Со (fi,7, М2), (*b7i. *1; I2J2, h) are now

evaluated for the IE in the form Eq. ( These will be used in subsequent sections for the

MSR of the integral operator ojKc- We now denote matrices with their dimensions (depending on resolutions)

In case of the IE ( with b(x) = 1, g(t) = w(t), the matrix elements (h, h), (h', h,j, fc),

f3ulc k: I2), ywKc (/1 • Jb hi: 12 J2 •> ^2) for ojK-c can be obtained by using their corresponding values for the equation through the following rule.

Evaluation of the integrals involving product of a(x), scale functions and wavelets

If we denote









In definitions (, Zj, I2 = 0, 1, .... К — 1.

Multiscale representation (regularization) of the operator ωKC in LMW basis

To obtain MSR of uilCc, we write


The integrals in ( are defined in the sense of CPV whenever x € supp Then the


of uiJCc in the basis (Фц, (,/_[)Ф) can be written in the form

where submatrices P^K(:, <*ыКс, > /Я^с > 7„Kc are given hi previous subsection.

Multiscale approximation of solution

For the IE (, we assume that the input function / € L2[0,1]. The solution u(x) of the IE ( is also L2[0,1] so that it has multiscale expansion ( in LMW basis. Using the MSRs ( for Aj and ( for udCc, the IE ( with b(x) = А иj(x),g(t) = 1 or h(x) = 1, g(t) = A Lo(t) (Л € R) can be recast into a system of linear algebraic equations given by

Components a0, (j_i)b of coefficients of multiscale expansion of the unknown function u(x) in the LMW basis at resolution J are unknown. Components Co, (,/_itd of coefficients of the multiscale expansion of f(x) in the same LMW basis at the same resolution J are known.

The unknown coefficient ao. (,/_i)b can be found as

The matrices A j in ( and ojICqj in ( are of order (2JК) x (2JK) and are sparse in general for a(x), ш(х) € L2([0,1]). Due to the inherent structure of multiresolution approximation of functions and operators, the coefficient matrix Aj + Л ojICq j involved in Eq. ( is in block form. In contrast to the Galerkin approximation in the basis of orthogonal functions with support [0,1], the algebraic equations for unknown coefficients a0, (j^i)b can be split into block by block. Since the matrices in the diagonal blocks7Q [j, j)+X'yUICc (j. j) (presented in ( and ( for j = 0,1,.., 2 ' — 1 are well-conditioned for а(х),ш(х)L2([0,1]) and Л € t. solutions ao. (J-i)b can be obtained by inverting the matrices pa +Л p„Kc , an [j,j) +A atulCc (j, j), {j,j)+А /ЗшКс (j, j),

7a Chi) + ^ 7uijcc (i;i) efficiently in place of inverting the full matrix (Aj + A ujIC^Ij).

At this point it may be worthy to mention that the rate of convergence of the approximate solution excluding the rapidly varying part xa or (1 — x)$ with respect to the resolution is of exponential order. Consequent ly, for physical problems, expected accuracy can be acheived by the resolution J < 5 so that the algebraic equations can be solved efficiently by using the numerical solver now available and thus J need not be chosen to be large.

The multiscale approximate solution of the integral equation ( at level J is given by

where the coefficients al0 0, fe(l = 0,1,......К — 1, к = 0,1,......2-*-1, j = 0,1,... J — 1) are

obtained from (

Estimate of Hölder exponent of u(x) at the boundaries

The behaviour of м(а;) at any point in the dyadic interval Ij,k (= зг-]) of [0,1] can be estimated

in terms of the wavelet coefficient fc by the following theorem.

Theorem 6.8. If Vj,k stands for the Holder exponent v of u(x) in Ij.k(= (з> then an

estimate for Vj at the endpoints is given by (Section 3.8.2]

Estimation of error

In the processes for obtaining numerical solution, the CTSIE ( has been reduced to a system of linear equations.

Here Uj,s. given by (, is the MSA of u(x). If we denote then the L2-en-or(a-posteriori) of the solution Ujls is given by

Although hl- fc’s, j > J appearing in the formula ( are not known, the sum in the RHS can be approximately obtained in terms of sum of squares of wavelet coefficients of previous two resolution, viz., J — 2, J — 1 in the formula

so that a qualitative estimate for ||<5u||^2[o,i] can be found as where

depending on bj fc’s, which are known.

Applications to problems in elasticity

To test the efficiency and domain of applicability of the approximation scheme in LMW basis developed here, we consider two specific examples from problems in elasticity theory.

Example 6.9. In contact problem

Here we have considered the contact problem corresponding to a rigid punch with sharp corners sliding on a half plane studied earlier by Miller and Keer (Miller and Keer, 1985) and Jin, Keer and Wang(Jin et ah, 2008). The equilibrium state of the system can be described by the IE

with the load condition (constraint)

Here, u(x) is the unknown function related to the stress intensity factor, //. is the coefficient of friction, cr is the Poisson ratio and L is the external load.

The exact solution to Eq. ( without the load condition ( can be found in the handbook (cf. (Polyanin and Manzhirov, 2008, pp. 765)) of integral equation compiled by Polyanin, and Manzhirov as


For numerical computation, ц and a are chosen such that ,5* = —0.34. Then the exponents of u(x) near x = 0 and 1 are found from ( as ±0.34 which are very close to the approximate values of the Holder exponents estimated by using the formula ( It is important to mention here that u°(x) ss a;-0'66(l — a:)-0'34 in /го-class (Karczmarek et al., 2006), satisfies the homogeneous Ecp ( So, the solution ( of Eq. ( with /3* = —0.34 is not unique in the sense that whenever ф is a solution of (, и + c uo is also a solution for any choice of c. The arbitrariness of c can be fixed by the load condition ( so as to provide the unique solution of ( in /го-class satisfying condition (

Since numerical solution of algebraic equations exists when such solution is unique, we choose MSA of unknown solution u(x) as

where ui(x) = (I — гг)" with —1 < a, /3 < 0 and гг + 8 = —1. Substitution of ( in

Ecp ( and ( gives system of linear algebraic equations ( and

In order to solve this system of К2J+1 + 1 equations in К2J+l variables which are components of ao, ,/_ib. one equation in ( has been replaced by ( and then the resulting system is solved by using standard procedure. The consistency of the solution has been tested by their substitution into the equation replaced. The difference has been found to be О 10-7 .

For L = 2 ; the value of the components aj are presented in Table 6.6 while the value of

the wavelet coefficients, i.e., the components of 0bT are found to be identically zero. Use of the results presented in Table 6.6 in ( shows that the approximate solution ( provides the exact solution

Table 6.6: Values of the coefficients of scale functions in the multiscale approximate solution ( of Ecp ( (in case of К = 4).





Example 6.10. Crack problem

Here we have considered the CSIEFK (Jin et ah, 2008; Theocaris and Ioakimidis, 1977) with the condition

appearing in the study of crack problems along the interface of two plane isotropic elastic media.

For solving the above IE numerically, we take here 7 in such a way such that u> = ^= 0.1. The forcing term f(x) is choosen as

Case 1. f(x) = 1,

Case 2. f(x) = 1 + (2x - l)3 + 3(2x - l)4,

Case 3. f(x) = 1 + (2x — l)8.

Using MSR Aj in ( and in ( for w(ai) = 1, J = 4 and К = 4, we get

( ^

the system of Eqs. ( for . Also the condition ( transformed into

< зЬт ,

Then the solution is obtained by solving the system of linear equations ( and Eq. (, i.e., one equation in ( is replaced by Eq. ( From the wavelet coefficients of the solution we calculate the Holder exponent near the end points by using the formula ( The estimated values v3 0 and 1/3.7 are presented in Table 6.7 for different choices of f(x). The values of и for all f(x) are converging to the number —5 + 0.1 i near x = 0 and —^ — 0.1 i near x = 1. This estimated results completely agree with the theoretical values of the Holder exponents. Now following the transformation

where ш(х) = xa(1 — x)@ with a = — + 0.1 i and ,3 = — | — 0.1 i, the IE ( is converted into

with the condition (cf. Eq. (

Table 6.7: Values of the Holder exponent of the solution of Ecj. (, condition ( for forcing terms f(x) presented in cases 1-3.

in Case 1

in Case 2

in Case 3

-0.499 + 0.1 i

-0.497 + 0.1 i

-0.497 + 0.1 i

-0.499-0.1 i

-0.497-0.1 i

-0.497-0.1 i

Comparing Eq. ( with Eq. ( corresponding to the second case of Eq. ( for unknown v, we get a(x) = —7 i oj{x), A = j. We have solved the system of linear algebraic equations ( corresponding to Eq. ( with additional equation

The approximate solution in LMW basis have been found by using numerical solution of algebraic equations ( and compared with the exact solutions given below.

Case 1.

Case 2.

Case 3.

This solutions have been obtained by using formula Eq. (

The absolute error is found to be zero in case 1. The pointwise absolute errors in approximate solution in cases 2 and 3 at resolution J = 4 and = 4) LMW basis have been presented in Fig. 6.6 and Fig. 6.7. The error estimate (in sup norm) for real as well as complex part of the approximate solution in cases 2 and 3 are presented in Table 6.9. The error is is found to be of order ()(10-4). It is observed that the absolute error decreases with increase in the resolution)J).

The formula

for the stress intensity factor used earlier by Theocaris and Ioakimidis (Theocaris and Ioakimidis, 1977) have been used to get their approximate value from the approximate solution obtained here. Those values have compared in Table 6.8 with the results obtained by other methods. From this table we conclude that the present method predicts the result with same accuracy as obtained by (Jin et al., 2008) but with better accuracy than the result given in (Theocaris and Ioakimidis, 1977).

Table 6.8: Comparison of the stress intensity factor.


Case 1

Case 2

Case 3


1 - 0.2 t

2.4069 - 0.8043 i

1.2488 - 0.3243 i

Method in

(Theocaris and Ioakimidis, 1977)

1 - 0.2 t

1.7035 - 0.51167 i

1.2488 - 0.324 i

Method in (Jin et al., 2008)

1 - 0.2 t

2.4069 - 0.8043 i

1.2488 - 0.324 i


1 - 0.2 t

2.4069 - 0.8043 i

1.2488 - 0.324 i

Plots of pointwise absolute error of the solution of Eq. ( for Case 2

Figure 6.6: Plots of pointwise absolute error of the solution of Eq. ( for Case 2: i) real part of the error (left) ii) imaginary part of the error (right).

Plots of pointwise absolute error of the solution of Eq. ( for Case 3

Figure 6.7: Plots of pointwise absolute error of the solution of Eq. ( for Case 3: i) real part of the error (left) ii) imaginary part of the error (right).

Table 6.9: Errors (in sup norm) corresponding to the approximate solution Vjls(x) of Eq. ( J = 2,3,4.


Case 2

Case 3


1.29 x 10-2

1.83 x 10-2

3.38 x 10-2

4.42 x 10-2


1.74 x 10-3

1.81 x 10-3

6.83 x 10-3

6.35 x 10-3


1.12 x 10-4

1.24 x 10-4

7.69 x 10-4

7.28 x 10-4

Example 6.11. We consider here the integral equation



The exact solution can be found as

We have evaluated ao, b,(0 < j < 2) by using the formula ( for this problem. It is found that all the wavelet coefficients (components of bj, j = 0,1, 2) are zero. The values of the components of ao are presented in Table 6.10 and found to coincide with the coefficients of the multiscale expansion of the exact solution given in (

Table 6.10: Values of the coefficients of scale functions in the multiscale approximate solution u^,s(x) of Eq. ( for К = 5.

  • 2
  • 3


  • 2
  • 3/5



Example 6.12. Variable coefficient with nonsmooth and logarithmic singular forcing term We consider IE



The exact solution can be found as

nonsmooth at x = |. The approximate solution in the LMW basis corresponding to К = 5 at resolution J = 3 have been obtained. The value of the coefficients of a0 and b0 in the multiscale approximation of the solution are presented in Table 6.11 while the values of the wavelet coefficients, i.e., the components of b; (1 < j < 3) are found identically zero. Therefore, the multiscale approximation of the solution in the LMW basis corresponding to К = 5 is reproducing the exact solution to Eq. (

Table 6.11: Coefficients of Фо, Фо.о hr the multiscale approximate solution u^s of Eq. (

  • 1
  • 2





  • 1
  • 16


  • 1
  • 8/ТУ





< Prev   CONTENTS   Source   Next >