Samarskii A.A., Vabishchevich P.N. Numerical Methods for Solving Inverse Problems of Mathematical Physics

Подождите немного. Документ загружается.

126 Chapter 4 Boundary value problems for parabolic equations

Exercise 4.12 Modify the program PROBLEM3 by implementing in it the three-layer

difference scheme (4.132). Conduct computational experiments on the use of the

scheme for approximate solution of the model problem for equation (4.117), (4.133).

Exercise 4.13 Construct an automodel solution of the traveling-wave type for equa-

tion (4.117) with coefﬁcients

k(x, u) = u

, f (x, t, u) = 0.

Conduct numerical experiments to model such regimes with the use of the program

PROBLEM3.

Exercise 4.14 Numerically examine the development of local perturbations in the

Cauchy problem for equation (4.117) with

k(x, u) = u

, f (x, t, u) = u

σ +1

Exercise 4.15 Write, around the program PROBLEM2, which numerically solves dif-

ference elliptic equations, a program that numerically solves the linear Dirichlet

boundary value problem for the two-dimensional parabolic equation

∂u

∂t



α=1

∂

∂x



k(x)

∂u

∂x



= f (x, t), x ∈ , t > 0,

in the rectangle

 ={x | x = (x

, x

), 0 < x

< l

,α= 1, 2}.

5 Solution methods for ill-posed problems

Basic approaches for the approximate solution of ill-posed problems use this or that

perturbation of the initial problem; such a perturbation allows one to pass to some

“close” yet well-posed problem. In this way, various regularization algorithms can

be obtained. With the example of the ﬁrst-order operator equation, below we con-

sider basic approaches to the solution of unstable problems. In particular, an ill-posed

problem can be replaced with a well-posed problem by perturbing the initial equation,

by passing to a variational problem, etc. Of primary importance here is an appropriate

choice of the regularization parameter, whose value needs to be properly matched with

input-data inaccuracies. Illustrative calculations performed for the ﬁrst-kind integral

equation are given.

5.1 Tikhonov regularization method

Following A. N. Tikhonov, consider the general approach of constructing stable com-

putational algorithms for the approximate solution of ill-posed problems. The method

is based on the passage from the initial ﬁrst-kind equation to a minimization problem

for a functional with an additional stabilizing term.

5.1.1 Problem statement

Normally, methods for the numerical solution of ill-posed problems are considered as

applied to the ﬁrst-kind linear operator equation

Au = f. (5.1)

The right-hand side and the solution itself belong to some metric spaces.

To not overload the consideration with technical details, we restrict ourselves to the

case of a linear operator A under conditions acting in a Hilbert space H (for simplicity,

here we assume that u ∈ H , f ∈ H, i.e., A : H → H ). In H , we introduce the scalar

product (u,v)and the norm u for elements u,v ∈ H .

Even more, we assume that the operator A is self-adjoint, positive and compact,

with eigenvalues λ

tending to zero as k →∞(λ

≥ λ

≥···≥λ

≥···> 0),

and the related eigenfunction system {w

}, w

∈ H , k = 1, 2,... is an orthonormal

system complete in H . Hence, for each v ∈ H we have

v =

∞



k=1

= (v, w

Typical of practical studies is a situation in which input data are given with some

inaccuracy. To model this general situation, consider a case in which the right-hand

128 Chapter 5 Solution methods for ill-posed problems

side of (5.1) is given accurate to some δ. Instead of f , we know f

such that

 f

− f ≤δ. (5.2)

In a more general case, we consider problems (5.1) in which not only the right-hand

side but also the operator A are known approximately.

We pose a problem in which it is required to ﬁnd an approximate solution of equa-

tion (5.1) with an approximately given right-hand side f

. This approximate solution

is denoted by u

, and the parameter α can be naturally related with the inaccuracy

level in the right-hand side, i.e., α = α(δ).

5.1.2 Variational method

The main idea behind the construction of stable methods for solving ill-posed problems

is based on the use of some a priori information about the input-data inaccuracy. Once

the right-hand side is given with an inaccuracy, no attempts to solve the equation

= f

(5.3)

exactly are necessary. We can try to compensate for the uncertainty in the right-hand

side, for instance, by passing to some another yet well-posed problem

= f

whose operator A

possesses improved properties compared to A.

In variational methods, instead of solving equation (5.3), they minimize the norm

of the difference r = Av − f

,orthediscrepancy functional

(v) =Av − f



There are many solutions to such a variational problem that satisfy equation (5.3) ac-

curate to some discrepancy δ. What is only needed is to wisely take into account all

available information about the inaccuracy in the right-hand side and, in this way, to

single out a most appropriate solution.

In the Tikhonov regularization method, they introduce some smoothing functional

(v) =Av − f



+ αv

. (5.4)

The approximate solution of the initial problem (5.1), (5.2) is the extremal of the func-

tional:

) = min

v∈H

(v). (5.5)

In (5.4), α>0istheregularization parameter, whose value must be matched

with the right-hand side inaccuracy δ. For a bounded solution to be separated out, the

discrepancy functional contains an additional stabilizing functional v

Section 5.1 Tikhonov regularization method 129

The key point in the theoretical consideration of approximate algorithms is related

with the proof of convergence of the approximate solution to the exact solution. It is

required to ﬁnd under which conditions the approximate solution u

found from (5.4),

(5.5) converges to the exact solution of problem (5.1). Ideally, it is required not only to

establish the fact that convergence takes place but also to ﬁnd the rate of convergence.

We represent the approximate solution in the operator form

= R(α) f

. (5.6)

If the approximate solution converges to the exact solution with the right-hand side

inaccuracy tending to zero, then they say that the operator R(α) is a regularizing op-

erator. With the chosen structure of R(α), it is also required to substantiate the choice

of regularization parameter α as a function of δ.

5.1.3 Convergence of the regularization method

In a consideration of conditionally well-posed problems, it is required to identify the

class of desired solutions and explicitly give a priori constraints on the solution. In

problem (5.1), (5.2), we are interested in bounded solutions and, hence, the a priori

constraint on the solution looks as

u≤M, (5.7)

where M = const > 0. The main result here can be formulated as follows.

Theorem 5.1 Let for the right-hand side inaccuracy the estimate (5.2) holds. Then,

the approximate solution u

found as the solution of problem (5.4), (5.5) with α(δ) →

0, δ

/α(δ) → 0 converges in H , as δ → 0, to the bounded exact solution u.

Proof. Under the indicated assumptions concerning the operator A, the exact solution

of (5.1) can be represented as

u =

∞



k=1

( f,w

. (5.8)

Suppose that

v =

∞



k=1

, c

= (v, w

Then, the functional J

(v) takes the form

(v) =

∞



k=1



(λ

− ( f

))

+ αc



130 Chapter 5 Solution methods for ill-posed problems

Condition (5.5) is equivalent to

∂ J

∂c

= 2λ

(λ

− ( f

)) + 2c

= 0, k = 1, 2,... .

From here, we obtain the following representation for the solution of problem (5.4),

(5.5):

∞



k=1

+ α

( f

. (5.9)

For the inaccuracy z = u

− u, we use the representation

z = z

+ z

, z

= u

− R(α) f, z

= R(α) f − u. (5.10)

In (5.10), R(α) f is the solution to the minimization problem for the smoothing func-

tional with the right-hand side given exactly.

In view of (5.8), (5.9), we obtain:

z



∞



k=1

(λ

+ α)

(( f

) − ( f,w

))

For non-negative x we have:

+ α

(

√

x −

√

α/x)

+ 2

√

≤

√

and, hence, under the assumptions of (5.2) there holds the two-sided inequality

z



≤

4α

∞



k=1

(( f

) − ( f,w

))

≤

4α

. (5.11)

Estimate (5.11) expresses the fact that the solution of problem (5.4), (5.5) is stable

with respect to weak right-hand side perturbations.

Let us estimate now z

in the representation (5.10) for the approximate-solution

inaccuracy. In this case, we deal with the conditions of closeness of the solution of

(5.1) to the solution of the minimization problem for the smoothing functional with

the same right-hand sides. From (5.8), (5.9) we have:

z



∞



k=1

(λ

+ α)

( f,w

)

We are going to establish the closeness of R(α) f to u for functions from class (5.7).

Let us show that for an arbitrary ε>0 we can choose α(ε) such that z



≤ ε for all

α from the interval 0 <α≤ α(ε). For all functions from (5.7), the series

u

∞



k=1

( f,w

)

Section 5.2 The rate of convergence in the regularization method 131

converges and, hence, we can found a number n(ε) such that

∞



k=n(ε)+1

( f,w

)

≤

Under the stated conditions, we obtain the inequality

z



≤

n(ε)



k=1

(λ

+ α)

( f,w

)

n(ε)



k=n(ε)+1

( f,w

)

At the expense of a sufﬁciently small chosen α(ε), for the ﬁrst term we have:

n(ε)



k=1

(λ

+ α)

( f,w

)

≤

The latter inequality holds for all α in the interval 0 <α≤ α(ε). Hence, we obtain

that s(α) =z

→0asα → 0.

Substitution into (5.10) yields

z≤z

+z

≤

√

+ s(α). (5.12)

By virtue of the above, if α(δ) → 0 and δ

/α(δ) → 0asδ → 0, then z→0.

The proved theorem remains valid under even more general conditions. A funda-

mental generalization here can be obtained considering problems (5.1), (5.2) whose

operator A is not a self-adjoint operator.

5.2 The rate of convergence in the regularization method

In the class of bounded solutions, we have established previously that the approximate

solution in the Tikhonov regularization method converges to the exact solution. Under

more tight constraints on the exact solution, one can evaluate the rate at which the

method converges.

5.2.1 Euler equation for the smoothing functional

Instead of the extremum problem (5.4), (5.5), we can consider a related Euler equation.

In the latter case, the approximate solution can be found as the solution of the following

second-kind equation:

∗

+ αu

= A

∗

. (5.13)

The transfer from the ill-posed problem (5.3) to the well-posed problem (5.13) can

be made passing to a problem with a self-adjoint operator A

∗

A. This can be done by

132 Chapter 5 Solution methods for ill-posed problems

multiplying equation (5.3) from the left by A

∗

, followed by a subsequent perturbation

of the resulting equation with the operator α E.

The mentioned equivalence between the minimization problem (5.4), (5.5) and

equation (5.13) provides us with some freedom in the computational realization of

the method. The variational approach is quite general, allowing one to consider var-

ious classes of problems from a uniﬁed methodological standpoint. The numerical

solution can most conveniently be found using the related Euler equation.

In the case of A = A

∗

≥ 0, we can restrict ourselves to perturbing the operator

itself:

+ αu

= f

. (5.14)

Problem (5.14) refers to the simpliﬁed regularization algorithm.

In fact, we can say that, in addition to variational solution methods for ill-posed

problems, a second class of approximate methods can be identiﬁed, which (see, for

instance, (5.13), (5.14)) features some perturbation of the operator in the initial or

transformed problem.

5.2.2 Classes of a priori constraints imposed on the solution

We have considered the Tikhonov regularization method under the assumption (5.7)

about the exact solution of the ill-posed problem (5.1), (5.2). For the solution inaccu-

racy we have obtained the estimate (5.12), in which s(α) → 0 provided that α → 0.

Note that the convergence was established in the same norm in which the a priori

constraints on the solution were formulated. We would also like to obtain an esti-

mate for the rate of convergence of the approximate solution to the exact solution with

α(δ) → 0 and δ → 0.

To clarify the situation, consider the problem of approximate calculation of a deriva-

tive. We use the difference relation

◦

u(x + h) − u(x − h)

, (5.15)

and assume that the function u(x) is differentiable at each x. We would like to know

how accurately the difference derivative (5.15) approximates the derivative du/dx at

some point x.

To prove the convergence and evaluate the rate with which the central difference

derivative (5.15) converges to du/dx(x), let us formulate more tight constraints on

u(x). If, for instance, u(x) is a twice differentiable function, then from (5.15) we

obtain

◦

+ O(h).

For a three times differentiable function we have

◦

+ O(h

Section 5.2 The rate of convergence in the regularization method 133

Thus, with reduced smoothness of differentiable functions the approximation inaccu-

racy decreases.

We encounter a similar situation when examining the convergence rate in the regu-

larization method. Instead of (5.17), we are going to formulate more tight constraints

on the exact solution of problem (5.1). It seems reasonable to associate the require-

ment for enhanced smoothness of the exact solution with the operator A. We assume

that the exact solution belongs to the class

A

−1

u≤M. (5.16)

Then, in the case of a self-adjoint positive operator A the intermediate (between (5.7)

and (5.16)) position belongs to the following type of classes of a priori constraints on

the solution:

u

−1

≤ M. (5.17)

Under conditions (5.16) and (5.17), one can try to render more concrete the depen-

dence s(α) in the inaccuracy estimate of type (5.12) in the regularization method.

5.2.3 Estimates of the rate of convergence

We now derive estimates for the solution inaccuracy in the Tikhonov regularization

method based on the a priori estimates of the operator equation (5.13) which, in view

of A = A

∗

> 0, assumes the form

+ αu

= Af

. (5.18)

Theorem 5.2 For the inaccuracy z = u

− u of the approximate solution of prob-

lem (5.1), (5.2) found from (5.13), there hold the a priori estimate

z

≤

2α

, (5.19)

for exact solutions from class (5.16) and the estimate

z

≤

√

(5.20)

for solutions from (5.17).

Proof. We subtract the equation

u = Af,

from (5.18); this yields the following equation for the inaccuracy z = u

− u:

z + αz = A( f

− f ) − αu.

134 Chapter 5 Solution methods for ill-posed problems

We scalarwise multiply this equation by z and obtain

Az 

+ αz

= (( f

− f ), Az ) − α(u, z).

Then, insertion of the inequality

(( f

− f ), Az ) ≤

Az 

 f

− f 

gives

Az 

+ αz

≤

 f

− f 

− α(u, z). (5.21)

Consider ﬁrst the case of a priori constraints (5.16). The use of the inequality

−α(u, z) =−α(A

−1

u, Az) ≤

Az 

A

−1

u

in (5.21) leads us to the estimate

αz

≤

 f

− f 

A

−1

u

With inequalities (5.2) and (5.16) taken into account, we arrive at the desired estimate

(5.19).

In the case of (5.17) the last term in the right-hand side of (5.21) can be estimated

−α(u, z) =−α(A

−1/2

u, A

1/2

z) ≤

√

α A

1/2

z

√

A

1/2

u

Taking the fact into account that

Az 

+ αz

−

√

α A

1/2

z



√

Az −

√



z

we obtain

z

≤

 f

− f 

√

A

1/2

u

In view of (5.2), (5.17), we have the estimate (5.20).

Tightened constraints on the smoothness of the exact solution result in an improved

rate of convergence of the approximations to the exact solution (see estimate (5.19),

(5.20)).

5.3 Choice of regularization parameter

Below several choices of the regularization parameter, whose value must be matched

with the input-data inaccuracy, are outlined.

Section 5.3 Choice of regularization parameter 135

5.3.1 The choice in the class of a priori constraints on the solution

In the theory of approximate solution methods for ill-posed problems, considerable at-

tention is paid to the choice of regularization parameter. The most widespread choices

here are the choices based on the discrepancy between the solutions, on the general-

ized discrepancy (which takes into account both the right-hand side inaccuracy and the

inaccuracy in the operator A), the quasi-optimal choice, etc. A good choice of regu-

larization parameter largely determines the efﬁciency of the computational algorithm.

The value of the regularization parameter α must be matched with the input-data

inaccuracy: the smaller is the inaccuracy, the smaller is the value of regularization

parameter that is to be chosen, i.e., α = α(δ). To comply with the structure of the

inaccuracy (see (5.12), (5.19), (5.20)), the regularization parameter cannot be cho-

sen too small: it is in the fact that, with decreased value of regularization parameter,

the inaccuracy also decreases, that the ill-posedness of the problem is manifested.

Hence, there exists an optimal value of regularization parameter that minimizes the

approximate-solution inaccuracy.

The optimal value of regularization parameter depends not only on the inaccuracy

in the right-hand side, but also on the class of a priori constraints on the exact solution.

For instance, in the case of bounded solutions (class (5.7)) the above estimate (5.12)

for the approximate-solution inaccuracy in the Tikhonov method does not allow one

to explicitly give the optimal value of regularization parameter.

On narrowing the class of exact solutions, it becomes possible to render concrete the

choice of regularization parameter. In the class of exact solutions (5.16) there holds the

a priori estimate (5.19) for the inaccuracy, and for the optimal value of regularization

parameter we obtain the expression

opt

, A

−1

u≤M. (5.22)

With this value of regularization parameter, a rate

z≤

√

Mδ

of convergence of the approximate solution to the exact solution can be achieved.

A similar consideration for the choice of regularization parameter in the class (5.17)

of a priori constraints on the exact solution leads us to

opt





2/3

, u

−1

≤ M (5.23)

with

z≤

√



To summarize, we can formulate the following statement.