V. I. Krylov N.S.Skoblya

A Handbook

of Methods

of Approximate Fourier Transformation and

Inversion

of the Laplace Transformation

Mir Publishers Moscow

B. MW. KPDIJIOR, H. C. CKOBJIA

METOJIDI NPHBJIMJKEHHOTO NPEOBPASOBAHHA DYPbE MW OBPAWEHHA TIPEOBPASOBAHHA JIATIJIACA

CIIPABOUHAH KHYTA

HanatenbctBo «Hayka» Mocxksa

V. I. Krylov and N. S. Skoblya

A HANDBOOK

OF METHODS

OF APPROXIMATE FOURIER TRANSFORMATION AND

INVERSION

OF THE LAPLACE TRANSFORMATION

Translated from the Russian by George Yankousky

MIR PUBLISHERS MOSCOW

First published 1977 Revised from the 1974 Russian edition

Ha antt2auticrom a3zvuKe

Contents

List of symbols ........ Be ie eh igh > Ge at o> ae ee 8 PIOlace: 2 oa: He ee ES Awe 9 Part One

INVERSION OF THE LAPLACE TRANSFORMATION

Chapter 1. Introduction 1.1 Basic concepts in the theory of the Laplace transforma-

DIODD Be hee oo oy eee A Ne a ote, lp AN Ta st BeBe ap! Ho Poe 15 1.2 Complex integral for computing inverse Laplace trans- POV: «5 a Si ae a ee ask Seg, oaks oe, RE oy Sek, ae SE Be 22

1.3 Representing functions by the Laplace integral . . . . 29 1.4 Ill-conditioning of the problem of computing inverse Laplace transforms ..........040+80808+8 30

Chapter 2. Some Analytical Methods for Computing Inverse Laplace Transforms

2.1 Finding the original function via the inversion formula 32 2.2 Expanding the original function into power series . . 36 2.3 Expanding the original function into generalized power

SORICS!: °° 54.. sic 4o ha are St Ge a OR we, Goa: SS cai 38

Chapter 3. Methods of Numerical Inversion of Laplace Transforms Based on the Use of Special Expansions

3.4 Computing inverse Laplace transforms by polynomials or-

thogonal on a finite interval ........2.2.2.2.. 41 3.2 Computing inverse Laplace transforms with the aid of the Fourier sine series ...........e448-4 67

3.3 Computing inverse Laplace transforms with the aid of series in terms of generalized Chebyshev-Laguerre po- lynomials .......2....0.8888888 8 70

Chapter 4. Methods of Computing the Mellin Integral with the Aid of Interpolation Quadrature Formulas

.1 The general theory of interpolation methods ..... .2 The equal-interval interpolation method ...... -3 The unequal-interval interpolation method .... . -4 Other interpolation methods. Using the truncated Tay- TOR BORICS «se Be ia cack ye IP oe Bek itn oe lee a ce -0 Some theorems on convergence of dekege opto te -6 Theorems on the convergence of interpolation methods Ol-INversion, * a e BO r S

Chapter 5. Methods of Numerical Inversion of Laplace Trans- forms via Quadrature Formulas of Highest Accuracy

5.4 The theory of quadratureformulas .......... 5.2 Orthogonal polynomials connected with the quadrature formula of highest accuracy ........... 5.3 Methods for computing the coefficients and points of a quadrature formula .........2.2.2+.2.2+86-

Chapter 6. Methods of Inverting Laplace Transforms via Quad- rature Formulas with Equal Coefficients

6.4 Constructing a computation formula ........ 6.2 Remark on the spacing of points .........

Part Two

FOURIER TRANSFORMS AND THEIR APPLICATION TO INVERSION OF LAPLACE TRANSFORMS

Chapter 7. Introduction

7.1 Fourier transforms .......4...+.+.4.6442-.4 7.2 Reducing integrals of the Mellin type to the Fourier transformation .......2..2.2+.+2+-++.44+4

Chapter 8. Inversion of Laplace Transforms by Means of the Fourier Series

8.1 The case of a rapidly decreasing original function f (z) 8.2 The case of rapid decrease of the modulus of the image TUNCUION, Cp) e0usse 4 Se a ech Ee ae le

118

145 148

150 157

159

CGhaptor 9. pe Formulas for Computing Fourier Integ- rals

9.1 Some preliminary POMALKS:. 6's ey bs se Ge: we Rk 9.2 Algebraic interpolation of the function f (z) 9.3. Interpolation by rational functions

Chapter 10. Highest-Accuracy Formulas for Computation

10.4 Introduction ..........2.20 2480088 6 10.2 Constructing a formula of highest degree of accuracy .

Part Three ~

ISOLATING SINGULARITIES OF A FUNCTION IN COMPUTATIONS

Chapter 11. Isolating Singularities of the Image Function F (p)

40.4 Introduction: se oe 86 cae ek a ae ve ee

11.2 Removing and weakening the singularities of the image MINCtION.F3(p)- 6 -i8o Gk ae. Se ce aa ee eS te

11.3 A remark on the increase in the rate of approach to zero of the image function F (p) ........2.2..-.

11.4 A table of image functions F (p) and the corresponding original functions f (x) for constructing the singular part of the image function IF, (p)

Chapter 12. Isolating Singularities of a Function in the Fourier Transformation

12.4 Removing discontinuities of the first kind ..... 12.2 Increasing the rate of approach to zero of the function undergoing transformation

e e e e e e e e e e ° e e e e e e s e e s e e °

164 166 201

233 236

245 247 252

204

299 263

266 268

f (t) F (p)

@ (p)

Pc (Pp)

Ps (p) Pr'(a) Pr (2) P7 (2) Tn (2)

Uz, (2)

& res f (a) 8

List of Symbols

original function

image function under the Laplace transfor- mation

complex Fourier transform of a function f Fourier cosine transform of a function f Fourier sine transform of a function f Jacobi polynomial of degree n of the para- meters a, B

ar Jacobi polynomial for the interval 10, 4

shifted Legendre polynomial for the inter- val [0, 4]

shifted Chebyshev polynomial of the first kind for the interval [0, 1]

shifted Chebyshev polynomial of the second kind for the interval [0, 4]

generalized Chebyshev-Laguerre polynomial of degree n

confluent hypergeometric function

the sum of all sets M; with respect to &

the intersection (common portion) of all sets M;

total variation of a function f on the inter- val (a, b) (the interval need not be indicated) = — lim g (f + 8) (6>0, 6—>0)

= lim g (t — 6) (6>0, 6 0)

residue of the function f (z) at the polez = a

Preface

The problem of approximate computation of inverse l.aplace transforms and, in particular, numerical inversion arose out of the need to obtain a numerical solution where existing tables of functions and their transforms (images, or image functions) do not enable one to obtain the original function from the image function or require excessive com- putation.

Many methods of inversion have been constructed over the past two decades but they are still scattered among specialized journals and books and, as a rule, are not gene- rally known. As far as the authors know, there are no books in the scientific literature that contain a systematic pre- sentation of all these methods.?

The authors’ purpose here has been to write a book with a full description of the present state of the inversion problem and to make the text useful for problem solving with the help of the Laplace transformation. The first aim was rela- tively easy to attain since the literature on the inversion problem is still slight and comparatively easy to survey.

1 In a way, reference [13] is an exception, but its purpose was morely to unite all known auxiliary numerical tables connected with computing the Mellin integral. The mathematical theory of inver- sion of Laplace transforms is given only briefly, in the form of ex- planations to the tables.

As to how useful this text will actually be, only time will show, but the authors believe that a few remarks are in order. Most important undoubtedly is the remark that the problem has not been investigated with sufficient thorough- ness and the results obtained are still small in number. For those readers who are not acquainted with the problem of inverting Laplace transforms, a few explanatory remarks should be made. What follows is a nonrigorous (but picto- rial) presentation.

The inversion problem is merely the problem of finding the solution f (z) of an integral equation of the first kind:

\ f (2) e-P* da = F (p) (+) 0

where F (p) is taken to be a known function of the complex argument p, the function being analytic in a certain half- plane of the form Re p > y (y < ov).

The equation (*) is the Laplace transformation of the function f (x) into F (p). The kernel of the integral, e-?*, is an entire analytic function of the arguments z and p with smooth variation, and the operation of integration which performs the averaging of f with weight e-?* can appreciably smooth the peculiarities in the behaviour of the function f being transformed.

In the inversion problem, one uses the image function F (p), which, due to its analyticity, varies very smoothly, to restore all possible irregularities in the behaviour of the original function f (x). It is therefore to be expected that the apparatus used to solve the inversion problem cannot be simple and rough but must be complicated and sensitive

10

even to small nuances in the behaviour of F (p) in order to perform the delicate work of recovering the original func- tion f (z).

Note yet another property of the inversion problem—the instability of the original function f relative to small varia- tions in the image F. This instability is evident from a glance at the transformation (+). Indeed, suppose f is the given original function and F is the corresponding image function. Subject f to any change, even a strong one, but over a very small interval. Call the new original function f,. Such a change will hardly affect the integral and will alter the image function F (p) only slightly, so that the new image function F, (p) will be close to F (p).

We could have constructed a new original function f, by changing f on many intervals (instead of one) of the half axis 0 < x < oo, but in a way that the total sum of the lengths of these intervals is a small quantity, for it will still be possible to make the new image function F, (p) as close to F (p) as desired so that to the close-lying image functions F (p) and F, (p) there will correspond the original functions f (z) and f, (z) that differ radically over many intervals of the half axisO < zt < oo.

We can now give a more complete explanation as to why and in what respect the content of this text is to be regarded as insufficient and we can attempt to indicate the reasons why we believe that these defects will hardly be removed in the near future.

At the present time we stand at the beginning of the construction of a theory of approximate inversion of the Laplace transform. In most published papers we find either

11

modifications of familiar methods of computation or newly devised computational schemes that have no counterparts in the past. Instances of the former kind are methods of computing the complex Mellin integral based on the idea of interpolation or the idea of Gauss of attaining the highest possible degree of accuracy, and others. Now the use of orthogonal polynomials for inversion appeared only in con- nection with the Laplace transformation and did not have closely related analogs in the past.

The construction of a computational method or any indica- tion of the possibility of construction is only the first step towards constructing a theory of the method. The steps that follow include determining the conditions of convergence of the rule and the rate of convergence, finding estimates of a priori and a posteriori errors, elaborating ways of improving convergence if it is not sufficiently fast, and so on. Such investigations of the inversion problem are particularly important due to its instability.

As the reader will see from the contents of this book, few results of this kind have been obtained, and then only in the simplest cases.

The main reasons that hamper a rapid solution of the problems at hand are those two properties of the inversion problem that were mentioned above: the unavoidable com- plexity of the mathematical apparatus involved in inverting Laplace transforms and the instability of the inversion problem relative to variations in the functions F, which instability must give rise to some form of instability in any computational process of the solution of this problem.

The fact that the rules for approximate inversion of

12

transforms have not been studied with sufficient thorough- ness makes this handbook incomplete in many respects. It is hard to expect a rapid improvement in our knowledge in these questions in the immediate future, and so for some years to come this book will most likely remain incomp- lete.

So far we have spoken only of inverting Laplace transforms and have not mentioned harmonic analysis at all. This is because in this book harmonic analysis is touched on only lo the extent that it is connected with the inversion problem; namely, we have taken from harmonic analysis only such problems involving computation of Fourier integrals as are closely related to inverting Laplace transforms; and comput- ing Fourier integrals is for us one of several possible ways of solving the inversion problem. For this reason, all that has been said about inversion refers also to computing Fourier integrals.

This text is aimed at a broad category of readers engaged in the theory of the Laplace transformation or its scientific and technical applications. For this reason, the authors did not strive for particular brevity in presentation and attempt- ed to make the text accessible to nonmathematicians as well. It is assumed the reader has a basic knowledge of analysis and the theory of functions of a complex variable as given in any extended college course of mathematics.

The authors

13

Part One

INVERSION OF THE LAPLACE TRANSFORMATION

Chapter 1

Introduction

1.4 Basic concepts in the theory of the Laplace transformation

The past few decades have seen particularly frequent and successful use of operational methods on the basis of the Laplace transformation applied in mathematics, mechanics and engineering. These methods have found extensive appli- cations in the theory of thermal conductivity, in electric and radio engineering, in the study of nonstationary pheno- mena in electric networks, in problems of the dynamics of nutomatic control systems, in the theory of linear differen- tial, integral, and difference equations, and in many other problems.

The operational method of problem solving may be divided into four stages:

(1) one passes from the desired original function f (¢) to the image function F (p);

19

(2) operations are performed on F (p) that correspond to operations on f (t); an equation in F (p) is then obtained which frequently turns out to be much simpler than the equation in the original function; for example, an ordinary differential equation is replaced by an algebraic equation, a partial differential equation is replaced by an ordinary differential equation, and so forth;

(3) the image equation thus obtained is solved for F (p);

(4) from the image function F (p) thus found one passes to the original function f(t), which is the desired func- tion.

In many cases, the most difficult stage is the fourth, i.e. finding the original function f (t) from the image func- tion F (p), or computing the inverse Laplace transform. Extensive tables of correspondences between original func- tions and image functions are available so that one can find the original function from a given image function. However, these tables do not by far encompass all cases that may be encountered in practice and, what is more, it often happens that the original function is expressed in terms of very complicated functions that are hard to compute and are not always tabulated. Then it is either impossible or unde@ sirable to find the exact original function. This necessitates the construction of approximate methods for computing inverse Laplace transforms that permit finding the original function in a broad class of cases. These methods are dis- cussed in Part One of this text.

Let us now recall certain familiar facts concerning the theory of the Laplace transformation.

Given on the half axis 0 << t< oo a function f (t) that is integrable! with its absolute value on any finite interval la, b] O<ax<b< ow).

1 Here and henceforth we assume the Riemann integral.

16

Let us introduce the complex parameter p = o ++ it und define the Laplace transform of the function f by

F (p) = \ f(t) e-P¢ dt (1.1.1) 0

Ilere, the value of the improper integral over the half axis [(0, oo) will be taken to mean the limit to which the integral over the finite interval [0, B] tends as B-»oo, so that

\ f(e-Ptdt=lim \ f (t) e-Pt dt (1.1.2) 0

B->0o

We say that-the Laplace transformation is applicable to the function f for the value p of the parameter if for that value of p the integral (1.1.1) converges.

It can be verified that if the transformation (1.1.1) is applicable to f for p = po = oo + itp, then it is applicable to f for any value p = o + it for which Re (p — py) = == 6 — 6), > 0. Indeed, consider the function @q (t) =

t

=| f (u) e-?o” du. If the integral (1.1.1) converges for 0 P = Po, then q (t) has a finite limit lim @g (¢) and conse- t-»0o quently is bounded on the half axis 0 < t < o: le |<@<©

To prove the convergence of the integral (1.1.1), we take advantage of the Bolzano-Cauchy criterion. Take three arbitrary positive numbers a, b, c and assume Re (p — py) = = 6—o,) >c. In the computations given below we made

2—0281 47

use of integration by parts:

a+b a+b | | f(penrtae]=| | e-@-rorde

a a- a

a+b =|0() en (PP | "+ (P—po) \ @(t)e-@-P dt | a

<e-.20+| p-- pol @ [ee dt =e-% (9 +1e—hl) ¢Q

The last term in this chain of relations does not depend on b and for any fixed p (as a increases) becomes less than any preassigned positive number. Therefore the Bolzano-Cauchy criterion holds true for the integral under the limit sign in (1.1.2), and the integral (1.1.1) converges.

Moreover, let p vary in a bounded closed region D lying inside the half-plane Re p > op. Clearly, for p € D there exist numbers c and M such that the following inequalities hold true: Re (p — po) = 6 — op > ¢ and | p — pp | < M. From the inequalities thus obtained it follows that the Bolzano-Cauchy criterion will hold uniformly with re$pect

b

to the parameter p in D. Since the integral \ f (t) e-?! dt

0 is Clearly an entire analytic function of p, it follows, from its uniform convergence in D to the integral (1.1.1), that F (p) is an analytic function regular in D and since D is any interior region of the half-plane Re p > op, the function F (p) is regular throughout this half-plane.

We now consider the set £ of all real values p = o of the parameter p for which the Laplace transform (1.1.1) is applicable to the function /, and denote by y the lower bound (the greatest lower bound of the values o) of this set:

y = info E

The value of y is determined on the basis of the following facts.

1. When y has a finite value, we can say that the Laplace transform (1.1.1) is applicable to f throughout the open half- plane Re p > y, and F (p) is a regular function of p at least in this half-plane and will not be applicable to f for a single value of p in the half-plane Re p < y.

2. When y = —oo, the Laplace transform (1.1.1) is upplicable to f for any value of p and F (p) is a regular function over the entire complex plane p.

3. When y = -++oo, the Laplace transform (1.1.1) is not upplicable to f for any value of p.

The number y may be called the boundary of the index of convergence, and the straight line Rep =o = y the boundary of the region of convergence of the Laplace trans- form.

In this connection, one ordinarily refines the concept of the original function. The function f is called the original function if it has the following properties:

(1) f is defined on the axis —oco < t < oo and is integrable with absolute value on every finite interval;

(2) the function f vanishes for f < 0;

(3) the Laplace transform is applicable to f for at least one value of p.

The following theorem holds true.

Theorem 1. To every original function { there corresponds a number y (—co < y < c©) such that the transform (1.1.1) is applicable to f for any p, where Re p = o > jy; here, F (p) is a regular function in the half-plane Rep =a >y. The transform (1.1.1) is not applicable to f for any value of p for which Rep < y.

The function F is called the Laplace image function of f. It is sometimes difficult to find the value of y. There is simple procedure for estimating the upper bound of y in many cases.

2 19

Theorem 2. If there are two numbers M (0 <M < oo) and @ (—co <a < oo) such that for any t > 0 the inequality

if (t) | < Mex (4.4.3)

holds true, then a > y. Proof. Let Rep =o>a

\If @ em | dt = (lA @ lene dt = | Me~@- dt" < © 0 0 0

Thus, for any real value p = o exceeding a, the integral (1.4.1) will be absolutely convergent, and since y is the lower bound of such values of p, it must be true that y < a, which completes the proof.

The condition (1.1.3) holds for a broad class of functions, in particular for-most functions encountered in applications; for this reason, the original functions f are frequently defined somewhat differently: namely, the first two properties are retained, while the third is replaced in the following fashion.

There exist two numbers M and a such that the inequality

lf () |< Me’, 0<t<oaw

is valid.

Also note the variation of the image function F (p) as the point p goes to infinity. A sufficient theorem for what follows is

Theorem 3. If the integral (1.1.1) converges absolutely for the value p = Po = 09 + ito, then F (p) tends to zero as p goes to infinity by any law as long as p remains in the half- plane Rep = 06 > Oo.

Proof. Let p = o + it and o > oy. By assumption, the integral (1.1.1) is absolutely convergent for py = oo + it, and so it will converge for the value of p that is taken. Take a positive number A whose value will be determined

20

below and split the integral (1.1.1) into two summands: 00 A co F (p)= \ f(t) e-Ptdt= \ f (t) e-P# dt + \ f(tje-P?edt=,+ 1, 0 0 A First estimate the second integral:

[To1< J [f(t e-Pet [e-(o-oo dex | f(t) e-Pet| dt A A

If we have an arbitrary positive number ¢, the number A can be chosen so that the last term of the inequalities is less than €/3. Then we have | J, | < e¢/3 for any value of p for which Re p = o > oy. Therefore it suffices to be sure that J, tends to zero aS p — oo.

The function f is taken to be absolutely integrable on any finite interval, say on [0, A], and therefore there is a function- g(t), defined and continuously differentiable on [0, Al, for which the following inequality holds true:

A \ LF —e let dt-< + 0

For the integral J, we have

A A =) @—g@lewdt+ | ged =1+l, 0 0 The integral J, may be estimated as follows:

A A al< J [f—e(t)le-stdt< | [ft)—g )|e-ot at< & 0 0 To estimate J,, first perform the integration by parts:

A

eee -t| I,= > & (tye? ;

A 1 C g(t) e-vt ty le enas

21

Then | Le1 <7 Ue (0) [+18 (A) fem]

A

1 ae +r J le’ le dt=—-—_- M

A

M =|g(0)|+]g(A)]e-%4+ \ | g’ (t) | e-%# dt 0

If |p|>3M/e, then |I,|<e/3 and

2 lhil<lsl+lnul< ze

Finally, if we take advantage of the estimate of J, obtained above, we can Say that when | p | > 3M/e the inequality

2 1

will hold for F (p), and since e is an arbitrary positive quantity, it follows that lim F(p) =0 (Re p = o> 0p).

poo

1.2 Complex integral for computing inverse Laplace transforms

We now give one of a number of possible general methods for computing the inverse Laplace transform when such an inversion is given with the aid of a complex integral in which the integration is carried out along some straight line parallel to the imaginary axis of the complex plane,

To obtain the rule for inversion we will proceed from the Fourier double integral. Given on the real axis ¢ (—co < <(t< oo) an arbitrary function g(t). We say that the function g is representable by a Fourier double complex integral if for all values of t (—co < t< oo) the following equation! holds true:

4 The conditions that are sufficient for representing a function by the Fourier integral are given in Chap. 7.

22

+le (t+0)+ ¢(t—0)]=— \ ett \ g(x)e-** dxdt (4.2.1)

This equation is called Fourier’s formula.

The inner integral with respect to the variable z is as- sumed to be absolutely convergent for all t, which is clearly equivalent to the absolute integrability of g(t). The outer integral is regarded in the sense of the principal value, that is, as the limit of the integral taken along the interval [—b, b] symmetric relative to the point t = 0 provided b — oo.

Let us establish a connection between the Fourier formu- la (1.2.1).and the inverse Laplace transform. Suppose that f(t) is an arbitrary original function and that for a certain value c the product f(t) e-*¢ is an absolutely integrable function on the half axis [0, oo). Suppose g (t) = f (t) e-*. Note that for negative ¢ the function f and, hence, g (t) are equal to zero:

We assume that g(t) can be represented as the Fourier integral (1.2.1). To simplify notation, we will assume that at the discontinuity points of g(t) the following equation holds:

g(t) = le (+0) +g (t—0))

This condition only slightly alters the problem, since the set of discontinuity points of f and g is of measure zero and altering their values on such a set will not affect the magni- tudes of the integrals. It clearly holds at all points of conti- nuity of g(t). The Fourier formula for the function g takes the form \

g(t)==— \ ere \ g(x) e-** dz dt —0o 0

or, if we revert to f,

23

co

fQ=—— \ ee+inyt \ f(x)e-+in=drdt (1.2.2) — 00 0 Let us consider the Laplace transform of the function f:

F (p)= | f(t)e-Ptdt 0

Due to the absolute integrability of f (¢) e-*t, on the straight line p = c+ it (—o < t< oc) the transform converges for all t and, hence, it will converge in the half-plane Rep =o-ctec. Besides, F (p) will be a regular analytic function for Rep =o>e.

The value of the image function F (p) is of interest when o =:

F (e+ it) = \ f (t) e-e+ine dt 0

The image function F (p) stands under the sign of the inner integral in (1.2.2) so that for f (¢) the equation

f@®= — F (c + it) eetivt dr holds true. This is exactly the representation of the original function f in terms of its image function F; it is preferably written with the aid of the complex variable p = o + it. On the straight line p = c + it (—co < tT < oo) we haye dp = idt, and the resulting equation becomes ctioo

f\=za J F(pyertap (1.2.3)

C—100 This is called Mellin’s formula.

The foregoing enables us to regard the following theorem aS proven.

24

Theorem 4. Let the original function f(t) be such that for any value c the function g(t) = f(t) e- satisfies the conditions:

(1) g(t) is absolutely integrable on the half axis [0, oo);

(2) g(t) is representable by the Fourier double integral.

Then the representation (1.2.3) of the original function f in terms of its image F is true; it is assumed here that the fol- lowing equation holds at the points of discontinuity of f:

f(t)=+ Uf E+0)+ fF (t—0)].

Thus the problem of computing the inverse Laplace trans- form reduces to the problem of computing the contour inte- gral (1.2.3) of a certain regular function. Since it is by far not always possible to find an exact expression of the inte- gral (1.2.3) in terms of known functions, one can attempt to set up rules for obtaining it numerically. But this is a hard problem since, firstly, the contour of integration in the integral (1.2.3) is infinite and, secondly, the integrand e?* oscillates on thedine of integration, these oscillations being the stronger, the greater the value assumed by the para- meter f¢.

On the other hand, the function F (p) under the integral sign is not an arbitrary function but is an image function with all the properties that were mentioned above. This somewhat simplifies computing the integral (1.2.3) since the properties of the image function F (p) may be taken into nccount beforehand when setting up rules for a numerical inversion of the Laplace transform.

Setting up rules for computing the Mellin integral (these (ake into account the properties of the image function) will be discussed in Chapters 4 to 6.

1.3 Representing functions by the Laplace integral

We list the conditions that are sufficient for a given ltinction of a complex variable, F (p), which is analytic in

20

the half-plane Re p > a, to serve as the image function of a certain original function; that is, so that it can be repre- sented by a convergent Laplace integral.

Let us first recall, without proof, certain facts from the theory of functions of a complex variable that we will need in the sequel.

Jordan’s lemma. /f on a certain sequence of arcs of circles CR: |p| =R,, Rep <c(R, > o, ¢ is fixed), the function F (p) tends to zero uniformly with respect to arg p, then for any positive t

lim \ F (p) ert dp =0 (1.3.1) Ti OO Ch If the same conditions hold on the sequence of arcs of circles CR, |p|=R,, Rep>c, then for any negative t

lim \ F (p) ert dp =0. (1.3.2)

Theorem 5 (Cauchy’s theorem), Jf a function f (z) is analy- tic in a simply connected domain D, then the integral of the function along any closed contour C lying in D is equal to zero:

J f(z) dz=0.

Theorem 6 (residue theorem). Let a single-valued function { (z) be continuous on the boundary C of a domain D and be everywhere analytic inside this domain, except for a finite number of singular points ay, A, .. +, A,- Then

\ f(2)dz=2ni 5\ res f (ax) (1.3.3) C

heat

where res f (a) is the residue of the function f at the singular point a.

26

The residue of the function f (z) at a pole of order n may be found from the formula

res f(a) =p lim Fy Mea)" Fla) (1.3.4)

For poles of the first order the formula (1.3.4) can be simpli- fied to res f (a) = lim [(z— a) f(z)] (1.3.5) za

Here, if in the neighbourhood of the point a the function f (z) is defined as the quotient of two analytic functions at this point,

_. § (8) and g (a) £0 while wp (z) has a zero of the first order at a (that is, p (a2) =0, w’ (a) =40), then (1.3.5) can be replaced by

— lim 2 (7-9) =) @(z) —_ (2) res f(a) = lim @—a)=lin Saree ew @ 43-6) z—a

Now let us return to the proof of the theorem on the repre- sentability of F (p) by the Laplace integral.

Theorem 7. Jf the function F (p) is analytic in the half- plane Re p > a and tends to zero as | p | > co in any half- plane Rep >c>a uniformly with respect to arg p, and the integral

0-+i00 F (p) dp (1.3.7) 0— too converges absolutely, then F (p) is the image function of c+ioo 4 f(t)=z | F(pyertdp (1.3.8) c—too

which means it can be represented by the convergent Laplace integral

27

F (p) = \ e-P#f (t) dt 0

for Rep > c; here the Laplace integral converges absolutely. Proof. Take a number py such that Re py > c; then from (1.3.8) it follows that

[ e-mey (t) dt=— [eons an ePt F (p) dp) dt (1.3.9) 0 0 Cc — 700

Consider the inner integral. In it, p =c + iy, dp =idy and hence it can be rewritten as

c+-i00 (ore)

\ ePt F' (p) dp == ie \ et F (c + iy) dy Let us estimate the last integral:

| J er (e+iy) dy|< | [F(et+ay)ldy (1.3.10)

By virtue of the conditions of the theorem, the integral (1.3.7) is absolutely convergent and so the integral on the left of (1.3.10) converges uniformly with respect to ¢ and, hence, we can interchange the order of integration in (1.3.9):

oo c-+ioo oo | e-Petf (t) dt= \ F (p) dp { e@-po dt 0 C—i00 0 , c-++i00 F(p) — : dp (1.3.11)

The last equation is true since the inner integral converges

due to the fact that Re(p — p,) <0 and t>0. Consider the arc CR: |p | = R, Rep >c. By the theo-

rem, on this arc max | F (p) | =ap—0 as Ro and

28

hence

| fh aot | CR

and this integral too tends to zero as R — oo. From this it follows that the straight-line path of integration in (1.3.11) can be replaced by the closed contour Cp made up of the arc Cp and the line segment [c + ib, c — ib] traversed downwards. Then the formula (1.3.11) can be written as:

j e~Petf (t) dt =

\ 0 Cp

F (P) J = dp (1.3.42) We drop the minus sign since we have changed the direction along the straight line. We compute the integral in the right member of (4 3.12) by using the residue theorem. The F (p)

function inside the contour Ce has only one singular

point—a pole of the first order at the point p = po. Its residue at this point can be computed from (1.3.6) and will be equal to F (py). ‘Then from (1.3.3) we find

co

\ e-Potf (t) dt = F (pp) (1.3:13) 0

And since py is any point of the half-plane Rep > c, it follows from (1.3.13) that F (p) is a convergent Laplace integral for all p for which Re p > c. Later on we will show that this integral is absolutely convergent as well.

We will now show that if the conditions of the theorem are fulfilled, then f (¢), which is represented by the integral (1.3.8), will have property (2) in the definition of an original function. Indeed, for t << 0, we have by the Jordan lemma

lim | e*F(p)dp=0 R--c ¥ CR 29

Hence, the straight-line path of integration in the formula

(1.3.8) may be replaced by the contour C, which was defined earlier. Then, for t< 0, we get, by the Cauchy theorem,

1 f(t) => | oF (p) ap =0

Cr since the integrand is analytic inside Cp. Hence, property (2) is fulfilled for f. Besides, we will show that a condition of the form (1.1.3) is fulfilled for the function f when a = c. Indeed, from (1.3.8) follows

' C-+i00

fl=|45 \ ert (p) dp|

6-100

‘

<—. e! \ | F(c+iy)|dy-= Met (4.3.14) so that the inequality (4.1.3) holds.

Let us now return to the integral (1.3.13) and show that it converges absolutely. True enough, let po = Xo + iYo, then by virtue of (1.3.14) and the fact that rz) >,

M

Ig—c

\ |e-Petf (¢) | dt< M \ e~(%o-elt dt = 0 | 0

1.4 Ill-conditioning of the problem of computing inverse Laplace transforms |

The problem of recovering the original function f (t) from the operational image function F (p) may be regarded as a problem in solving the first-order integral equation

\ e-Ptf (t) dt =F (p) (1.4.4) 0 30

which belongs to the class of what are called ill-conditioned problems. These problems have two properties that greatly complicate their solution: they are not solvable for all values of the numerical or functional parameters defining their solution, and small variations in these parameters may he associated with a large change in the solution (the insta- bility of such problems was mentioned in the preface and we are compelled to repeat some of the things said there). We now give a more detailed explanation of these properties for the problem of computing inverse Laplace transforms and for the associated integral equations.

For the sake of definiteness consider the case where the image function F (p) is known on the real half axis p >a and the argument p is a real variable. The equation (1.4.1) does not have a solution for all functions F (p) that are continuous or even smooth when p >a. In particular, it is unsolvable if F (p) isnot an analytic function when p >a. Now suppose that F (p) is the image function of some func- tion f and the equation (1.4.1) is, hence, solvable. Replace /(t) by some other function f/f, (¢) that differs from f (¢) by a perturbation of large magnitude on a rather small interval and is coincident with f (¢) on the remaining portion of the half axis [0, oo). The new original function f, will be associated with the image function F, (p) that differs but slightly from F (p) for arbitrary p >a. Hence, to a small change in the right-hand member of (1.4.1) there can corre- spond an arbitrarily large variation in the solution f in a uniform métric. It can be demonstrated that a similar situation will obtain in other metrics.

The ill-conditioning of the problem of solving equation (1.4.1), which is equivalent to computing the inverse Laplace transform, can complicate the numerical solution but does not make it impossible. To the solution of equation (1.4.1) we can apply the regularization methods developed by A. N. Tikhonov, V. K. Ivanov and others. We will not dwell on these methods here.

31

Chapter 2

Some Analytical Methods for Computing Inverse Laplace

Transforms

2.1 Finding the original function via the inversion formula

It was shown in Chap. 1 that if f (¢) is the original function and F (p) is the Laplace image function, then to compute the original function from the image we can take advantage of the complex integral

c+ioo

fQ=s> \ Flee ap (2.1.1)

C—ico

where c is the abscissa in the half-plane of the absolute convergence of the Laplace integral. True, it is rather difficult to utilize (2.1.1) for a direct computation of f (¢), since it requires a knowledge of the function F(p) for com- plex values p = c+ iy (—o < y < oc) and the integral is improper with an oscillating kernel. But since (2.1.1) is an integral of an analytic function taken along a contour in the complex plane, it can be transformed by applying familiar methods of the theory of functions of a complex variable, for example, by changing the path integration, computing residues, and the like. In certain cases such transformations permit obtaining an expression of practical convenience for the original function, from which we can