On direct numerical treatment of hypersingular integral equations arising in mechanics and acoustics

19 Pages·0.14 MB·English

by  G. Iovane

On direct numerical treatment of hypersingular
integral equations arising in mechanics and
3
0
0
2 acoustics

n
a
J

6 Gerardo Iovane
2

1 D.I.I.M.A., University of Salerno, 84084 Fisciano (SA), Italy
v
4 e-mail: [email protected]
3
0
1
I.K.Lifanov
0
3
0 The Zhukovsky Aeroforce Engineering Academy,
/
h
p Leningradskii prospect 40, Moscow 125190, Russia
-
h
t e-mail: lifanov [email protected]
a
m
: M.A. Sumbatyan
v
i
X
Faculty of Mechanics and Mathematics, Rostov
r
a
State University, Zorge 5, Rostov-on-Don 344090, Russia
e-mail: [email protected]
Abstract
In this paper we present a treatment of hypersingular integral equa-
tions, which have relevant applications in many problems of wave dynam-
ics, elasticity and fluid mechanics with mixed boundary conditions. The
main goal of the present work is the development of an efficient direct nu-
merical collocation method. The paper is completed with two examples
taken from crack theory and acoustics: the study of a single crack in a
1
linear isotropic elastic medium, and diffraction of a plane acoustic wave by
a thin rigid screen.
1 Introduction
Many problems of fluid mechanics, elasticity, and wave dynamics (acoustics) with
mixed boundary conditions can be reduced to hypersingular equations of the
following type
1
1
ϕ(t) +K (x,t) dt = f (x), x < 1 , (1.1)
"(x t)2 0 # | |
Z
−1 −
where K (x,t) is a regular ( in some sense) part of the kernel . Direct numerical
0
treatment of eq.(1.1) is not easy. Many authors applied the approach, which may
be currently considered as classical. If one defines a solution bounded on the
both endpoints x = 1 ( that will be shown below to vanish as x 1 at
± | ± |
q
x 1 ), then one may search it in the form of the series
→ ±
∞
ϕ(x) = √1 x2 ϕ U (x), x < 1 , (1.2)
j j
− | |
j=0
X
where U (x) are the Chebyshev’s polynomials of the second kind [1]. Then, using
j
the relation
1 1
√1 t2 d √1 t2
− U (t)dt = − U (t)dt =
(x t)2 j −dx x t j
−Z1 − −Z1 −
d
= π [T (x)] = π(j +1) U (x)
j+1 j
− dx −
( T (x) is the Chebyshev’s polynomial of the first kind), applying the scalar
j
product of eq.(1.1) with the functions √1 x2 U (x), one can reduce initial
i
−
integral equation (1.1) to the infinite linear algebraic system of the second kind:
π2 ∞
(i+1) ϕ + k ϕ = f , i = 0,1,..., (1.3)
i ij j i
− 2
j=0
X
1
k = √1 x2 √1 t2 U (x) U (t) K (x,t)dxdt . (1.4)
ij i j 0
− −
ZZ−1
2
We have used here the well known orthogonality relation [1]
1
π
√1 x2 U (x) U (x) = δ , i,j = 0,1,...,
i j ij
− 2
Z
−1
where δ is the Dirac delta.
ij
The study of qualitative properties of the infinite system (1.3) is not the aim
of the present work, and we only note that, when the regular kernel K (x,t) is
0
of a convolution type [K (x,t) = K (x t)], this approach is often applied in a
0 0
−
different way.
We rewrite eq. (1.1) in equivalent form
1 ∞
1
K(x t) ϕ(t)dt = f (x) , x < 1 , K(x) = L(s) e−isx ds, (1.5)
− | | 2π
Z Z
−1 −∞
where L(s) is the Fourier transform of K(x). Let L(s) be integrable for finite s,
and for s it has the following asymptotics
→ ∞
1
L(s) = A s +L (s) , L (s) = O , (δ > 0) , (1.6)
| | 0 0 s1+δ
(cid:18) (cid:19)
then
2A 1 ∞
K(x) = +K (x), K (x) = L (s)e−isxds ,
−x2 0 0 2π 0
Z−∞
where K (x) is at least continuous, and the generalized value of the integral [2]
0
has been used, that is,
∞ ∞
2
s e−isxds = 2 lim e−ǫsscos(sx)ds = . (1.7)
| | ǫ→+0 − x2
Z Z
−∞ 0
Equation (1.5) together with (1.6) is evidently equivalent to eq.(1.1). If one
substitutes the series expansion (1.2) into eqs.(1.5) and then applies a scalar
product with functions √1 x2 U (x), one finally arrives at the infinite linear
i
−
algebraic system [1]
1
∞
l ϕ = f , (i = 0,1,...), l = (1 x2)(1 t2)U (x)U (t)dxdt
ij j i ij i j
− −
jX=0 Z−Z1 q
(1.8)
∞ ∞
J (s)J (s)L(s)ds
L(s)e−is(x−t)ds i+1 j+1 .
× ∼ s2
Z Z
−∞ −∞
3
equivalent to system (1.3). Note that thelast integralis finite foralli,j = 0,1,...,
since J (s) O(sn), s 0; J (s) O(1/√s), s , which, on taking
n n
∼ → ∼ → ∞
into account the asymptotic property (1.6), is evident.
The integral equation (1.1)with a hypersingular kernel canthus be reduced to
an infinite system of linear algebraic equations. This implies numerical treatment
ofdoubleintegralsin(1.4)orsingleintegralsoverinfiniteintervalofstronglyoscil-
lating functions, like in (1.8). When solving integral equations with more regular
kernels, a powerful direct collocation technique, in the framwork of the Boundary
Element Method [3], can successfully be applied. For the Cauchy-type integral
equations, a direct collocation method has been developed by Belotserkovsky and
Lifanov [4, 5]. Hence, the main goal of the present paper is to develop an efficient
direct numerical collocation method.
2 Some properties of hypersingular integrals
First,oneshouldclarifyinwhichsensehypersingular integralsasgivenbyeq.(1.1)
may be treated, since they do not exist either as improper integrals of the first
kind or as Cauchy-type singular integrals. At least three different definitions of
hypersingular integrals are known (see [6]):
1. The integral is a derivative of the Cauchy principal value
b b
ϕ(t) d ϕ(t)
dt = dt. (2.1)
(x t)2 −dx x t
Za − Zα −
2. The integral is treated as a Hadamar principal value [4,5]
b x−ε b
ϕ(t) ϕ(t)dt 2ϕ(x)
dt = lim + . (2.2)
(x t)2 ε→+0  (x t)2 − ε 
Zα − Zα x+Zε −
  
3. The integral is a residue value, in the sense of generalized functions [2],
that is is an analytical continuation of the integral
b
x t αϕ(t)dt (2.3)
| − |
Zα
4
from domain, where it exists in the classical sense, to the value α = 2 .
−
When ϕ(x) 1, the three different approaches give the same result:
≡
b b
dt d dt d b x a b
= = ln − = − .
(x t)2 − dx x t dx x a! (x a)(b x)
Za − Zα − − − −
b
dt 1 1 1 1 2 a b
= lim + + + = − .
(x t)2 ε→+0 − x a ε ε x b − ε (x a)(b x)
Za − (cid:20)(cid:18) − − (cid:19) (cid:21) − −
b x b α+1 α+1
(x a) (b x)
x t αdt = (x t)αdt+ (t x)αdt = − + − , (2.4)
| − | − − α+1 α+1
Za Za Zx
when Re(α) > 1. The analytical continuation of (2.4) from the half-plane
−
Re(α) > 1 to the value α = 2 gives
− −
b dt (x a)α+1 (b x)α+1 a b
= lim − + − = − .
(x t)2 α→−2" α+1 α+1 # (x a)(b x)
Za − − −
All three definitions are equivalent to each other, if the density ϕ(x) is ana-
lytical (i.e. differentiable) on the open interval (a,b). Really, let ϕ(x) be analytic
and Re(α) > 1, then
b x−ε b x−ε b
d ϕ(t)dt d ϕ(t)dt ϕ(t)dt
= lim + = lim +
−dx x t −ε→+0 dx   x t ε→+0  (x t)2−
Za − Za x+Zε − Za x+Zε −
   
b x−ε b
2ϕ(x)
, x t αϕ(t)dt = lim + x t αϕ(t)dt =
− ε # | − | ε→+0 | − |
Za Za x+Zε
 
x−ε b
= lim (x t)αϕ(t)dt+ (t x)αϕ(t)dt =
ε→+0 − − 
Za xZ+ε
 
d x−ε(x t)α+1 b (t x)α+1
= lim − ϕ(t)dt+ − ϕ(t)dt ,
dx ε→+0− α+1 α+1 
Za xZ+ε
 
so, by applying analytical continuation to the last relation, one can see (with
the use of a standard (ε,δ) formalism) that the right-hand side results in (2.1).
Therefore, equivalence of 1 to 2 and 3 is evident. Let us prove that if ϕ(x)
∈
C (a,b), then a finite value of the limit in expression (2.2) exists, and so for
2
5
x (a,b) the integral bϕ(t)/(x t)2dt is finite in any sense. Really, expression
∈ a −
in the square bracketsRin eq.(2.2) is
x−ε b
ϕ(t) ϕ(x) ϕ′(x)(t x) a b b x
+ − − − dt+ϕ(x) − +ϕ′(x)ln − ,
  (x t)2 (x a)(b x) x a
Za x+Zε − − − −
 
which has a finite limit at ε +0.
→
3 Integral equation with characteristic hyper-
singular kernel
Consider hypersingular equation with the characteristic kernel
b
g(t)dt
= f′(x), x (a,b), f (x) C (a,b) . (3.1)
(x t)2 ∈ ∈ 2
Za −
Let us prove that a bounded solution of Eq.(3.1) is unique and is given as follows
b
(x a)(b x) f (t)dt
g(x) = − − .
q π2 (t a)(b t) (x t)
Za − − −
q
Really, according to the definition (2.1), equation (3.1) is equivalent to
b b
d g(t)dt g(t)dt
= f′(x), = f (x)+C, x (a,b) ,
dx (x t) − ∼ x t − ∈
Za − Za −
where C is an arbitrary constant. Now an inversion formula for the Cauchy
characteristic integral operator [7] determines the bounded solution as
b
(x a)(b x) f (t)
g(x) = − − dt, x (a,b) , (3.2)
q π2 (t a)(b t) (x t) ∈
Za − − −
q
and the constant C as
b
1 f (t)
C = dt .
π2 (t a)(b t)
Za − −
q
Thus, as indicated in the Introduction, any bounded solution of Eq.(3.1) vanishes
at x a,b.
→
6
Toconstruct adirect collocationtechnique tosolveeq.(3.1) forarbitraryright-
handside, wedividetheinterval(a,b)intonsmallequalsubintervalsofthelength
h = (b a)/n, by the nodes a = t , t , t ,...,t t = b , t = a + jh, j =
0 1 2 n−1, n j
−
0,1...,n. The central points of each sub-interval (t ,t ) are denoted by x , thus
i−1 i i
x = a+(i 1/2)h, i = 1,...,n.
i
−
If we try to arrange an approximation of integral in (3.1) by using finite sum,
like in regular cases, then for x = x we have
i
b g(t)dt n tj dt h/2 dt
g(t ) = g(t ) +
Za (xi −t)2 ≈ jX=0 j tjZ−1 (xi −t)2 i Z−h/2 t2 (3.3)
1 1 n 1 1
+ g(t ) = g(t ) ,
j j
Xj6=i xi −tj − xi −tj−1! jX=1 xi −tj − xi −tj−1!
where a value of the above hypersingular integral of the function 1/t2 has been
used. So, we try to approximate eq.(3.1) by the linear algebraic system
n 1 1
g(t ) = f′(x ), i = 1,...,n. (3.4)
j i
jX=1 xi −tj − xi −tj−1!
Further considerations are similar to those in [4,5]. Let us prove that, by
assuming x = x (a,b) is fixed, the difference between solution g(x ) of the
l l
∈
system (3.4) and analytical solution (3.2) tends to zero, when n . Obviously,
→ ∞
the system (3.4) can be rewritten as
n 1 1
g(t ) = f′(x ), i = 1,...,n , (3.5)
j i
jX=1 xj −ti − xj −ti−1!
and its principal determnant is
1 1 1 1
...
x t − x t x t − x t
(cid:12) (cid:18) 1 − 1 1 − 0(cid:19) (cid:18) n − 1 n − 0(cid:19) (cid:12)
(cid:12) (cid:12)
∆ = (cid:12) ... ... ... (cid:12). (3.6)
(cid:12) (cid:12)
(cid:12) 1 1 1 1 (cid:12)
(cid:12) (cid:12)
(cid:12) ... (cid:12)
(cid:12)(cid:12) x1 −tn − x1 −tn−1! xn −tn − xn −tn−1! (cid:12)(cid:12)
(cid:12) (cid:12)
(cid:12) (cid:12)
Let us re(cid:12)write the i th row (i = 2,...,n) of ∆ as a sum of al(cid:12)l other rows,
−
7
from k = 1 to k = i:
t t t t
1 0 1 0
− ... − n
(x t )(x t ) (x t )(x t ) (t t )
(cid:12) 1 − 1 1 − 0 n − 1 n − 0 (cid:12) j − 0
(cid:12) (cid:12) j=1
∆ = (cid:12) ... ... ... (cid:12) =
(cid:12) (cid:12) Qn ×
(cid:12) t t t t (cid:12) (x t )
(cid:12)(cid:12) n − 0 ... n − 0 (cid:12)(cid:12) i=1 i − 0
(cid:12)(cid:12)(cid:12) (x1 −tn)(x1 −t0) (xn −tn)(xn −t0) (cid:12)(cid:12)(cid:12) Q (3.7)
(cid:12) (cid:12)
(cid:12) 1 1 (cid:12)
...
x t x t (t t ) (t t ) (x x )
(cid:12)(cid:12) 1 − 1 n − 1 (cid:12)(cid:12) j j − 0 q<p q − p q<p p − q
(cid:12) ... ... ... (cid:12) = ,
×(cid:12) (cid:12) Q(x t ) QQ (xQQt )
(cid:12)(cid:12)(cid:12) 1 ... 1 (cid:12)(cid:12)(cid:12) i i − 0 q, p q − p
(cid:12) x t x t (cid:12) Q QQ
(cid:12) 1 − n n − n (cid:12)
(cid:12) (cid:12)
where(cid:12) the lower limit in all p(cid:12)roducts is 1 and the upper is n. A known value of
(cid:12) (cid:12)
the last determinant has been taken from [4,5].
According to the Cramer’s rule, one needs to calculate the determinant ∆
l
where a column with the right-hand side (3.5) is substituted instead of the l th
−
column of ∆ (3.6):
1 1 1 1
... f′(x ) ...
1
x t − x t x t − x t
(cid:12)(cid:12) (cid:18) 1 −1 1 1 −1 0(cid:19) (cid:18) n −1 1 n −1 0(cid:19) (cid:12)(cid:12)
(cid:12) ... f′(x ) ... (cid:12)
(cid:12) 2 (cid:12)
∆l = (cid:12)(cid:12) (cid:18)x1 −t2 − x1 −t1(cid:19) (cid:18)xn −t2 − xn −t1(cid:19) (cid:12)(cid:12) =
(cid:12) ... ... ... ... ... (cid:12)
(cid:12) (cid:12)
(cid:12) (cid:12)
(cid:12) 1 1 1 1 (cid:12)
(cid:12) ... f′(x ) ... (cid:12)
(cid:12) n (cid:12)
(cid:12) x1 tn − x1 tn−1! xn tn − xn tn−1! (cid:12)
(cid:12) − − − − (cid:12)
(cid:12) (cid:12)
(cid:12) t t t t (cid:12)
(cid:12) 1 − 0 ... f′(x ) ... 1 − 0 (cid:12)
1
(x t )(x t ) (x t )(x t )
(cid:12)(cid:12) 1 −t1 t1 − 0 n −t1 tn − 0 (cid:12)(cid:12)
(cid:12) 2 − 0 ... f′(x )+f′(x ) ... 2 − 0 (cid:12)
(cid:12) (x t )(x t ) 1 2 (x t )(x t ) (cid:12)
= (cid:12) 1 2 1 0 n 2 n 0 (cid:12),
(cid:12) − − − − (cid:12)
(cid:12) ... ... ... ... ... (cid:12)
(cid:12) (cid:12)
(cid:12) (cid:12)
(cid:12) t t n t t (cid:12)
(cid:12) n − 0 ... f′(x ) ... n − 0 (cid:12)
(cid:12) k (cid:12)
(cid:12)(cid:12) (x1 −tn)(x1 −t0) kX=1 (xn −tn)(xn −t0) (cid:12)(cid:12)
(cid:12) (cid:12)
where w(cid:12)e have used the same summation of rows as in the case of ∆. (cid:12)
(cid:12) (cid:12)
The last determinant may be calculated arranging expansion by elements of
the l th column as follows (see also [4,5])
−
m
(t t ) f′(x ) (t t ) (x x )
∆ = j j − 0 n k=1 k ( 1)m+l q<p; q,p6=m q − p q<p; q,p6=l p − q ,
l Q(x t ) Pt t − Q Q (xQ tQ)
i6=l i − 0 mX=1 m − 0 q6=lp6=m q − p
Q Q Q
8
so
m
f′(x ) (x t ) (x t )
g(x ) = ∆l = (x t ) n k=1 k p l − p q q − m . (3.8)
l ∆ l − 0 (t Pt )(x tQ) (t Qt ) (x x )
mX=1 m − 0 l − m p6=m m − p q6=l q − l
Q Q
The last expression at h 0 n , under the condition x (a,b) is
l
→ ∼ → ∞ ∈
fixed, can be estimated as follows:
n l−1 n
(x t ) (x t ) (x t )
l p l p l p
− − −
p=1 p=1 p=l+1
= (x t ) =
Q(x x ) l − l l−Q1 × Qn
qQ6=l q − l q=1(xq −xl) g=l+1(xq −xl) (3.9)
Q Q
( 1)nh l−1 1 l−1 1 h√b x
= − 1 . 1+ ( 1)n − l , n ,
2 pY=1 − 2p! pY=1 2p! ∼ − π√xl −a → ∞
if x belongs to the open interval (a,b) (i.e. l = n, l = 1). Here we have used
6 6
the asymptotic estimate [4,5]
n β nβ
1+ = +O nβ−1 , n ,
m! Γ(1+β) → ∞
mY=1 (cid:16) (cid:17)
Further, by analogy
(x t )
q q − m = ( 1)n(t x )m−1 1+ 1 n−m 1 1 ( 1)n h√tm −a.
pQ6=m(tm −tp) − m − m qY=1 2q! qY=1 − 2q! ∼ − π√b−tm
Q
Other terms in eq. (3.8) can be simplified as follows (h 0):
→
m tm
(x t ) (x a), (t t ) (t a), h f′(t ) f′(t)dt = f(t ),
l 0 l m 0 m k m
− → − − → − →
kX=1 Za
hence expression (3.8) with h 0 tends to
→
h n f (t )
m
g(x ) (x a)(b x )
l ∼ π2 l − − l (t )(b t ) (x t ) ∼
q mX=1 m−a − m l− m
q
b
(xl a)(b xl) f (t)dt
− − .
∼ q π2 (t a)(b t) (x t)
Za − − l −
q
An example on application of the proposed numerical method is shown in
Fig.1, where (a,b) = ( 1,1).
−
9
4 Full hypersingular equation
Consider the full equation
b
1
+K (x,t) g(t)dt = f′(x), x (a,b) , (4.1a)
"(x t)2 0 # ∈
Za −
where
∂K (x,t)
1
K (x,t) = . (4.1b)
0
∂x
Itsboundedsolutioncanbeconstructedbyapplyinginversionofthecharacteristic
part, that reduces eq.(4.1a) to a second-kind Fredholm integral equation
b
g(x)+ N (x,t)g(t)dt = f (x), x (a,b) , (4.2)
1 1
∈
Za
where
b
(x a)(b x) K (τ,t)dτ
N (x,t) = − − 1 ,
1 q π2
(τ a)(b τ) (x τ)
Za − − −
q
b
(x a)(b x) f (τ)dτ
f (x) = − − .
1 q π2 (τ a)(b τ) (x τ)
Za − − −
q
Let us prove that, if f (x) C (a,b); K (x,t) C [(a,b) (a,b)], then for
1 1 1
∈ ∈ ×
any x (a,b) the difference between solution g(x) of the linear algebraic system
∈
n 1 1
+hK (x ,t ) g(t ) = f′(x ), i = 1,...,n
0 i j j i
jX=1"xi −tj − xi −tj−1 #
and the bounded solution of eq.(4.1) tends to zero when h 0 (i.e. n ).
→ → ∞
Really, if one transfers the terms, related to the regular kernel, to the right-hand
side and solves so written linear algebraic system with the characteristic matrix
1/(x t ) 1/(x t ), then one arrives at a finite-difference approximation
i j i j−1
− − −
of the eq.(4.2). The proof is finally completed, if one applies classical results on
numerical solution of the second-kind Fredholm integral equation.
An example for K (x,t) = A(x t), f′(x) = π f (x) = πx is
0
− − ∼ −
shown in Fig.2, where numerical solution is compared with the exact one g(x) =
8√1 x2(4+Ax)/(32+A2) in the case A = 3, (a,b) = ( 1,1).
− −
10