9512.net

# Numerical solution of Fredholm integral equations

1. Introduction

Numerical solution of Fredholm integral equations using simple wavelet expansions
Starting with a wavelet expansion of the unknown function f , we exhibit a scaling function expansion and then report on numerical experiments done with both the Daubechies D6 and D4 scaling functions. We begin by describing the D6 case (the D4 case is similar; see 1] and 2]). We rst introduce the Daubechies D6 scaling function and its corresponding W6 wavelet w. D6 is de ned recursively by a two-scale equation with 6 terms, that is, (x) =
5 X

B. BERTRAM

We consider the numerical solution of Fredholm integral equations of the second kind of the form Zb f (x) = g(x) + k(x; t)f (t)dt; a x b: (1)
a

i=0

ci (2x ? i);

where c0 = 0:44300824740057; c1 = 1:0587400346807; c2 = 0:59544707975904; c3 = ?0:13601649480113; c4 = ?0:28845532715960, and c5 = 7:7276460120481D ? 02. This scaling function and the ensuing wavelet are continuously di erentiable functions. In addition, the coe cients have the following properties: (i) orthogonality: 1 n X 2 if m = 0; ck ck?2m = 0 otherwise. (ii) smoothness:
k=?1

1 X
k=?1

(?1)k km ck = 0; 0 m 2:

The two-scale equation is then used to compute at a family of dyadic points i=2j . The W6 wavelet is de ned to be

w(x) =
We wish to express f as

1 X
i=?1

(?1)i c1?i (2x ? i):

f (x) =

1 1 X X
j =?1 k=?1

Ajk wjk (x);
1

where

jx ? wjk (x) = w((2 ? k) k) if j = 0,1. x if j ?

2. The scaling function expansion
f (x) =

To nd an approximation for the solution f using the D6 wavelets, we use an idea similar to one found in 2]. We assume that a wavelet expansion for f is
1 1 X X
j =?1 k=?1

Ajk w(2j x ? k):

(We remark that this expansion agrees with the previous wavelet expansion for f when the j = ?1 term is appropriately de ned as above.) It can be shown that the partial sums of the above series can be expressed in terms of a single series involving only the scaling function. Lemma 1. If the nth partial sum wavelet series approximation of f is de ned to be

Sn (x) =
then

n 1 X X j =?1 k=?1

Ajk w(2j x ? k);

Sn?1 (x) =

1 X
k=?1

Bnk (2n x ? k):

Proof. The result follows by the use of induction, a shifting of indices, the de nition of wavelet and interchanging the order of summation.
We consider equations of the form (1) where, for testing purposes, we use solution functions on the interval 0,1]. We expand f (as indicated above) in terms of the scaling function (either D4 or D6), that is,

3. The computational technique

fn (x) =

X

j

aj (2n x ? j );

the projection of f on to the subspace Vj of the multiresolution analysis de ned by (see 2]). A close inspection of the argument 2n x ? j and anconsideration of the support of indicate that the bounds for j are ?2 j 2 ? 1 for the D4 case, and ?4 j 2n ? 1 for the D6 case. We now discretise the interval of integration into 2n intervals. If the kernel of the integral equation is continuous, we use Simpson's rule to approximate the integral. If the kernel is weakly singular, we use a linear interpolation of the non-singular part of the kernel times , and 2

use product integration to produce the approximate integral (see 3]). We do this for the required number of j values and produce a system which is then solved for the coe cients aj . Since the number of j values is greater than 2n , we choose the necessary additional points at the next level of dyadic points. It is worth noting that this choice of evaluation nodes coupled with the use of the nth level approximation of the unknown function using the scaling function alone necessitates knowledge of only at integer and half-integer points. This provides a savings in storage needed and computational time used in the construction of the matrix and reconstruction of the solution from the coe cients. We present some of the computational results in a table, listing the kernel, solution function, type of wavelet, level of discretisation and average absolute error. The programs were written in FORTRAN and run using double precision on a Sun Sparc ELC workstation. In addition, routines from LINPACK were used to solve the systems involved. Kernel xt xt (xt)9 (xt)9 jx ? tj?1=2 jx ? tj?1=2 ext ext Solution x x4 1: 1: x(1 ? x) x(1 ? x) x(1 ? x) x(1 ? x) Wavelet D4 D4 D4 D6 D6 D4 D4 D6 Level (n) 6 6 9 9 9 9 9 7 Average Error exact 1.50 10?8 4.33 10?12 4.32 10?12 1.35 10?4 3.02 10?3 7.0 10?8 2.24 10?9

4. Computational results

Table 1. Computational results. The D4 wavelet solution compares favourably with that obtained by means of the Nystrom-Simpson method. Both methods take less than a second and yield identical results for the n = 6 case for the (xt)9 kernel and solution f = 1:0. The results remain identical for the n = 9 case, but the wavelet solution takes about 4 times longer to compute. Computational results also indicate that the D6 wavelet is more computationally unstable than the D4 wavelet, yielding condition numbers for the resultant matrices which are several orders of magnitude larger. For the exponential kernel, the D6 results are comparable to Nystrom-Simpson up through n = 7, where roundo error and possible oversampling begin a ecting the results. The D4 results are less accurate, but consistently decrease at successive levels of discretisation.

5. Error Theorem 1. The scaling function expansion of f 2 L2(R) converges in sup norm

to f pointwise almost everywhere. Proof. From Lemma 1, we see that our representation of the scaling function expansion of f is equivalent (at any given level of discretisation) to the wavelet expansion of f ; hence, Theorem 2.1 in 5] can be applied to obtain the result. 3

The error in the wavelet approximation is O(2?n) for the D4 case and O(2?2n) for the D6 case 5]. Therefore, in the cases considered, the error will be the larger of this error and that from the numerical integration, that is, O(2?2n) for trapezoid integration (linear interpolation of the kernel) and O(2?4n ) for Simpson integration (quadratic interpolation of the kernel). We remark that the numerical results in Section 4, substantiate the statement in 5] that convergence rates of wavelet expansions are sensitive to both the smoothness of the wavelet and that of the function being expanded.

6. Conclusions

We have demonstrated that the expansion of the solution function in terms of scaling functions is a viable alternative to existing methods for solution of second kind integral equations. The use of the wavelet expansion makes the pointwise convergence of the discretised integral operator to the continuous one automatic except on a set of measure zero (and such points are never included in the set of nodes used). In several instances, this method compares not only in accuracy, but also in e ciency. Since the wavelet coe cients are obtained as a by-product, the solution can be obtained anywhere and the coe cients are available for further work. Further research will include (i) the use of norm-reducing matrix transformations to lower the size of the condition numbers involved (see 6] for a similar idea applied to wavelet-like bases), and (ii) an exploration into the use of the order of convergence to classify unknown solution functions (one idea is to use veri cation schemes (see 7]) instead of estimations for such a study). In addition, we conjecture that we have oversampled the function in the D6 case by the time that n > 7. Acknowledgements. The author wishes to thank Professors O.G. Ruehr and Archie G. Gibson for their time, expertise and helpful comments regarding this work, Robin L. Hocken for her help in the author's beginning exploration of wavelets as a computational tool, and Professor Alex Stone of the University of New Mexico for his generous computing support during the author's sabbatical leave.

References
1. I. Daubechies, Ten lectures on wavelets, SIAM, 1992. 2. J. Smith, The design and analysis of parallel algorithms, Oxford University Press, Oxford, 1993. 3. L.M. Delves and J.L. Mohamed, Computational methods for integral equations, Cambridge University Press, Cambridge, 1985. 4. S. Bertoluzza, Some error estimates for wavelet expansion, Math. Models Methods Appl. Sci. 2 (1992), 489{506. 5. S. E. Kelly, M.A. Kon and L.A. Raphael, Pointwise convergence of wavelet expansions, J. Amer. Math. Soc. 30 (1994), 87{94. 6. B. Alpert, G. Beylkin, R. Coifman and V. Rokhlin, Wavelet-like bases for the fast solution of second kind integral equations, SIAM J. Sci. Statist. Comput. 14 (1993), 159{184. 4

7. C. Falco Korn, B. Hornmann and C.P. Ullrich, Veri cation may be better than estimation, SIAM J. Sci. Statist. Comput. 17 (1996), 1013{1017. Department of Mathematical Sciences, Michigan Technological University, Houghton, Michigan 49931, USA

5

the boundary integral equations discrete into algebraic equations, and then using numerical methods, to obtain the solution of the boundary integral equation...

The numerical accuracy of the second-kind Fredholm integral equations are ...One of the most commonly used methods in solving integral equations is the...
Tikhonov吉洪诺夫正则化
SIAM, 9, 387-392 Phillips DL, 1962, A technique for the numerical solution of certain integral equations of the first kind, J Assoc Comput Mach, 9,...

N: Nearly singular approximations of CPV integrals with end- and corner-singularities for the numerical solution of hypersingular boundary integral equations. ...

numerical method solving algebraic equations, thus ...Keywords: Boundary Element Method; Integral Equation...1903 年由 Erik Ivar Fredholm 奠基了积分方程理论...

Fredholm Equations of the First Kind, Pitman, ...19. D. L. Colton, Solution of boundary value ...Volume 4: Integral Equations and Numerical Methods...