Convergence and adaptive discretization of the IRGNM Tikhonov and the IRGNM Ivanov method under a tangential cone condition in Banach space.

In this paper we consider the iteratively regularized Gauss-Newton method (IRGNM) in its classical Tikhonov version as well as two further-Ivanov type and Morozov type-versions. In these two alternative versions, regularization is achieved by imposing bounds on the solution or by minimizing some regularization functional under a constraint on the data misfit, respectively. We do so in a general Banach space setting and under a tangential cone condition, while convergence (without source conditions, thus without rates) has so far only been proven under stronger restrictions on the nonlinearity of the operator and/or on the spaces. Moreover, we provide a convergence result for the discretized problem with an appropriate control on the error and show how to provide the required error bounds by goal oriented weighted dual residual estimators. The results are illustrated for an inverse source problem for a nonlinear elliptic boundary value problem, for the cases of a measure valued and of an L∞ source. For the latter, we also provide numerical results with the Ivanov type IRGNM.

Introduction

In this paper we consider a nonlinear ill-posed operator equationwhere the possibly nonlinear operator with domain maps between real Banach spaces X and Y. We are interested in the ill-posed situation, i.e., F fails to be continuously invertible, and the data are contaminated with noise, thus regularization has to be applied (see, e.g., [8, 27], and references therein).

Throughout this paper we will assume that an exact solution of (1) exists, i.e., , and that the noise level in the (deterministic) estimateis known.

Partially we will also refer to the formulation of the inverse problem as a system of model and observation equation Here and are the model and observation operator, so that with the parameter-to-state map satisfying and , (1) is equivalent to the all-at-once formulation (3), (4).

Newton type methods for the solution of nonlinear ill-posed problems (1) have been extensively studied in Hilbert spaces (see, e.g., [2, 20] and the references therein) and more recently also in a in Banach space setting. In particular, the iteratively regularized Gauss–Newton method [1] can be generalized to a Banach space setting by calculating iterates in a Tikhonov type variational form assee, e.g., [11, 16, 17, 21, 28] where , is a sequence of regularization parameters, and is some nonnegative regularization functional. Alternatively, one might introduce regularization by imposing some bound on the norm of x, or, again, generally, on a regularization functional of xwhich corresponds to Ivanov regularization or the method of quasi solutions, see, e.g., [7, 13–15, 22, 24, 26]. A third way of incorporating regularization in a Newton type iteration is Morozov regularization, also called the method of the residuals, see, e.g., [9, 22, 23]for some , where the choice of the bound in the inequality constraint is very much inspired by the inexact Newton type regularization parameter choice in [10].

We restrict ourselves to the norm in Y as a measure of the data misfit, but the analysis could as well be extended to more general functionals satisfying certain conditions, as e.g., in [11, 28]. Here is a set (possibly chosen with convenient properties for carrying out the minimization) containing and being contained in , such that F satisfies additional conditions on see (8), (10) below. If F is defined on all of X, then the minimization problem (5) can be posed in an unconstrained way .

As a restriction on the nonlinearity of the forward operator F we impose the tangential cone condition(also called Scherzer condition, cf. [25]) for some constant . Here, for any ,is a sublevel set of the regularization functional and R will be specified in the convergence result Theorem 1.

Note that the convergence conditions imposed in [11, 16, 17, 21, 28] in the situation without source condition, namely local invariance of the range of , are slightly stronger, since this adjoint range invariance is sufficient for (8). However, most probably the gap is not very large, as in those application examples where (8) has been verified, the proof of (8) is actually often done via adjoint range invariance. In (5), (6), (7), the bounded linear operator is not necessarily a Gâteaux or Fréchet derivative of F, but just some local linearization (in the sense of (8)), satisfying additionally the weak closedness conditionIn here, and are topologies on X and Y (e.g., just the weak or weak* topologies) such that bounded sets in Y are -compact and the norm in Y is -lower semicontinuous.

The remainder of this paper is organized as follows. In Sect. 2 we state and prove convergence results in the continuous and discretized setting. Section 3 shows how to actually obtain the required discretization error estimates by a goal oriented weighted dual residual approach and Sect. 4 illustrates the theoretical findings by an inverse souce problem for a nonlinear PDE. In Sect. 5 we provide some numerical results for this model problem and Sect. 6 concludes with some remarks.

Convergence

In this section we will study convergence of the IRGNM iterates first of all in a continuous setting, then in the situation of having discretized for computational purposes.

The regularization parameters , , are chosen a priori(note that for ),andwith as in (14), and the iteration is stopped according to the discrepancy principlewith some fixed chosen sufficiently large but independent of .

Theorem 1

Let be proper, convex and lower semicontinuous with and let, for all in case of (5), or for all in case of (6), or for in case of (7), the sublevel set (9) be compact with respect to the topology on X.

Moroever, let F satisfy (8), (10).

Finally, let the family of data satisfy (2).

Then for fixed , , the iterates according to (5)–(7) are well-defined and satisfy for all , which denotes the stopping index according to the discrepancy principle (14) with sufficiently large, and this stopping indes is finite.

Moreover, for both methods we have -subsequential convergence as i.e., has a -convergent subsequence and the limit of every -convergent subsequence solves (1). If the solution of (1) is unique in , then as .

Additionally, satisfies the asymptotics .

Proof

Existence of minimizers of (5)–(7) for fixed k, and follows by the direct method of calculus of variations: In all three cases, the cost functionalis bounded from below and the admissible setis nonempty (for (6) this follows from and for (7) from (8), (14) and (13), see (16) below). Hence, there exists a minimizing sequence forwith bounded linearized residuals forand . By -compactness of , the sequence has a -convergent subsequence with limit . Moroever, -compactness of norm bounded sets in Y together with --closedness of and lower semicontinuity of the norm in Y, implies that in all three cases and , hence is a minimizer.

Note that (ii) follows from (i) by standard arguments and our assumption on -compactness of . Thus it remains to prove (i) and (iii) for the three versions (5), (6), (7) of the IRGNM.

For this purpose we are going to show that for every , there exists such that , and the stopping criterion according to the discrepancy principle is satisfied. For (5), we also need to show that for , whereas in (6) this automatically holds by (12). The same holds true for (7): If , then by (8), (14) and (13) we haveso is admissible, hence , i.e., .

We start with the Tikhonov version (5) and carry out an induction proof of the following statement: For all where for some fixed small . We will require , which by definition of q (19) is achievable for sufficiently small, due to , cf. (11). By Lemma 2 (see the “Appendix”) the right hand side estimate in (17) impliesUsing the minimality of and (2), (8) together with , we haveFrom (21) and (14) we inferUsing this and again (14) in (22) yieldsOn the other hand, since we have established , we can apply (8) to the left hand side of (22) to obtainsee also [11, Lemma 5.2] and [21, proof of Theorem 3].

To handle the power p we make use of the following inequalities that can be proven by solving extremal value problems, see the “Appendix”for all and , where for the right hand inequality to hold, additionally is needed.

Hence, in case the following general estimate holdsfor

So in order for this recursion to yield geometric decay of , we need to ensurefor a proper choice of . To obtain the largest possible (and therefore least restrictive) bound on , we rewrite the requirement above asas can be found out by evaluating the derivative of Thus we will furtheron set and assume that is sufficiently small so that (26) holds with , i.e., (19). Then, using (11), estimate (25) can be written aswhich we first of all regard as a recursive estimate for .

To derive a similar estimate also in the complementary case , we use that fact that, for as in (18), this inequality just meansand, using (22) and the left hand part of (24),hence after addition we again get (27) (even with a slightly smaller value of ).

Thus in both cases, using Lemma 2 we can conclude thatThis finishes the induction proof of (17) for all .

We next show that the discrepancy stopping criterion from (14), i.e., for , will be satisfied after finitely many, namely , steps. For this purpose, note that , provided is chosen sufficiently large, which we assume to be done. Thus, indeed, using (11), (21), we havewhere the right hand side falls below as soon asThus we get the upper estimate .

For the Ivanov version (6), it only remains to show finiteness of the stopping index, as boundedness of the values by holds by definition. Applying the minimality argument with being admissible [cf. (12)] to (6) leads to the special case , in (25)Our notation becomeswhich givesand by induction, one can concludewhere the right hand side is smaller than (with ) for allso that we can again conclude .

Finally we consider (7), where boundedness of the values by holds by minimality and the fact that is admissible, cf. (16). Geometric decay of the residuals follows by the estimateand (13), i.e.,withso that similarly to above we end up with a logarithmic estimate for .

Remark 1

Convergence of to as holds along the convergent subsequence according to Theorem 1 (ii), first of all for the Morozov and the Ivanov version of the IRGNM, with the choice for the latter, since in both cases holds for all and is lower semicontinuous. The same holds true also for the Tikhonov version with the alternative choice of such thatfor some constants , satisfying in place of (11), as can be seen directly from (22). If is defined by the norm on a space with the Kadets-Klee property, and is the weak topology of this space, then this implies norm convergence of to along the same subsequence.

Remark 2

The fact that stays in [cf. (15)] is crucial for the applicability of the tangential cone condition (8) in these iterates. If the functional quantifies some distance to an a priori guess , (e.g., for some norm and some ), then with small R means closeness of x to in a certain sense. Thus, the smaller R is, the better (8) might get achievable with some . On the other hand, making R according to (15) small means closeness of to . Thus we deal with local convergence, as typical for Newton type methods.

Now we consider the appearance of discretization errors in the numerical solution of (5), (6) arising from restriction of the minimization to finite dimensional subspaces and leading to discretized iterates and an approximate version of the forward operator i.e., we consider the discretized version of Tikhonov-IRGNM (5)of Ivanov-IRGNM (6)and of Morozov-IRGNM (7)respectively. Moreover, also in the discrepancy principle, the residual is replaced by its actually computable discretized versionWe define the auxiliary continuous iterates andrespectively in order to be able to use minimality, i.e., compare with the continuous exact solution . For an illustration we refer to [18, Figure 1].

First of all, we assess how large the discretization errors can be allowed to still enable convergence. Later on, in Sect. 3, we will describe how to really obtain such estimates a posteriori and to achieve the prescribed accuracy by adaptive discretization.

Corollary 1

Let the assumptions of Theorem 1 be satisfied and assume that the discretization error estimates (note that no absolute value is needed in (38), (40); moreover, (40) is only be needed for (5) and (7)) hold withfor all and constants sufficiently small, .

Then the assertions of Theorem 1 remain valid for in place of with (34) in place of (14) and (42) in place of (23).

For the Tikhonov version (31), in order to inductively estimate , given , note that from (43) with replaced by k, we get like in (23) thatwherefor , which are chosen small enough so that . As before, from the minimality of and (2), (8) as well as , we havethen using (38), (40),Hence, with the same technique as in the proof of Theorem 1, using (24) with , we haveusing (41). From this, by induction we concludeHence, by (39), (41), we have the following estimatewhere the right hand side falls below as soon asfor . Note that , provided is chosen sufficiently large, which we assume to be done. That is, we have shown that the discrepancy stopping criterion from (34) will be satisfied after finitely many, namely , steps.

On the other hand, the continuous discrepancy at the iterate defined by the discretized discrepancy principle (34) by (39), (41) satisfiesTo estimate , note that according to our notation, from (43), we get, like in (23), that for all Now we show finiteness of the stopping index for the discretized Ivanov-IRGNM (32). By minimality of and (38), for this problem we havewhich withby induction, (39) and (41) giveswhere the right hand side is smaller than for allwith , so that we can again conclude .

It remains to show finiteness of the stopping index for the discretized Morozov-IRGNM (33). By minimality of we have (30) with replaced by , thus the inequalities (38) and (41) yieldthen, by (39) and inductionwhereand the right hand side of (44) falls below for allwhere , and we can again conclude .

Boundedness of the values for (33) by follows like in the proof of Theorem 1 together with (40), (41).

Error estimators for adaptive discretization

The error estimators , and can be quantified, e.g., by means of a goal oriented dual weighted residual (DWR) approach [3], applied to the minimization problems(note that the last constraint is added in order to enable computation of below)andwhich are equivalent to (5), (6), and (7), respectively, withas quantities of interest [where is only needed for (5) and (7)]. We assume that and the norms can be evaluated without discretization error, so the discretized versions of only arise due to discreteness of the arguments. Indeed, it is easy to see that the left hand sides of (38) and (39) can be bounded (at least approximately) by combinations of and , using the triangle inequality: where we will neglect .

It is important to note that is not equal to , see [18].

The computation of the a posteriori error estimators is done as in [18]. These error estimators can be used within the following adaptive algorithm for error control and mesh refinement: We start on a coarse mesh, solve the discretized optimization problem and evaluate the error estimator. Thereafter, we refine the current mesh using local information obtained from the error estimator, reducing the error with respect to the quantity of interest. This procedure is iterated until the value of the error estimator is below the given tolerance (41), cf. [3].

In this case, all the variables are subject to a new discretization. For better readability we will partially omit the iteration index k and the discretization index h. The previous iterate is fixed and not subject to a new discretization.

Consider now the cost functional for (45)and define the Langrangian functionalneglecting for simplicity (cf. Remark 2) the constraints defined by . The first-order necessary optimality conditions for (45) are given by stationarity for the Lagrangian L. Setting , they readand for the discretized problem,To derive a posteriori error estimators for the error with respect to the quantities of interest (), we introduce auxiliary functionals :Let and be continuous and discrete stationary points of satisfyingrespectively. Then, are continuous and discrete stationary points of L and there holds . Thus the z part, as computed already during the numerical solution of the minimization problem (45) (or (46)) remains fixed for all . Moreover, after computing the discrete stationary point for L (e.g., by applying Newton’s method), it requires only one more Newton step to compute the coordinate of the stationary point for M fromAccording to [3], there holdswith a remainder term R of order that is therefore neglected. Thus we usewhere is an operator defined such that approximates the interpolation error as in [18], typically defined by local averaging, to define the estimators , , according to the rulecf. (48), (49). The estimators obtained by this procedure can be used to trigger local mesh refinement until the requirements (41) are met cf. [3].

Explictly, for (for simplicity) such a stationary point can be computed by solving the following system of equations (analogously for the discrete stationary point of L) Note that (58) is decoupled from the other equations and that if is injective, Eq. (54) implies .

Summarizing, since we have a convex minimization problem, after solving a nonlinear system of seven equations to find the minimizer, we need only one more Newton step to compute the error estimators to check whether we need a refinement on the mesh or not.

Regarding the problem (46) related to the Ivanov-IRGNM, we have the Lagrangian functional (50) with the cost functional defined byand the indicator functional takes the role of a regularization functional. The resulting optimality system is the same as above, cf. (52)-(58), just with (52) replaced bySimilarly for (47) for Morozov-IRGNM, with the cost functionwe end up with an optimality system by setting and replacing (53), (555657) in (52)–(58) by respectively.

Note that the bound on only appears—via (51)—in connection to the assumption , for in (41). This may be satisfied in practice without refining explicitly with respect to , but simply by refining with respect to the other error estimators (and in the Tikhonov or Morozov case). The fact that and only differ in the discretization level, motivates the assumption that for small h, we have and . Thefore, the algorithm used in actual computations will be built neglecting and hence skipping the constraint in (45), (46), (47), which implies a modification of the Lagrangian (50) accordingly. Therefore, the corresponding optimality systems for in the Tikhonov case is given by Note that Eq. (66) is decoupled from the others. Therefore, the strategy is to solve (66) for first, then solve the linear system (62), (63), (65) for , and finally compute via the linear equation (64). Here, the system (62), (63), (65) can be interpreted as the optimality conditions for the following problemFor the Ivanov case, we have to solve (63)–(66) within place of (62), hence again (66) is decoupled from the other equations, (64) is linear with respect to , once has been computed, and the remaining system for can be interpreted as the optimality conditions for the following problemThe Morozov case requires solution of (62) (with ), (65), (66), (60), (61). Thus again, we first solve (66) for , then the system (62) (with ), (65), (60), which is the first order optimality condition forwith Lagrange multiplier for the equality constraint, and finally the (now possibly nonlinear) inclusion (61) for .

Remark 3

Since DWR estimators are based on residuals which are computed in the optimization process, the additional costs for estimation are very low, which makes this approach attractive for our purposes. However, although these error estimators are known to work efficiently in practice (see [3]), they are not reliable, i.e., the conditions , can not be guaranteed in a strict sense in the computations, since we neglect the remainder term R and use an approximation for . As our analysis in Theorem 1 is kept rather general, it is not restricted to DWR estimators and would also work with different (e.g., reliable) error estimators.

Model examples

We present a model example to illustrate the abstract setting from the previous section. Consider the following inverse source problem for a semilinear elliptic PDE, where the model and observation equations are given by where denotes the extension by zero of a function on to a function on all of . We first of all consider Tikhonov regularization and, aiming for a sparsely supported source, therefore use the space of Radon measures as a preimage space X. Thus we define the operators , , and the injection , , where is a bounded domain in with or 3, with Lipschitz boundary and are the control domain and the observation domain, respectively.

A monotonicity argument yields well posedness of the above semilinear boundary value problem, i.e., well-definedness of as a solution to the elliptic boundary value problem (68), (69), as long as we can guarantee that for any , i.e., the embeddings and are continuous for some , which (by duality) is the case iff embeds continuously both into and . By Sobolev’s Embedding Theorem, this boils down to the inequalitieswhich by elementary computations turns out to be equivalent towhere the left hand side is larger than one and the denominator on the right hand side is positive due to the fact that for we have . Taking the extremal bounds for —note that the lower bound is increasing and the upper bound is decreasing with q—in (72) we getThus, as a by-product, we get that for any there exists such that continuously embeds into , withFor the regularization functional , the IRGNM-Tikhonov minimization step is given by (ignoring h in the notation)Here and below and denote the integrals with respect to the Lebesgue measure and with respect to the measure , respectively.

Therefore, to compute this Gauss–Newton step, one first needs to solve the nonlinear equationfor , then solve the following optimality system with respect to (written in a strong formulation)which can be interpreted as the optimality system for the minimization problemwith Lagrange multiplier for the equality constraint, and finally, compute by solvingFor carrying out the IRGNM iteration, is not required, but we need it for evaluating the error estimators.

For the Ivanov case, we consider the same model and observation equations (68), (69), (70) but now we intend to regularize by imposing bounds and thus use the slightly different function space setting, , , and the injection .

The IRGNM-Ivanov minimization step with the regularization functional is given byFor the Gauss–Newton step, one needs to first solve the nonlinear equation (74) for , and then solve the following optimality system with respect to which can be interpreted as the optimality system for the minimization problemwith Lagrange multiplier for the equality constraint. Finally, is computed from (76).

For the IRGNM-Morozov case, using for simplicity the regularization functional , and leaving the rest of the setting as in the IRGNM-Ivanov case, the step is defined bySo again we first solve (74) for , then the minimization problemor actually its first order optimality systemfor , and finally,for .

For numerically efficient methods to solve the minimization problems (75) and (77) we refer to e.g., [4-6] and the references therein.

We finally check the tangential cone condition in case and, for simplicity also , in both settings(where we will have to restrict ourselves to ) andFor this purpose, we use the fact that with the notation , , and , the functions satisfy the homogeneous Dirichlet boundary value problems for the equationsUsing an Aubin-Nitsche type duality trick, we can estimate the norm of w via the adjoint state , which solveswith homogeneous Dirichlet boundary conditions, so that by Hölder’s inequalitywhere we aim at choosing , such that indeedand the embeddings , , are continuous. If we succeed in doing so, we can bound by some constant , which will be small provided and hence is small. Thus, the numbers n, m are limited by the requirements, i.e., by duality,and the fact that should be contained in for , and , which via Hölder’s inequality inand duality leads to the requirementsIn case with and , (78) will not work out, since according to (73), m cannot be chosen larger or equal to four.

In case with and , we can choose, e.g., , to satisfy (78), (79), (80) as well as as in (73).

The same choice is possible in case with .

Numerical tests

In this section, we provide some numerical illustration of the IRGNM Ivanov method applied to the example from Sect. 4, i.e., each Newton step consists of solving (74) and subsequently (77). For the numerical solution of (74) we apply a damped Newton iteration to the equation wherewhich is stopped as soon as has been reduced by a factor of 1.e−4. The sources and states u are discretized by piecewise linear finite elements, hence after elimination of the state via the linear equality constraint, (77) becomes a box constrained quadratic program for the dicretized version of , which we solve with the method from [12] using the Matlab code mkr_box provided to us by Philipp Hungerländer, Alpen-Adria Universität Klagenfurt. All implementations have been done in Matlab.

We performed test computations on a 2-d domain , on a regular computational finite element grid consisting of triangles, with . We first of all consider (below we will also show results with ) and the piecewise constant exact source functionwhere B is the ball of radius 0.2 around cf. Fig. 1, and correspondingly set . In order to avoid an inverse crime, we generated the synthetic data on a finer grid and, after projection of onto the computational grid, we added normally distributed random noise of levels to obtain synthetic data .

Fig. 1

Left: exact source ; right: locations of spots for testing weak * convergence

In all tests we start with the constant function with value zero for . Moreover, we always set . According to our convergence result Theorem 1 with , we can expect weak * convergence in here. Thus we computed the errors in certain spots within the two homogeneous regions and on their interface,cf. Fig. 1, more precisely, on squares located at these spots, corresponding to the piecewise constant functions with these supports in order to exemplarily test weak * convergence. Additionally we computed errors.

Table 1 provides an illustration of convergence as decreases. For this purpose, we performed five runs on each noise level for each example and list the average errors.

Table 1

Convergence as : averaged errors of five test runs with uniform noise

δ	err_spot_1	err_spot_2	err_spot_3	err_L^1(Ω)
0.1000	0	4.0818	8.0043	0.0627
0.0667	0.1558	3.6454	7.8451	0.0541
0.0333	0	3.0442	6.5726	0.0370
0.0100	0	0	3.9091	0.0188

In Fig. 2 we plot the reconstructions for and . For , the noise levels correspond to a percentage of of the deviation of the exact state from the background state . In case of , where the background state is the corresponding percentages are . For an illustration of the noisy data as compared to the exact ones, see Figs. 3 and 4. Indeed, the box constraints enable to cope with relatively large noise levels, even in the rather nonlinear regime with .

Fig. 2

Reconstructions from noisy data with (top to bottom) for (left) and (right)

Fig. 3

Exact and noisy data () for

Fig. 4

Exact and noisy data () for

Conclusions and remarks

In this paper we have studied convergence of the Tikhonov type, the Ivanov type, and the Morozov type IRGNM with a stopping rule based on the discrepancy principle type. To the best of our knowledge, the Ivanov and Morozov IRGNMs have not been studied so far and in all three Tikhonov, Ivanov, and Morozov type IRGNMs, convergence results without source conditions so far use stronger assumptions than the tangential cone condition used here. We also consider discretized versions of the methods and provide discretization error bounds that still guarantee convergence. Moroever, we discuss goal oriented dual weighted residual error estimators that can be used in an adaptive discretization scheme for controlling these discretization error bounds. An inverse source problem for a nonlinear elliptic boundary value problems illustrates our theoretical findings in the special situations of measure valued and sources. We also provide some computational results with the Ivanov IRGNM for the case of an source. Numerical implementations and tests for a measure valued source, together with adaptive discretization is subject of ongoing work, based on the approaches from [4–6, 18, 19]. Future research in this context will be concerend with convergence rates results for the Ivanov and Morozov IRGNMs under source conditions.