Probl. Anal. Issues Anal. Vol. 8 (26), No 2, 2019, pp. 51-66 51
DOI: 10.15393/j3.art.2019.5570
UDC 519.632.4
P. N. IVANSHIN, E. A. SHIROKOVA
THE SOLUTION OF A MIXED BOUNDARY VALUE PROBLEM FOR THE LAPLACE EQUATION IN A MULTIPLY CONNECTED DOMAIN
Abstract. Here we apply the Cauchy integral method for the Laplace equation in multiply connected domains when the data on each boundary component has the form of the Dirichlet condition or the form of the Neumann condition. This analytic method gives highly accurate results. We give examples of applications of the method.
Key words: Cauchy integral, Laplace equation, mixed boundary value problem, multiply connected domain, approximate solution 2010 Mathematical Subject Classification: 35J57, 30E10
1. Introduction. The article extends the results of [1], where the approximate analytical solution of the Dirichlet problem for the Laplace equation in simply connected and doubly connected domains with smooth boundaries was reduced to systems of linear algebraic equations and the solution had the form of the real part of a Cauchy integral. The Laplace equation arises in different areas, such as electrostatics (where it describes the electrostatic potential), stationary potential incompressible fluid flows, and steady-state heat conduction [12]. We consider the solution of the Dirichlet problem as the base for constructing solutions of the Neumann and the mixed problems. There exist three main methods of solving the Dirichlet problem: the Fourier method and the Green-function method for domains of special types; the method of double layer logarithmic potential with the density that is the solution of a Fredholm integral equation [5,6]; numerical methods (the finite difference method, the finite element method, the boundary element method [14]). The numerical finite difference method of solving of Dirichlet problems for irregular domains meets difficulties connected with construction of adequate difference schemes for such domains [7]. The solution of a Dirichlet problem by the
© Petrozavodsk State University, 2019
numerical finite element method for irregular domains meets the difficulties connected with an adequate discretization of the domain. However, recently some boundary value problems for the Laplace equation have been solved by the finite element method applying modern computer approches of automatical domain discretization. Unfortunately, these solutions converge only in Lp spaces [8-11].
Here we present the Cauchy integral method applied to the mixed boundary value problem in the multiply connected domain. We reduce a mixed boundary problem to a linear integral equation. In order to do this, we apply the technique of the Cauchy integral and analytic functions. A number of articles are devoted to the study of these equations from the perspective of their reduction to infinite linear systems [16-19]. We apply the method based on the Fourier expansions of the kernels of integral operators. The similar technique was applied in [2,3] for construction of the conformal mapping. The main idea of the proposed solution is to reconstruct the boundary values of a conjugate harmonic function. We give an analytical approximate solution of the mixed problem through the Fourier coefficients of the conjugate harmonic functions. The constructed approximate solution of the mixed problem is the real part of an analytic function i (i. e., the real part of a Cauchy integral plus certain logarithmic summands), so it is possible to investigate this solution via differentiation. 2. The main results. Consider a multiply connected domain D,
n
dD = U Ck. Assume that ED. Let the boundary components of the k=0
domain be given as Ck: zk (t) = xk (t)+iyk (t), t E [0, 2n], k = 0, ...,n. Here C0 is the exterior boundary component. Assume that zk(t) is differentiable and zk(t) = 0, t E [0, 2n], Vk E {0,..., n}. We pass each component of the boundary dD increasing the parameter t so that the domain D remains on the left.
The mixed problem,. Let us find the function / such that
Af(x,y) = 0, (x,y) E D,
given the boundary values /1Ck = Uk(t), t E [0, 2n], for k = 1,..., m, and
df
the normal derivatives —— |C, = Vk(t), t E [0, 2n], for k = m + 1,..., n.
dN k
2n
Here the compatibility conditions J Vk(t)|zk(t)|dt = 0, k = m + 1,..., n,
0k
should hold for the boundary components where the Neumann conditions
are given [20, Ch.12, Example 5, p. 618]. Note that the boundary condition on the exterior boundary component L0 may be either of the Dirichlet form or of the Neumann form.
Definition 1. We call the mixed boundary value problem with the Dirichlet condition on the exterior contour the mixed boundary value problem of the first type (MBVP1).
Definition 2. We call the mixed boundary value problem with the Neumann condition on the exterior contour the mixed boundary value problem of the second type (MBVP2).
Theorem 1. The mixed boundary value problems of both types for a domain with the smooth boundary components and the boundary data of the class C2a have unique solutions that can be approximately obtained via the Cauchy integral method.
Proof. We give the constructive proof of the theorem. The differential class of the boundary data ensures convergence of the approximate solution to the exact one.
The uniquenece of the solution of the mixed problem is easily obtained as the consequence of the Stokes formula.
Structure of the solution depends on the type of the boundary data on the exterior component Co.
We base our solution on reconstruction of the boundary values of the harmonic functions conjugate to the functions defined by the boundary conditions on the contours Cj, j = 0,... ,n. The solution f(x,y) is the real part of a single-valued Cauchy integral and a sum of logarithms as in [4]. In order to construct the solution, we fix the points z* inside the domains bounded by the contours Cj, j = 1,... ,m.
Let us consider the MBVP1. Consider the harmonic in D function
m
p(x, y) = f (x,y) - E Ajln |x + iy - z*|. Here j=i
p(x,y)lcl = pi(t) = Ui (t) Aj ln | zi(t) - z* |, l = 0 ,...,m.
j=i
Denote by q(x, y) the harmonic in D function conjugate to the function
t
m
p(x,y). Then q(x,y)lol = qi(t) = Vi(t)lzl(t)ldt+BiAjarg(zг(t)-zm),
o j=i
l = m + 1,..., n. Later we will find the real numbers Aj, j = 1,..., m, l = m + 1,..., n.
For the MBVP2, consider the harmonic in D function
m—1 m— 1
j •
p(x y) = / (x y) - X] ln |x ± ^- j±ln |x ± ^- zm ^ • j=i j=i
m— 1 m—1
Here p(x,y)|c; = Pi (t) = u (t) Aj ln |zi(t) - z*| +ln |zi (t) - zm | ^ Aj,
j=i j=i
l = 0,..., m.
Denote by q(x, y) the harmonic in D function conjugate to the function p(x,y). Then
q(x,y)h = qi(t) =
t i i
/m—1 m—1
VKt)|zi'(t)|dt + Bi - ^ Ajarg(zi(t) - z*) + arg(zi(t) - zm) ^ Aj, o j=1 j=1
l = m+1,..., n, 0. Later we will find the real numbers Aj, j = 1,..., m-1, Bi, l = m +1,...,n,0.
In both cases p(t) + iq(t) represents the boundary value of an analytic in D function if and only if [4]
1 i p(t ) + *q(T)
P(t)± iq(t) = 1/Pfi±iTtT dz (t )
ni j z(t) - z(t)
z(t) - z(t)
dD
We separate either the real or the imaginary unknown part of these boundary values on the boundary components and get the following Fredholm system of equations.
For the MBVP1, we consider the imaginary part of equation (1) for t e Ck, k = 0,... ,m:
(n
- ^ iPj(t)Re[dr ln(Zj(t) - Zk(t))]dr± j=00
n 2n x
± / qj(T)Im[dr ln(Zj (t) - Zk(t))]dr). j=o 0 /
Similarly, for t G Ck, k = m + 1,...,n, we consider the real part of equation (1)
1
2n
Pk(t) = - J] qj(r)Re[5r ln(zj(r) - zk(t))]dr+
3=0 0
n 2n x
+ It /Pj(r)Im[dr ln(zj(r) - zk(t))]drV
j=0 0
For the MBVP2, the equation for t G C0 is the equation of the form
/ n 2n
Po(t) = - ( ^ / 3(r)Re[dr ln(zj(r) - zo(t))]dr+
j=0 0
n 2n \
+ It /Pj(r)Im[dT ln(zj(r) - zo(t))]dr).
j=0 0
These equations constitute system (1) of n + 1 Fredholm equations of the second kind.
For the MBVP1, we search for the solution in the form of the vector-function (qo(t),... ,qm(t),Pm+i(t),... ,pn(t)):
L0,n \ L1,n
(I + Ko, o Ko, i ... Ko, m Lo,m+i
Ki, o I + Ki, i ... Ki,m Lo,m+i
Km, o Km, i . . . I + K m, m Lm, m+i
Lm+1, o Lm+i, i ... Lm+i, m I — Km+i,m+i
\ Ln, o Ln, i ... Ln, m Kn, m+i
1 qo \
Lm, n Km+i,n
I Kn, n )
x
qm pm+i
\ Pn J
Wo, 0 \ m I Wj, 0
+ E Aj
Wo,
j=i
W
j, n
n
For the MBVP2, we search for the vector-function (q^t),..., qm(t), Pm+i(t),... ,Pn(t),Pb(t)) :
\
/ I + Ki, i . . Ki,m L0,m+i . . Li,n Li, 0
Km, i . 1 + K m, m Lm,m+i . . Lm, n Lm, 0
Lm+i, i . . Lm+i,m 1 — Km+i, m+i . . Km+i,n — Km+i, 0
Ln, i . . Ln, m Kn, m+i . 1 — K — n, n —Kn, 0
V L0, i . . L0, m — K0, m+i . . —K0,n 1 — K0, 0
x
x
qi
qm
'Pm+l Pn
V Po )
( wo, l \
Wo,n V wo, 0 )
mi
( Wj, i \
+ £ Aj
j=i
w
J,n
\ Wj, 0 /
The right-hand sides of equations (2) and (3) contain the linear combinations of the known functions that are either reconstructed from the boundary conditions or derived from the real or imaginary parts of ln(zj (t) — —zk(t)), j = 0,..., m, k = 1,..., n, where j = k.
Here we have the operators
2n
g(t) = Im[5r ln(zj (t ) — zfc (t))]g(r)dr,
2n
Lj,fcg(t) = J Re[5r ln(zj(t) — zfc(t))]g(T)dr. o
The second operator is a principle value integral defined for the given density functions of a Holder class. Systems (2) and (3) can be solved in the space (L2[0, 2n])n+i by reducing them to infinite systems over Fourier coefficients of each unknown function. If the unknown function has the
form rj (t) = Dj + Y1 Oj, k cos(kt) + k sin(kt), j = 0,..., n, then the fc=i
values Oj,k, ,5j,k, j = 0,... ,n, k G N, can be obtained from the infinite
system. We find the values of Dj later. For the MBVP1, we fix D0 = 0,
2n
and for the MBVP2, we choose Bo, so that f qo(t)dt = 0.
o
Approximate solution of the infinite system over the Fourier coefficients is the solution of the truncated system.
Assume that each contour Cj is given as the polynomial
Nj
Mj
,(t) = 5] Chk eikt + £ j -k e-ik\ j = 0,...,n
k=l
k=l
Note that
. t - t\ .T + t
ln[zj(t) - zj(t)] = ln(2i) + ln ( sin +
So,
Nj
+ ln £ Cj,keiktY, ei'
k-l
ikt\^ eil(r-t)
k=l
l=0
Mj k-l
E Cj,-ke-ikr E eiltyT-t) k=l l=0
Red ln(zj(t) - zj(t))] = 2cot (+
+Re
Nj
dT ln £ c3,keiktj2
k-l
ikt\^ eil(T-t)
k=l
l=0
Mj
E
k=l
k-l
Cj,-ke-ikT E eil(T-t) l=0
and
Im[dT ln(zj (t) - zj (t))]
—+ Im 2
Nj
k-l
6t ln £ Cj,keiktj2
il(T-t)
k=l
l=0
Mj
- Cj,-ke k=l
kl
-ikT\^ ea(T-t)
l=0
Solution of systems (2), (3) reduces to solution of the finite systems over the coefficients of the Fourier polynomial expansions of the unknown functions qo,... ,qm, Pm+i, ...,Pn or qi,... ,qm, Pm+i,... ,Pn,Po for the MBVP1 or MBVP2, respectively. It is possible to replace the infinite system with the finite one, due to the lemma similar to the following statement (see also Lemma 1 of [2] for the (n + 1)-dimensional case).
z
2
T
Lemma 1. Let Zj (t) = eipt be such that zj (t) = 0, j = 0, ...,n,
p=-t
t G [0, 2n], and let Uj (t), Vj (t) G C2,a[0, 2n], j G {0,1,...,n}. Then the approximate solution of the uniquely resolvable Fredholm system of integral equations of the second kind (2) or (3) can be reduced to solution of a finite linear system of rank F with the error estimated by O(1/F2).
Proof. The solution y- (t) of this system of equations can be represented
in the form of the trigonometric series: y-(t) = ak, j cos kt + j sin kt.
k=i
2n
Restrictions J yj (t)dt = 0, j = 0,...,n, ensure the uniqueness of the 0
solution of the system (1). Now this system of equations can be reduced to the infinite system of equations over the finite number of the unknown coefficients ak, j, , j = 0,..., n, k = 1,..., ro; in the equivalent form
Y = PY + Q, (4)
where Y = (ai,o,ai, 1,..., 0,^1,1,... ,^i,n,«2,o,...) G /2, the infi-
nite matrix P consists of the elements
2n 2n
"2 f (mt)dt / g(PT)(arg[zk(t) - zj (t)])Tor
2n 2n
-^J f (mi)d^g(pr)(ln |zfc(t) - Zj(t)|)Tdr 0 0
with f (x), g(x) equal to cos x or sinx. These elements are the Fourier coefficients of the double Fourier series of (arg(zk(r) — Zj (t)))T or ln |zk (r) — —Zj (t)|T, Q constitutes the sequence of the corresponding Fourier coefficients of the functions Uj (t), j = 0,..., m, and Vj(t), j = m + 1,..., n, in Case 1, and Uj (t), j = 1,..., m, and Vj(t), j = m + 1,..., n, 0, in Case 2. We need to construct the approximate solution y (t) of system of equa-
M
tions (4) in the trigonometric polynomial form y (t) = ak,j cos kt +
fc=i
+ sin kt in order to apply the truncated linear system approach as in [2]. So, we have to find the vector YF with zero coordinates starting with the (2M(n + 1) + 1)-th one, which approximates the infinite vector Y. Further on, we identify the vector-function Y, the integral operator P,
the vector-function Q with the sequence Y, the infinite matrix P and the sequence Q, respectively.
Evidently, the kernels (arg(zj(r) - zk(t)))'T of system (1) are infinitely differentiable for k = j. Due to the Cauchy Theorem, we have arg(zj (r)--zj(t)) = argz'(t + 9(r - t)) with 9 G (0,1); so this function is well-defined for t,r G [0, 2n] and limarg(zj(r) - zj(t))'T = k(t)|z'(t)|/2, where
T ^t
k(t) is the curvature of the boundary curve at the corresponding point. It can be easily verified that the kernel (arg(zj(r) - zj(t)))'T is at least twice differentiable with respect to both variables. So the double complex Fourier coefficients of (arg(zk(r) - zj(t)))'T have the following estimates:
I I U \Ck,j,l,t\ < |l|2|t|2 .
The singular kernels (ln \zj(r) - zj(t)\)T are the kernels of the integrals with the known Holder densities. So, these integrals are contained in the system of equations as the free summands.
For F = (n + 1)2M, integral system of equations (1) reduces to infinite linear system (4), which can be presented as follows:
(If - Pf S \ (Yi\ = (Qi\
V R I~ - y)\Y2j \Q2) .
(P0,0 ... P0,n\
......... I, M x M
Pn, 0 ... Pn, n)
matrices Pj,k correspond to the integral summands of (1), j,k = 0,..., n, S is an F x ^ matrix, R is an to x F matrix, V is an to x to matrix, IF, and I^ are the identity matrices of the relative sizes. Each of the vectors Qi and Yi has F coordinates, the vectors Q2 and Y2 have the infinite number of coordinates. The Fourier coefficients of the smooth functions tend to zero as their numbers tend to infinity; so the coefficients of the matrices S, R and V together with the coordinates of Q2 decrease rapidly as F ^ to. Due to the Lemma 1 assumptions and the Fourier coefficients speed of convergence to zero, the matrix norm of V and the vector norm of Y2 tend to zero as F ^ to.
Let us prove that there exists the number T G N, such that the matrix operator IF - PF is invertible VF > T, since the limit for PF integral operator P is compact and the operator I - P is invertible, due to the lemma assumption. Note that we do not distinguish a finitely dimensional
vector and the Fourier polynomial with the corresponding finite set of coordinates in our proof. Recall first that due to chapter VI, paragraph 1 of [21], ||P — PF|| — 0 if F — to. The operator norm that we deal with here is the usual operator norm for the Hilbert space mappings. Let us assume that VT e F there exists sj > T, such that the spectrum of PSl contains 1. Then there exists an infinite sequence )1eF C L2, such that ||vSl || = 1 and PSl= . Let us prove that then there exists at least one limit point for the sequence }1eF. Since the operator P is compact, there exist both a subsequence {vsfc. }jeF and an element wo e L2, so that PvSK. — wo as j —y to. Then
||PsK Vsk, — Wo || = ||Psfc. VSK. — + — Wo || ^
Kj Kj Kj Kj Kj Kj
< 11 PSKj — PvSKj 11 + 11 PvSKj — Wo 11 ^ ^ ||Psk,. — P|| + ||PVsk,. — Wo| — 0, j — TO.
Thus, ||vsKj — Wo|| = ||PsKjVsKj — Wo|| — 0, j — TO. Hence, vSK. — Wo, j — to. Note that, since ||vSK. || = 1, Vj e F, the element Wo is nondegenerate. Let us show now that the relation Pwo = Wo holds. Indeed, we have
||PWo — Wo || = ||Pwo + PvsKj — PvsKj — Wo || ^
^ ||P HH Wo — vSKj 11 + 11 PvSKj + PSKj vSKj — PSKj — Wo 11 ^
^ ||P ||||wo — vskj ! + |PSKj — Pvskj ! + ||wo — Pskj VsJ ^
^ ||P|||wo + VsKj || + ||P — PsKj || + ||Wo — VsKj ||— 0, j — TO.
Hence, the spectrum of P contains 1. This contradicts with the uniquiness of the corresponding mixed problem.
We now take the number F so that ||V|| < 1 and the matrix /F — PF has the inverse. Now we have the relation
(If — Pf )[/f — (If — Pf )-1S (/^ — V )-1R]Y1 = Q — S (/^ — V )-1Q2.
Obviously, one can choose the value of F so large that ||S(/^ — V)-1R|| = = O(1/F2) = r, where r < 1 is arbitrary small due to estimates of the corresponding Fourier coefficients. Now we estimate the norm of the difference between the solution Y1 and the solution Y1 of the truncated system (/f — Pf )Y = Q1:
||Y1 — Y11| ^ -^r ||(/f — Pf )-1||||S (/« — V )-1|nQ2n +
1r
+ II (If - pf)-
1 — r
Consider the first summand on the right-hand side of the last inequality. Recall the Jackson's inequality [22]: if f : [0,2n] ^ C is a p times differentiable periodic function such that |f(p)(x)| ^ 1, 0 ^ x ^ 2n, then, for every positive integer n, there exists a trigonometric polynomial
Tn-1 of degree at most n — 1, such that |f(x) — Tn-1(x)| ^ C(p~ for any
x G [0, 2n], where C(p) depends only on p. So, the vector norm of Q2 can be estimated by this inequality as K/F2. The second summand also behaves as O(1/F2), due to presence of the coefficient r. So, the error due to the series tail is O(1/F2). □
The matrix of the truncated system consists of the Fourier coefficients of 5i, j ± Ki, j or Li, j.
Each solution of the truncated system has the form of Fourier polynomial and is linear dependent on the real constants Aj, Bl and Dj.
Since p(t)+iq(t) is the boundary value of the analytic in D function, the corresponding Cauchy integral vanishes at the points that do not belong to D. So, we obtain the equations that constitute the linear system over 2n unknown real constants Aj, Bl, and Dj. The constants Aj, Bl together with the free summands Dk of the function Fourier series expansions, are the solutions of the complex system
p(r) + iq(T)
——-d(z(t)) = 0, l = 1,... ,n.
z(t) — z\ y v " ' ' '
dD
Here z*, I = m + 1,... ,n, are the points inside the domains bounded by the contours Cl, respectively. Note that we have n real constants Aj, Bl. Finally, for the MBVP1 the solution equals
f (x,y) = Re
m
» + E A ln(x +" j)
dD
For the MBVP2 the solution has the form
m— 1
f (x,y) = Re
1 t p(T) + iq(T) d(z(T)) + y-; Aj ln(x + iy - z;)-
2m J z(t) - x - iy y y " j
dD j=1
m-1
—in(x+¿y— oX! . j=1
The proof is completed. □
2. Examples. In all of the examples we take F = 200, M = 100.
Figure 1: The contour lines of the solution in Example 1 (left), the contour lines of the solution in Example 2 (center) and the extremal electric field strength at the domain boundary in Example 2 (right).
Example 1. Consider the problem for the symmetric domain with the elliptic boundary components zo(t) = 6ert—e-it, z1(t) = 0.8eit+0.2e-it+3i, z2(t) = 0.8ert + 0.2e-it — 3i. Here we solve the mixed boundary value problem with the following data: the unknown function possesses zero normal derivative at the exterior boundary component and equals 1 and 2 at the inner upper and lower components, respectively. It is the MBVP2.
Figure 1 (left) demostrates the contour lines of the solution in this case.
Example 2. Consider the problem for the asymmetric domain with the elliptic boundary components zo(t) = 6ert — e-it, z1(t) = 0.8ei(i+n/4) + + 0.2e-i(t+n/4) + 3i, z2(t) = 0.8ei4 + 0.2e-it — 3i. Here we solve the mixed boundary value problem with the following data: the function has zero normal derivative at the exterior large elliptic boundary component and equals 1 and 2 at the inner upper and lower elliptic components, respectively. Figure 1 (center) shows the contour lines of the solution.
Figure 1 (right) shows the extremal points of the electric field distribution at the domain boudary. The letter A denotes the points of the minimal electric field strength and the letter B shows the point of the maximal electic strength.
Example 3. Consider the problem for the symmetric domain with the elliptic boundary components z0(t) = 6elt—e~lt, Z\{t) = 0.8eli+0.2e_li+3i, z2(t) = 0.8eu + 0.2e~u - 3i, z3(t) = 0.8eu - 0.3e~u. Here we solve the following boundary value problem: the unknown function has zero normal derivative at the exterior boundary component and equals 1, 1.5, and 2 at the inner upper, middle, and lower components, respectively. Figure 2 demonstrates the solution graph.
Figure 2: Solution graph for Example 3.
Example 4. Consider the problem for the symmetric domain with the elliptic boundary components z0(t) = Qelt — e~lt, Z\{t) = 0.8ert + + 0.26-^ + 3?;, z2{t) = 0.8^ + 0.26-^-3?;, z3(t) = 0.8eu — 0.3e_ii. Here we solve the boundary value problem with the following data: the function possesses zero normal derivative at the exterior and the middle boundary components and equals 1 and 2 at the inner upper and lower components, respectively. Figure 3 (left) demonstrates the solution graph.
The letter A denotes the points of the minimal electric field strength and the letter B shows the point of the maximal electic field strength (Fig. 3, right).
Example 5. Consider the domain with the boundary given by the equations z0(t) = 5.5elt — 0.5e~lt and Zi(t) = 0.5ert — 0.1e_rt (Fig. 4, left).
Figure 3: The solution graph for Example 4 (left) and the extremal electric field strength at the domain boundary in Example 4 (right).
We approximate the exact solution x2 — y2 in this domain. Assume that the function values are known at the exterior boundary component and the normal derivatives are given at the interior boundary component. This is the MBVP1. We have no constants of type A. Difference between the constructed solution and the exact one inside the domain is less than 0.0006 (Fig. 4, right).
Figure 4: The domain in (left) and the difference between the approximate and exact solutions (right) in Example 5.
3. Conclusion. The Cauchy integral method can be applied not
only to Dirichlet problems but also to problems with mixed boundary conditions on multiply connected domains. The method is easily programmable. The solution is a smooth function that allows deeper insights into the solution properties. The process of solution does not require the domain mesh, nor construction and analysis of the weak solution.
References
[1] Shirokova E. A., El-Shenawy A. A Cauchy integral method of the solution of the 2D Dirichlet problem for simply or doubly connected domains Numerical Methods for Partial Differential Equations, vol. 34, 2018, pp. 2267-2278. DOI: https://doi.org/10.1002/num.22290.
[2] Shirokova E. A., Ivanshin P. N. Approximate Conformal Mappings and Elasticity Theory, Journal of Complex Analysis, vol. 2016, Article ID 4367205, pp. 1-8, 2016.
[3] Shirokova E. A. On approximate conformal mapping of the unit disk on an simply connected domain Russian Mathematics (Iz VUZ), 2014, vol. 58, N. 3, pp. 47-56.
[4] Gakhov F. D. Boundary Value Problems. Courier Corporation, 1990.
[5] Jaswon M. Integral equation methods in potential theory, I Proc. Roy. Soc. Ser. A 275, 1963, pp. 23-32.
[6] G. Symm, Integral equation methods in potential theory, II, Proc. Roy. Soc. Ser. A, vol. 275 (1963) 33-46.
[7] Mikhailov S. E. On an integral equation of some boundary value problems for harmonic functions in plane multiply connected domains with nonreg-ular boundary Mat. Sb. (N.S.) vol. 121(163), 1983, pp. 533-544.
[8] Axelsson O., Barker V. A. Finite Element Solution of Boundary Value Problems. Academic Press, Orlando, FL, 1984.
[9] Ciarlet P. G. The Finite Element Method for Elliptic Problems. North-Holland, New York, NY, 1978.
[10] Johnson C. Numerical Solution of Partial Differential Equations by the Finite Element Method. Cambridge University Press, Cambridge, MA, 1987.
[11] Babuska I.,. Tempone R, Zouraris G. E. Galerkin finite element approximations of stochastic elliptic partial differential equations, SIAM J. Numer. Anal. 42 (2) (2004) 800-825.
[12] Cialdea A., Leonessa V., Malaspina A. On the Dirichlet problem for the Stokes system in multiply connected domains Abstr. Appl. Anal. vol. 2013, 2013, pp. 1-8.
[13] Costabel M., Riva M. D., Dauge M., Musolino P. Converging expansions for Lipschitz self similar perforations of a plane sector Integr. Equ. Oper. Theory, vol. 88, 2017, pp. 401-449.
[14] Hsieh C., Kassab A. J. Complex variable boundary element methods for the solution of potential problems in simply and multiply connected domains Comput. Methods Appl. Mech. Eng. vol. 86, 1991, pp. 189-213.
[15] Liu C. A highly accurate collocation Trefftz method for solving the Laplace equation in the doubly connected domains Numer. Methods Partial Differ. Equ. vol. 24, 2008, pp. 179-192.
[16] Anselone P. M., Moore R. H. Approximate solutions of integral and operator equations J. Math. Anal. Appl. vol. 9, 1964, pp. 268-277.
[17] Buckner H. Numerical methods for integral equations in Survey of Numerical Analysis, John Todd, ed., McGraw-Hill, New York, 1962.
[18] Anselone P., Gonzalez-Fernandez J. Uniformly convergent approximate solutions of Fredholm integral equations J. Math. Anal. Appl. vol. 10 (3), 1965, pp. 519-536.
[19] Ray S. S., Sahu P. K. Numerical methods for solving Fredholm integral equations of second kind Abstr. Appl. Anal. vol. 2013, 2013, 17 p.
[20] Asmar N. H., Partial Differential Equations with Fourier Series and Boundary Value Problems, Dover Publications, 2016.
[21] Lyusternik L. A., Sobolev V. I. The elements of functional analysis. Gosu-darstv. Izdat. Tehn.-Teor. Lit., Moscow-Leningrad, 1951.
[22] Achieser N.I. Theory of Approximation. New York: Frederick Ungar Publishing Co. 1956.
Received June 27, 2018.
In revised form,, December 20, 2018.
Accepted December 20, 2018.
Published online January 18, 2019.
Kazan Federal University 18 Kremlyovskaya str., Kazan 420008, Russia E-mail: [email protected], [email protected]