Regularization of inverse problems by two-point gradient methods with convex constraints
Min Zhong Wei Wang Qinian Jin



In this paper, we propose and analyze a two-point gradient  method for solving inverse problems in Banach spaces which
is based on the Landweber iteration and an extrapolation strategy.  The method allows to use non-smooth penalty terms,
including the L^1 and the total variation-like penalty functionals, which are significant in reconstructing special
features of solutions such as sparsity and piecewise constancy in practical applications. The design of the method involves
the choices of the step sizes and the combination parameters which are carefully discussed. Numerical simulations are presented
to illustrate the effectiveness of the proposed method.



Magnetically tunable multi-band near-field radiative heat transfer
between two graphene sheets
     Ying Wu
    December 30, 2023
==============================================================================================






§ INTRODUCTION



In this paper we are interested in solving inverse problems of the form

    F(x)=y ,

where F:(F)⊂→ is a  Fréchet differentiable operator between two
Banach spaces  and . Throughout this paper we will assume (<ref>) has a solution,
which is not necessarily unique. Such inverse problems are ill-posed in the sense of unstable dependence
of solutions on small perturbations of the data. Instead of exact data y, in practice we are
given only noisy data y^$̣ satisfying

    y-y^≤.

Consequently, it is necessary to apply regularization methods to solve (<ref>) approximately (<cit.>).

Landweber iteration is one of the most prominent regularization methods for solving inverse problems
formulated in Hilbert spaces. A complete analysis on this method for linear problems as well as nonlinear
problems can be found in <cit.>. This method has received tremendous attention due to
its simple implementation and robustness with respect to noise.

The classical Landweber iteration in Hilbert spaces, however, has the tendency to over-smooth solutions, which makes it
difficult to capture special features of the sought solution such as sparsity and discontinuity. To overcome this drawback,
various reformulations of Landweber iteration either in Banach spaces or in a manner of incorporating general non-smooth
convex penalty functionals have been proposed, see <cit.>. In <cit.>, using
a gradient method for solving

    min1/sF(x)-y^^s,

the authors proposed the Landweber iteration of the form

    ξ_n+1^= ξ_n^-̣μ_n^F̣'(x_n^)̣^*J_s^(F(x_n^)̣-y^)̣ ,
    
    x_n+1^= J_q^^*(ξ_n+1^)̣ ,

for solving linear as well as nonlinear inverse problems in Banach spaces, assuming suitable smoothness and convexity onand, whereF'(x)andF'(x)^*denote the Fréchet derivative ofFatxand its adjoint,μ_n^$̣ denotes the step size,
and J_s^:→^* and J_q^^*:^*→ with 1<s,q<∞ denote the duality mappings with gauge functions
t→ t^s-1 and t→ t^q-1 respectively. This formulation of Landweber iteration, however, exclude the use of
the L^1 and the total variation like penalty functionals. A Landweber-type iteration with general uniform convex penalty functionals was
introduced in <cit.> for solving linear inverse problems and was extended in <cit.> for solving nonlinear inverse problems. Let Θ:→(-∞,∞] be a proper lower semi-continuous uniformly convex functional, the method in <cit.>
can be formulated as

    ξ_n+1^= ξ_n^-̣μ_n^F̣'(x_n^)̣^*J_s^(F(x_n^)̣-y^)̣,
    
    x_n+1^= argmin_x∈{Θ(x)-⟨ξ_n+1^,̣x⟩}.

The advantage of this method is the freedom on the choice of Θ so that it can be utilized in detecting different
features of the sought solution.

It is well known that Landweber iteration is a slowly convergent method. As alternatives to Landweber iteration, one may consider the
second order iterative methods, such as the Levenberg-Maquardt method <cit.>, the iteratively regularized Gauss-Newton
method <cit.>, or the nonstationary iterated Tikhonov regularization <cit.>. The advantage of these methods is that
they require less number of iterations to satisfy the respective stopping rule than the Landweber iteration, however they always require
to spend more computational time in dealing with each iteration step. Therefore, it becomes more desirable to accelerate Landweber
iteration by preserving its simple implementation feature.

For linear inverse problems in Hilbert spaces, a family of accelerated Landweber iterations were proposed in <cit.> using
the orthogonal polynomials and the spectral theory of self-adjoint operators. The acceleration strategy using orthogonal polynomials
is no longer available for Landweber iteration in Banach spaces with general convex penalty functionals. In <cit.> an acceleration of
Landweber iteration in Banach spaces was considered based on choosing optimal step size in each iteration step. In <cit.>
the sequential subspace optimization strategy was employed to accelerate the Landweber iteration.

The Nesterov's strategy was proposed in <cit.> to accelerate gradient method. It has played an important role on the development of
fast first order methods for solving well-posed convex optimization problems <cit.>. Recently, an accelerated version of Landweber iteration
based on Nesterov's strategy was proposed  in <cit.> which includes the following method

    z_n^= x_n^+̣n/n+α(x_n^-̣x_n-1^)̣,
    
    x_n+1^= z_n^+̣μ_n^F̣'(z_n^)̣^*(y^-̣F(z_n^)̣)

with ≥ 3 as a special case for solving ill-posed inverse problems (<ref>) in Hilbert spaces, where x_-1^=̣ x_0^=̣x_0
is an initial guess. Although no convergence analysis for (<ref>) could be given, the numerical results presented in <cit.>
clearly demonstrate its usefulness and acceleration effect. By replacing n/(n+) in (<ref>) by a general connection
parameters λ_n^$̣, a so called two-point gradient method was proposed in <cit.> and a convergence result was proved under
a suitable choice of{λ_n^}̣. Furthermore, based on the assumption of local convexity of the residual functional around the sought
solution, a weak convergence result on (<ref>) was proved in <cit.> recently.

In this paper, by incorporating an extrapolation step into the iteration scheme (<ref>), we propose a two-point gradient method
for solving inverse problems in Banach spaces with a general uniformly convex penalty termΘ, which takes the form

    ζ_n^= ξ_n^+̣λ_n^(̣ξ_n^-̣ξ_n-1^)̣,
    
    z_n^ = min_z∈{Θ(z) - łζ_n^,̣ z},
    ξ_n+1^ = ζ_n^-̣μ_n^F̣'(z_n^)̣^* J_s^ (F(z_n^)̣-y^)̣

with suitable step sizesμ_n^$̣ and combination parameters λ_n^$̣; after terminated by a discrepancy principle,
we then usex_n^:̣= min_x∈{Θ(x) - łξ_n^,̣ x}as an approximate solution.  We note that whenλ_n^=̣0, our method becomes the Landweber iteration of the form (<ref>)
and whenλ_n^=̣ n/(n+)it becomes a refined version of the Nesterov acceleration of Landweber iteration proposed in <cit.>.
We note also that, when bothandare Hilbert spaces andΘ(x)= x^2/2, our method becomes the two-point gradient
methods introduced in <cit.> for solving inverse problems in Hilbert spaces. Unlike the method in <cit.>, our method not only
works for inverse problems in Banach spaces, but also allows the use of general convex penalty functions including theL^1and the total
variation like functions.  Due to the possible non-smoothness ofΘand the non-Hilbertian structures ofand, we need to
use tools from convex analysis and subdifferential calculus to carry out the convergence analysis. Under certain conditions on the
combination parameters{λ_n^}̣, we obtain a convergence result on our method. In order to find nontrivialλ_n^$̣,
we adapt the discrete backtracking search (DBTS) algorithm in <cit.> to our situation and provide a complete convergence analysis of
the corresponding method by showing a uniform convergence result for the noise-free counterpart with the combination parameters chosen in
a certain range. Our analysis in fact improves the convergence result in <cit.> by removing the technical conditions on {λ_n^}̣
chosen by the DBTS algorithm.

The paper is organized as follows, In section 2, we give some preliminaries from convex analysis. In section 3, we
formulate our two-point gradient method with a general uniformly convex penalty term and present
the detailed convergence analysis. We also discuss the choices of the combination parameters by a discrete
backtracking search algorithm. Finally in section 4, numerical simulations are given to test the performance of the method.




§ PRELIMINARIES



In this section, we introduce some necessary concepts and properties related to Banach space and convex analysis,
we refer to <cit.> for more details.

Let  be a Banach space whose norm is denoted by ·, we use ^* to denote its dual space. Given x∈
and ξ∈^*, we write łξ, x=̊ξ (x) for the duality pairing. For a bounded linear operator A: → between
two Banach spaces  and , we use A^*: ^*→^* to denote its adjoint.

Given a convex function Θ: → (-∞, ∞], we use Θ(x) to denote the subdifferential
of Θ  at x∈, i.e.

    Θ(x) := {ξ∈^*: Θ(x̅) -Θ(x) -łξ, x̅-x≥̊0 x̅∈}.

Let (Θ): = {x∈:Θ(x)<∞} be its effective domain and let

    (Θ): = {x∈(Θ):Θ(x)≠∅}.

 The Bregman distance induced by Θ at x
in the direction ξ∈Θ(x)  is defined by

    D_ξΘ(x̅,x):=Θ(x̅)-Θ(x)-łξ, x̅-x,̊    ∀x̅∈

which is always nonnegative and satisfies the identity

    D_ξ_2Θ(x,x_2)-D_ξ_1Θ(x, x_1) =D_ξ_2Θ(x_1,x_2) +łξ_2-ξ_1, x_1-x

for all x∈(Θ), x_1, x_2∈(Θ), and ξ_1∈Θ(x_1), ξ_2∈Θ(x_2).


A proper convex function Θ: → (-∞, ∞] is called uniformly convex if there exists
a strictly increasing function h:[0,∞)→ [0,∞) with h(0) = 0 such that

    Θ(x̅+(1-)x) + (1-)h(x-x̅) ≤Θ(x̅) +(1-) Θ(x)

for all x̅,x∈ and ∈(0,1). If h(t) = c_0t^p for some c_0>0 and p≥ 2 in (<ref>),
then Θ is called p-convex. It can be shown that Θ is p-convex if and only if

    D_ξΘ(x̅,x)≥ c_0x-x̅^p,  ∀x̅∈,  x∈(Θ),  ξ∈Θ(x) .


For a proper lower semi-continuous convex function Θ: → (-∞, ∞], its Legendre-Fenchel conjugate
is defined by

    Θ^*(ξ):=sup_x∈{łξ, x-̊Θ(x)},   ξ∈^*

which is also proper, lower semi-continuous, and convex. If  is reflexive, then

    ξ∈Θ(x) ⟺ x∈Θ^*(ξ) ⟺Θ(x) +Θ^*(ξ) =łξ, x.̊

Moreover, if Θ is p-convex with p≥ 2 then it follows from <cit.> that (Θ^*)=^*,
Θ^* is Fréchet differentiable and its gradient ∇Θ^*:  X^*→ satisfies

    ∇Θ^*(ξ_1)-∇Θ^*(ξ_2) ≤(ξ_1-ξ_2/2c_0)^1/p-1,
      ∀ξ_1, ξ_2∈^*.

Consequently, it follows from (<ref>) that

    x=∇Θ^*(ξ) ⟺ξ∈Θ(x) ⟺ x =min_z∈{Θ(z) -łξ, z}.



If Θ is p-convex with p≥ 2, then for any pairs (x,ξ) and (x̅,ξ̅)
with x, x̅∈(Θ), ξ∈Θ(x), ξ̅∈Θ(x̅), we have

    D_ξΘ(x̅,x) ≤1/p^*(2c_0)^p*-1ξ-ξ̅^p* ,

where p^*: = p/(p-1).



Applying (<ref>),  x̅ = ∇Θ^*(ξ̅) and (<ref>), it follows that

    D_ξΘ(x̅,x)
        = Θ^*(ξ)-Θ^*(ξ̅)-łξ-ξ̅, ∇Θ^*(ξ̅)[0]
        = ∫_0^1 łξ -ξ̅, ∇Θ^* (ξ̅+t (ξ-ξ̅)) -∇Θ^*(ξ̅)d̊ t [0]
       ≤ξ-ξ̅∫_0^1  ∇Θ^* (ξ̅+t (ξ-ξ̅)) -∇Θ^*(ξ̅) dt
       ≤1/p^*(2c_0)^p*-1ξ-ξ̅^p*

which shows the result.


On a Banach space , we consider for 1<s<∞ the convex function x→x^s/s.
Its subgradient at x is given by

    J_s^(x):={ξ∈^*: ξ=x^s-1łξ, x=̊x^s}

which gives the duality mapping J_s^: → 2^^* of  with gauge function t→ t^s-1.
If  is uniformly smooth in the sense that its modulus of smoothness

    ρ_(t) := sup{x̅+ x+x̅-x- 2 : x̅ = 1,  x≤ t}

satisfies lim_t↘ 0ρ_(t)/t =0, then the duality mapping J_s^, for each 1<s<∞,
is single valued and uniformly continuous on bounded sets. There are many examples of uniformly smooth Banach spaces,
e.g., sequence space ℓ^s, Lebesgue space L^s, Sobolev space W^k,s and Besov space B^q,s with 1<s<∞.




§ THE TWO-POINT GRADIENT METHOD



We consider

    F(x) =y,

where F: (F)⊂→ is an operator between two Banach spaces  and . Throughout this paper,
we will always assume that  is reflexive,  is uniformly smooth, and (<ref>) has a solutions.
In order to capture the special feature of the sought solution, we will use a general convex function Θ: → (-∞, ∞]
as a penalty term. We will need a few assumptions concerning Θ and F.


Θ:→ (-∞, ∞] is proper, weakly lower semi-continuous, and p-convex with p≥ 2 such that the condition
(<ref>) is satisfied for some c_0>0.







  (a) There exists ρ > 0, x_0∈ and ξ_0∈Θ(x_0) such that B_3ρ(x_0)⊂(F)
and (<ref>) has a solution x^*∈(Θ) with

    D_ξ_0Θ(x^*, x_0) ≤ c_0 ρ^p,

where B_ρ(x_0) denotes the closed ball around x_0 with radius ρ.


  (b) The operator F is weakly closed on (F).


  (c)  There exists a family of bounded linear operators {L(x): →}_x∈ B_3ρ(x_0) such that
x→ L(x) is continuous on B_3ρ(x_0) and there is 0≤η <1 such that

    F(x)-F(x̅) -L(x̅) (x-x̅)≤ηF(x) -F(x̅)

for all x, x̅∈ B_3ρ(x_0). Moreover, there is a constant C_0>0 such that

    L(x) _→≤ C_0,   ∀ x ∈ B_3ρ(x_0).




All the conditions in Assumption 3.2 are standard. The condition (c) is called the tangential condition and is widely used in the
analysis of iterative regularization methods for nonlinear ill-posed inverse problems <cit.>. The weak closedness of F in
condition (b) means that for any sequence {x_n} in (F) satisfying x_n⇀ x and F(x_n)⇀ v,
then x∈(F) and F(x) = v, where we use “⇀" to denote the weak convergence.

Using the p-convex function Θ specified in Assumption <ref>, we may pick among solutions of (<ref>) the one
with the desired feature. We define x^† to be a solution of (<ref>) with the property

    D_ξ_0Θ(x^†,x_0) := min_x∈(Θ)∩(F) {D_ξ_0Θ(x,x_0) : F(x) = y } .

When  is reflexive, by using the weak closedness of F and the weak lower semi-continuity of Θ, it is standard
to show that x^† exists. According to Assumption <ref> (a), we always have

    D_ξ_0Θ(x^†, x_0) ≤ c_0 ρ^p

which together with Assumption <ref> implies that x^†-x_0≤ρ. The following lemma shows that x^† is
uniquely defined.

 
Under Assumption <ref> and Assumption <ref>, the solution x^† of (<ref>) satisfying (<ref>) is uniquely
defined.



This is essentially proved in <cit.>.


In order to construct an approximate solution to (<ref>), we will formulate a two-point gradient method with
penalty term induced by the uniformly convex function Θ. Let τ>1 be a given number. By picking
x_-1^=̣ x_0^:̣=x_0∈ and ξ_-1^=̣ξ_0^:̣=ξ_0∈∂Θ(x_0) as initial guess, for n≥ 0 we define

    ζ_n^= ξ_n^+̣λ_n^(̣ξ_n^-̣ξ_n-1^)̣,
    
    z_n^= ∇Θ^*(ζ_n^)̣,
    ξ_n+1^=ζ_n^-̣μ_n^Ḷ(z_n^)̣^* J_s^(r_n^)̣,
    
    x_n+1^= ∇Θ^*(ξ_n+1^)̣,

where r_n^=̣ F(z_n^)̣-y^$̣,λ_n^≥̣0is the combination parameter,μ_n^$̣ is the step sizes defined by

    μ_n^=̣{[ min{μ̅_0 r_n^^p(s-1)/L(z_n^)̣^*J_s^(r_n^)̣^p, μ̅_1}r_n^^p-s                                               r_n^>τ[2ex]
    0                                                      r_n^≤τ ].

for some positive constants μ̅_0 and μ̅_1, and the mapping J_s^: →^* with 1<s<∞ denotes the
duality mapping of  with gauge function t→ t^s-1, which is single-valued and continuous because
 is assumed to be uniformly smooth. We remark that when λ_n^=̣0 for all n, the method (<ref>) reduces to the
Landweber iteration considered in <cit.> where a detailed convergence analysis has been carried out. When λ_n^=̣ n/(n+)
with ≥ 3 for all n, the method (<ref>) becomes a refined version of the Nesterov acceleration of Landweber iteration
proposed in <cit.>; although there is no convergence theory available, numerical simulations in <cit.> demonstrate its
usefulness and acceleration effect. In this paper we will consider (<ref>) with λ_n^$̣ satisfying suitable conditions
to be specified later. Note that our method (<ref>) requires the use of the previous two iterations at every iteration step,
which follows the spirit from <cit.>; on the other hand, our method allows the use of a generalp-convex penalty
functionΘ, which could be non-smooth, to reconstruct solutions with special features such as sparsity and discontinuities.




 §.§ Convergence


In order to use our TPG-Θmethod (<ref>) to produce a useful approximate solution to (<ref>),
the iteration must be terminated properly. We will use the discrepancy principle with respect toz_n^$̣,
i.e., for a given τ>1, we will terminate the iteration after n_$̣ steps, wheren_:̣= n(,̣y^)̣is the integer such that

    F(z_n_^)̣-y^≤τ<̣F(z_n^)̣-y^ ,   0≤ n<n_,̣

and usex_n_^$̣ as an approximated solution. To carry out the convergence analysis of x_n_^$̣ tox^†as→̣0, we are going to show that, for any solutionx̂of (<ref>) inB_2ρ(x_0)∩(Θ), the Bregman distanceD_ξ_n^Θ(x̂, x_n^)̣,0≤n ≤n_$̣, is
monotonically decreasing with respect to n. To this end, we introduce

    _n: = D_ξ_n^Θ(x̂, x_n^)̣ - D_ξ_n-1^Θ(x̂, x_n-1^)̣.

We will show that, under suitable choice of {λ_n^}̣, there holds Δ_n≤ 0 for 0≤ n ≤ n_$̣.


Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Then, for any solution x̂ of (<ref>) in  B_2ρ(x_0)⋂(Θ), there holds

    D_ζ_n^Θ(x̂, z_n^)̣ -D_ξ_n^Θ(x̂, x_n^)̣
               ≤λ_n^_n+1/p^*(2c_0)^p^*-1(λ_n^+̣(λ_n^)^p^*)
    ξ_n^-̣ξ_n-1^^p^*.

If z_n^∈̣B_3ρ(x_0) then

    D_ξ_n+1^Θ(x̂,x_n+1^)̣ - D_ζ_n^Θ(x̂,z_n^)̣   ≤ -(1-η -1/p^*(μ̅_0/2 c_0)^p^*-1) μ_n^F(z_n^)̣-y^^s
            +(1+η)μ_n^F(z_n^)̣-y^^s-1.̣

By using the identity (<ref>), Lemma <ref> and the definition of ζ_n^$̣, we can obtain

    D_ζ_n^Θ(x̂, z_n^)̣ -D_ξ_n^Θ(x̂, x_n^)̣   = D_ζ_n^Θ(x_n^,̣ z_n^)̣+łζ_n^-̣ξ_n^,̣x_n^-̣x̂   ≤1/p^*(2c_0)^p^*-1ζ_n^-̣ξ_n^^p^*+ łζ_n^-̣ξ_n^,̣x_n^-̣x̂   = 1/p^*(2c_0)^p^*-1(λ_n^)̣^p^*ξ_n-1^-̣ξ_n^^p^*+ łζ_n^-̣ξ_n^,̣x_n^-̣x̂.

By using again the definition ofζ_n^$̣, (<ref>) and Lemma <ref>, and referring to
the definition of Δ_n and λ_n^≥̣0, we have

    łζ_n^-̣ξ_n^,̣ x_n^-̣x̂=λ_n^ł̣ξ_n^-̣ξ_n-1^,̣x_n^-̣x̂
       =λ_n^(D_ξ_n^Θ(x̂, x_n^)̣-D_ξ_n-1^Θ(x̂, x_n-1^)̣+D_ξ_n-1^Θ(x_n^,̣x_n-1^)̣)
       ≤λ_n^_n + 1/p^*(2c_0)^p^*-1λ_n^ξ_n^-̣ξ_n-1^^p^* .

The combination of the above two estimates yields (<ref>).

To derive (<ref>), we first use the identity (<ref>) and Lemma <ref> to obtain

    D_ξ_n+1^Θ(x̂,x_n+1^)̣ - D_ζ_n^Θ(x̂,z_n^)̣
            = D_ξ_n+1^Θ(z_n^,̣x_n+1^)̣+łξ_n+1^-̣ζ_n^,̣z_n^-̣x̂
           ≤1/p^*(2c_0)^p^*-1ξ_n+1^-̣ζ_n^^p^* + łξ_n+1^-̣ζ_n^,̣z_n^-̣x̂.̊

Recall the definition of ξ_n+1^$̣ in (<ref>), we have

    ξ_n+1^-̣ζ_n^^p^*   = (μ_n^)̣^p^*L(z_n^)̣^*J_s^𝒴(F(z_n^)̣-y^)̣^p^*.

According to the definition (<ref>) ofμ_n^$̣, one can see that
μ_n^≤̣μ̅_0 F(z_n^)̣-y^^s(p-1)/L(z_n^)̣^* J_s^ (F(z_n^)̣-y^)̣^p,

which implies that

    (μ_n^)̣^p^*-1L(z_n^)̣^*J_s^𝒴(F(z_n^)̣-y^)̣^p^*   ≤μ̅_0^p^*-1F(z_n^)̣-y^^s(p-1)(p^*-1)
        = μ̅_0^p^*-1F(z_n^)̣-y^^s.

Therefore, the first term on the right hand side of (<ref>) can be estimated as

    1/p^*(2c_0)^p^*-1ξ_n+1^-̣ζ_n^^p^*≤1/p^*(μ̅_0/2c_0)^p^*-1μ_n^F(z_n^)̣-y^^s.

For the second term  on the right hand side of (<ref>), we may use the definition of ξ_n+1^$̣,
the property ofJ_s^, and the definition ofμ_n^$̣ to derive that

    łξ_n+1^-̣ζ_n^,̣z_n^-̣x̂ = -μ_n^ł̣J_s^𝒴(F(z_n^)̣-y^)̣,L(z_n^)̣(z_n^-̣x̂)    = -μ_n^ł̣J_s^𝒴(F(z_n^)̣-y^)̣, y^-̣F(z_n^)̣-L(z_n^)̣(x̂-z_n^)̣       - μ_n^F(z_n^)̣-y^^s
       ≤μ_n^F(z_n^)̣-y^^s-1(δ+y-F(z_n^)̣-L(z_n^)̣ (x̂-z_n^)̣) 
           -μ_n^F(z_n^)̣-y^^s.

Recall that z_n^∈̣B_3ρ(x_0), we may use Assumption <ref> (c) to further obtain

    łξ_n+1^-̣ζ_n^,̣z_n^-̣x̂
       ≤μ_n^F(z_n^)̣-y^^s-1(+̣ηF(z_n^)̣-y) -μ_n^F(z_n^)̣-y^^s 
       ≤ (1+η) μ_n^F(z_n^)̣-y^^s-1-̣ (1-η) μ_n^F(z_n^)̣-y^^s.

The combination of above two estimates (<ref>) and (<ref>) with (<ref>) yields (<ref>).

Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Assume that τ>1 and μ̅_0>0 are chosen such that

    c_1:=1-η-1+η/τ -1/p^*(μ̅_0/2 c_0)^p^*-1>0.

If z_n^∈̣B_3ρ(x_0), then, for any solution x̂ of (<ref>) in B_2ρ(x_0)⋂(Θ), there holds

    _n+1   ≤λ_n^_n+1/p^*(2c_0)^p^*-1(λ_n^+̣(λ_n^)̣^p^*) ξ_n^-̣ξ_n-1^^p^*
           -c_1 μ_n^F(z_n^)̣-y^^s,

where Δ_n is defined by (<ref>).

By using the definition of μ_n^$̣ it is easily seen thatμ_n^̣̣≤μ_n^F(z_n^)̣-y^/τ.
It then follows from (<ref>) that

    D_ξ_n+1^Θ(x̂, x_n+1^)̣ - D_ζ_n^Θ(x̂, z_n^)̣≤ -c_1 μ_n^F(z_n^)̣-y^^s.

Combining this estimate with (<ref>) yields (<ref>).



We will use Lemma <ref> and Lemma <ref> to show thatz_n^∈̣B_3ρ(x_0)andΔ_n ≤ 0for alln≥ 0and that the integern_$̣ determined by (<ref>) is finite. To this end, we need to place
conditions on {λ_n^}̣. We assume that {λ_n^}̣ is chosen such that

    1/p^*(2c_0)^p^*-1((λ_n^)̣^p^*+λ_n^)ξ_n^-̣ξ_n-1^^p^*≤ c_0ρ^p

and

    1/p^*(2c_0)^p^*-1((λ_n^)̣^p^*+λ_n^) ξ_n^-̣ξ_n-1^^p^*
    -c_1/νμ_n^F(z_n^)̣-y^^s≤ 0

for all n≥ 0, where ν>1 is a constant independent of $̣ andn. We will discuss how to choose{λ_n^}̣to satisfy (<ref>) and (<ref>) shortly.


Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Let τ>1 and μ̅_0>0 be chosen such that (<ref>) holds.  If {λ_n^}̣ is chosen such that (<ref>) and (<ref>) hold, then

    z_n^∈̣B_3ρ(x_0)      x_n^∈̣B_2ρ(x_0)    n≥ 0.

Moreover, for any solution x̂ of (<ref>) in B_2ρ(x_0)⋂(Θ) there hold

    D_ξ_n^Θ(x̂, x_n^)̣≤ D_ξ_n-1^Θ(x̂, x_n-1^)̣

and

    ∑_m=0^nμ_m^F(z_m^)̣-y^^s ≤ν/(ν-1)c_1 D_ξ_0Θ(x̂, x_0)

for all n≥ 0.  Let n_$̣ be chosen by the discrepancy principle (<ref>), thenn_$̣ must be a finite integer.

We will show (<ref>) and (<ref>) by induction. Since x_-1^=̣ x_0^=̣ x_0, ξ_-1^=̣ξ_0^=̣ξ_0,
and z_0^=̣∇Θ^*(ζ_0^)̣ = ∇Θ^*(ξ_0) = x_0, they are trivial for n=0. Now we assume that
(<ref>) and (<ref>) hold for all 0≤ n≤ m for some integer m≥ 0, we will show that they are also true
for n = m+1. By the induction hypotheses z_m^∈̣B_3ρ(x_0), we may use Lemma <ref> and
(<ref>) to derive that

    Δ_m+1≤λ_m^Δ̣_m - (1-1/ν) c_1μ_m^F(z_m^)̣-y^^s.

Since λ_m^≥̣0 and ν>1, this together with the induction hypothesis Δ_m ≤ 0 implies that

    Δ_m+1≤ - (1-1/ν) c_1 μ_m^F(z_m^)̣-y^^s ≤ 0

which shows (<ref>) for n=m+1. Consequently, by taking x̂= x^† and using Assumption <ref> (a), we have

D_ξ_m+1^Θ(x^†, x_m+1^)̣≤ D_ξ_m^Θ(x^†, x_m^)̣≤⋯≤ D_ξ_0Θ(x^†, x_0)≤ c_0 ρ^p.

By virtue of Assumption <ref>, we then have c_0x_m+1^-̣x^†^p≤ c_0 ρ^p which together with x^†∈ B_ρ(x_0)
implies that x_m+1^∈̣B_2ρ(x_0). Now we may use (<ref>) in Lemma <ref>, (<ref>) and Δ_m+1≤ 0
to derive that

    D_ζ_m+1^Θ(x^†, z_m+1^)̣   ≤ D_ξ_m+1^Θ(x^†, x_m+1^)̣+ λ_m+1^_m+1 +c_0ρ^p 
       ≤ D_ξ_m+1^Θ(x^†, x_m+1^)̣ +c_0ρ^p
       ≤ D_ξ_0Θ(x^†, x_0)+c_0ρ^p 
       ≤ 2c_0ρ^p .

This together with Assumption <ref> yields x^†-z_m+1^≤ 2^1/pρ≤ 2ρ, and consequently
z_m+1^∈̣B_3ρ(x_0). We therefore complete the proof of (<ref>) and (<ref>).

Since (<ref>) and (<ref>) are valid, the inequality (<ref>) holds for all m≥ 0.
Thus

    (1-1/ν) c_1 μ_m^F(z_m^)̣-y^^s
    ≤ D_ξ_m^Θ(x̂, x_m^)̣ - D_ξ_m+1^Θ(x̂, x_m+1^)̣

for m≥ 0. Hence, for any integer n ≥ 0 we have

    (1-1/ν)c_1∑_m=0^n μ_m^y^-̣F(z_m^)̣^s
       ≤ D_ξ_0Θ(x̂, x_0)- D_ξ_n+1^Θ(x̂, x_n+1^)̣
       ≤ D_ξ_0Θ(x̂, x_0)

which shows (<ref>).

If n_$̣ is not finite, thenF(z_m^)̣-y^>τ$̣ for all integers m and consequently, by using
L(x)≤ C_0 from Assumption <ref> (c) and the property of J_s^, we have

    μ_m^= min{μ̅_0 F(z_m^)̣-y^^p(s-1)/L(z_m^)̣^* J_s^(F(z_m^)̣-y^)̣^p, μ̅_1}F(z_m^)̣-y^^p-s
       ≥min{μ̅_0/C_0^p, μ̅_1}F(z_m^)̣-y^^p-s.

Therefore, it follows from (<ref>) that

    ν/(ν-1)c_1 D_ξ_0Θ(x̂, x_0)
       ≥min{μ̅_0/C_0^p, μ̅_1}∑_m=0^n F(z_m^)̣-y^^p
       ≥min{μ̅_0/C_0^p, μ̅_1} (n+1) τ^p ^̣p

for all n≥ 0. By taking n→∞ we derive a contradiction. Therefore, n_$̣ must be finite.
 In the proof of Proposition <ref>, the condition (<ref>) plays a crucial role. Note that, by the
definition of our method (<ref>), z_n^$̣ depends onλ_n^$̣. Therefore, it is not immediately clear how to choose
λ_n^$̣ to make (<ref>) satisfied. One may ask if there existsλ_n^$̣ such that (<ref>) holds.
Obviously λ_n^=̣0 satisfies the inequality, which correspond to the Landweber iteration. In order to achieve acceleration,
it is necessary to find nontrivial λ_n^$̣. Note that whenF(z_n^)̣-y^≤τ$̣ occurs, (<ref>)
forces λ_n^=̣0 because μ_n^=̣0. Therefore we only need to consider the case F(z_n^)̣-y^ >τ$̣. By using
(<ref>) we can derive a sufficient condition

    1/p^*(2c_0)^p^*-1(λ_n^+̣(λ_n^)̣^p^*)ξ_n^-̣ξ_n-1^^p^*≤ M τ^p ^̣p,

whereM := c_1/νmin{μ̅_0/C_0^p, μ̅_1}.Considering the particular case whenp=2, this thus leads to the choice

    λ_n^:̣ =min{-1/2+√(1/4+4c_0 Mτ^2δ^2/ξ_n^-̣ξ_n-1^^2), n/n+} ,

where≥ 3is a given number. Note that in the above formula forλ_n^$̣, inside the “min" the second argument is
taken to be n/(n+) which is the combination parameter used in Nesterov's acceleration strategy; in case the first argument is large,
this formula may lead to λ_n^=̣ n/(n+) and consequently the acceleration effect of Nesterov can be utilized.
For general p>1, by placing the requirement 0≤λ_n^≤̣n/(n+) ≤ 1, one may choose λ_n^$̣ to satisfy

    2 λ_n^/p^*(2 c_0)^p^*-1ξ_n^-̣ξ_n-1^^p^*≤ Mτ^p ^̣p

which leads to the choice

    λ_n^=̣min{γ_0 ^̣p/ξ_n^-̣ξ_n-1^^p^*, n/n+},
      γ_0 := 1/2 (2 c_0)^p^*-1 p^* Mτ^p.

We remark that the choices ofλ_n^$̣ given in (<ref>) and (<ref>) may decrease to 0 as →̣0,
consequently the acceleration effect could also decrease for →̣0. Since for small values of $̣
the acceleration is needed most, other strategies should be explored. We will give a further consideration on the choice ofλ_n^$̣ in the next subsection.


In order to establish the regularization property of the method (<ref>), we need to consider its noise-free counterpart.
By dropping the superscript$̣ in all the quantities involved in (<ref>), it leads to the following formulation
of the two-point gradient method for the noise-free case:

    ζ_n    = ξ_n +λ_n (ξ_n-ξ_n-1),
    
    z_n    = ∇Θ^*(ζ_n),
    ξ_n+1   = ζ_n-μ_n L(z_n)^* J_s^(r_n),
    
    x_n+1   = ∇Θ^*(ξ_n+1),

where r_n = F(z_n)-y, λ_n≥ 0 is the combination parameter, and μ_n is the step size given by

    μ_n = {[ min{μ̅_0 r_n^p(s-1)/L(z_n)^*J_s^(r_n)^p, μ̅_1}r_n^p-s                                            F(z_n)  y,;                                                     0                                           F(z_n) = y. ].

We will first establish a convergence result for (<ref>). The following result plays a crucial role
in the argument.


Consider the equation (<ref>) for which Assumption <ref> holds.
Let Θ: → (-∞, ∞] be a proper, lower semi-continuous and uniformly convex function.
Let {x_n}⊂ B_2ρ(x_0) and {ξ_n}⊂^* be such that


  (i)ξ_n∈Θ(x_n) for all n;


  (ii) for any solution x̂ of (<ref>) in B_2ρ(x_0)∩(Θ) the sequence
{D_ξ_nΘ(x̂, x_n)} is monotonically decreasing;


  (iii)lim_n→∞F(x_n)-y=0.


  (iv) there is a subsequence {n_k} with n_k→∞ such that for any solution
x̂ of (<ref>) in B_2ρ(x_0)∩(Θ) there holds

    lim_l→∞sup_k≥ l |łξ_n_k-ξ_n_l, x_n_k-x̂|̊ = 0.

Then there exists a solution x_* of (<ref>) in B_2ρ(x_0)∩(Θ) such that

    lim_n→∞ D_ξ_nΘ(x_*, x_n)=0.

If, in addition, x^†∈ B_ρ(x_0)∩(Θ) and ξ_n+1-ξ_n ∈ℛ(L(x^†)^*)
for all n, then x_*=x^†.

This result essentially follows from <cit.> and its proof.

Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Assume that μ̅_0>0 is chosen such that

1-η - 1/p^*(μ̅_0/2 c_0)^p^*-1 >0

and the combination parameters {λ_n} are chosen to satisfy the counterparts of (<ref>)
and (<ref>) with =̣ 0  and

    ∑_n=0^∞λ_nξ_n-ξ_n-1<∞.

Then, there exists a solution x_* of (<ref>) in B_2ρ(x_0)∩(Θ) such that

    lim_n→∞x_n-x_*=0    lim_n→∞ D_ξ_nΘ(x_*, x_n)=0.

If in addition (L(x^†))⊂(L(x)) for all x∈ B_3ρ(x_0), then x_*=x^†.
 We will use Proposition <ref> to prove the result. By the definition x_n = ∇Θ^*(ξ_n)
we have ξ_n ∈Θ(x_n) which shows (i) in Proposition <ref>. By using the same argument for
proving Proposition <ref> we can show that z_n ∈ B_3ρ(x_0) and x_n ∈ B_2ρ(x_0) for all n with

    D_ξ_n+1Θ(x̂, x_n+1) ≤ D_ξ_nΘ(x̂, x_n)

for any solution x̂ of (<ref>) in B_2ρ(x_0) ∩(Θ) and

    ∑_n=0^∞μ_n F(z_n)-y^s <∞.

From (<ref>) it follows that (ii) in Proposition <ref> holds. Moreover, by using the definition of
μ_n and the similar derivation for (<ref>) we have
min{μ̅_0/C_0^p, μ̅_1}F(z_n)-y^p
≤μ_n F(z_n)-y^s ≤μ̅_1 F(z_n)-y^p.

Thus it follows from (<ref>) that
∑_n=0^∞F(z_n)-y^p <∞.

Consequently

    lim_n→∞F(z_n)-y= 0.

By using Assumption <ref> (c), (<ref>), (<ref>) and (<ref>) with =̣0, we have

    F(x_n)-F(z_n)   ≤1/1-ηL(z_n) (x_n-z_n)≤C_0/1-ηx_n-z_n[0]
        = C_0/1-η∇Θ^*(ξ_n)-∇Θ^*(ζ_n)
       ≤C_0/(1-η)(2c_0)^p^*-1ξ_n-ζ_n^p^*-1[0]
        = C_0/(1-η)(2c_0)^p^*-1λ_n^p^*-1ξ_n-ξ_n-1^p^*-1[0]
       ≤C_0/1-η(c_1 p^*/2 c_0 ν)^1/p (μ_nF(z_n)-y^s)^1/p[0]
       ≤C_0/1-η(c_1 p^*μ̅_1/2 c_0 ν)^1/pF(z_n)-y.

The combination of (<ref>) and (<ref>) implies that F(x_n)-y→ 0 as n→∞
which shows (iii) in Proposition <ref>.

In order to establish the convergence result, it remains only to show (iv) in Proposition <ref>.
To this end, we consider F(z_n)-y. It is known that F(z_n)-y→ 0 as n→∞.
If F(z_n)-y =0 for some n, then (<ref>) with =̣0 forces λ_n (ξ_n-ξ_n-1) =0.
Thus ζ_n = ξ_n. On the other hand, we also have μ_n=0 and hence ξ_n+1 = ζ_n. Consequently
ξ_n+1 = ζ_n =ξ_n and
ζ_n+1 = ξ_n+1 + λ_n+1(ξ_n+1-ξ_n) = ξ_n+1 = ζ_n.

Thus

z_n+1 = ∇Θ^*(ζ_n+1) = ∇Θ^*(ζ_n) = z_n

which implies that F(z_n+1)= F(z_n) =y. By repeating the argument one can see
that F(z_m)=y for all m≥ n. Therefore we can choose a strictly increasing sequence
{n_k} of integers by letting n_0=0 and for each k≥ 1, letting n_k be the first integer satisfying

    n_k≥ n_k-1+1   and   F(z_n_k)-y≤F(z_n_k-1)-y.

For such chosen strictly increasing sequence {n_k}, it is easily seen that

    F(z_n_k)-y≤F(z_n)-y,   0≤ n<n_k.

For any integers 0≤ l<k<∞, we consider
łξ_n_k-ξ_n_l,x_n_k-x̂=̊∑_n=n_l^n_k-1łξ_n+1-ξ_n, x_n_k -x̂.̊
By using the definition of ξ_n+1 we have
ξ_n+1-ξ_n = λ_n (ξ_n-ξ_n-1) -μ_n L(z_n)^* J_s^(F(z_n)-y).

Therefore, by using the property of J_s^, we have

    |łξ_n_k-ξ_n_l,x_n_k-x̂|
       ≤∑_n=n_l^n_k-1λ_n |łξ_n-ξ_n-1, x_n_k-x̂|̊
           + ∑_n=n_l^n_k-1μ_n |ł J_s^ (F(z_n)-y), L(z_n) (x_n_k-x̂)|̊[0]
       ≤∑_n=n_l^n_k-1λ_nξ_n-ξ_n-1x_n_k-x̂
           + ∑_n=n_l^n_k-1μ_nF(z_n)-y^s-1L(z_n)(x_n_k-x̂).

By using Assumption <ref> (c) and (<ref>), we obtain for n<n_k that

    L(z_n)(x_n_k-x̂)   ≤L(z_n)(x_n_k-z_n)+L(z_n)(z_n-x̂)
       ≤ (1+η)(F(x_n_k)-F(z_n) + F(z_n)-y)
       ≤ 2(1+η)F(z_n)-y + (1+η)F(x_n_k)-y
       ≤ 2(1+η)F(z_n)-y
           + (1+η)(F(x_n_k)-F(z_n_k)+F(z_n_k)-y)
       ≤ 3(1+η) F(z_n)-y + (1+η) F(x_n_k)-F(z_n_k).

By using (<ref>) and (<ref>), we have for n<n_k that

    F(x_n_k)-F(z_n_k)   ≤C_0/1-η(c_1 p^*μ̅_1/2 c_0 ν)^1/pF(z_n_k)-y
       ≤C_0/1-η(c_1 p^*μ̅_1/2 c_0 ν)^1/pF(z_n)-y.

Therefore

    L(z_n)(x_n_k-x̂)≤ C_1 F(z_n)-y

for n<n_k, where C_1 := 3(1+η) + (1+η)C_0/1-η(c_1 p^*μ̅_1/2 c_0 ν)^1/p.
Combining (<ref>)
with (<ref>) and using x_n_k∈ B_2ρ(x_0) we obtain

    |łξ_n_k-ξ_n_l,x_n_k-x̂|
       ≤∑_n=n_l^n_k-1λ_nξ_n-ξ_n-1x_n_k-x̂+ C_1 ∑_n=n_l^n_k-1μ_n F(z_n)-y^s
       ≤ 4ρ∑_n=n_l^n_k-1λ_nξ_n-ξ_n-1+ C_1 ∑_n=n_l^n_k-1μ_n F(z_n)-y^s.

By making use of (<ref>), we obtain, with C_2:= ν C_1/((ν-1) c_1),

    |łξ_n_k-ξ_n_l,x_n_k-x̂| 
       ≤ 4 ρ∑_n=n_l^n_k-1λ_nξ_n-ξ_n-1
    + C_2 ∑_n=n_l^n_k-1(D_ξ_nΘ(x̂,x_n)-D_ξ_n+1Θ(x̂,x_n+1)) [0]
        = 4 ρ∑_n=n_l^n_k-1λ_nξ_n-ξ_n-1
    +C_2 (D_ξ_n_lΘ(x̂,x_n_l)-D_ξ_n_kΘ(x̂,x_n_k)).

Let γ := lim_n→∞ D_ξ_nΘ(x̂, x_n) whose existence is guaranteed by the monotonicity
of {D_ξ_nΘ(x̂, x_n)}. Then

    sup_k≥ l|łξ_n_k-ξ_n_l,x_n_k-x̂|
       ≤ 4 ρ∑_n=n_l^∞λ_nξ_n-ξ_n-1 + C_2 (D_ξ_n_lΘ(x̂,x_n_l)-γ).

Thus it follows from (<ref>) that

    lim_l→∞sup_k≥ l|łξ_n_k-ξ_n_l,x_n_k-x̂|
    ≤ C_2 (lim_l→∞ D_ξ_n_lΘ(x̂,x_n_l)-γ)=0

which verifies (iv) in Proposition <ref>.

To show x_* = x^† under the additional condition (L(x^†))⊂(L(x)) for
all x∈ B_3ρ(x_0), we observe from (<ref>) and ξ_0-ξ_-1 =0 that

    ξ_n+1-ξ_n    = -μ_n L(z_n)^*(F(z_n)-y) + λ_n(ξ_n-ξ_n-1)
        = -∑_k=0^n(∏_i=k+1^nλ_i )μ_k L(z_k)^*(F(z_k)-y).

Since  is reflexive and (L(x^†))⊂(L(x)), we have (L(x)^*)⊂(L(x^†)^*) for all x∈ B_3ρ(x_0). Recall that z_k∈ B_3ρ(x_0). It thus
follows from the above formula that ξ_n+1-ξ_n ∈ℛ(L(x^†)^*). Therefore
we may use the second part of Proposition <ref> to conclude the proof.


Next, we are going to show that, using the discrepancy principle (<ref>) as a stopping rule, our method (<ref>)
becomes a convergent regularization method, if we additionally assume that λ_n^$̣ depends continuously on$̣
in the sense that λ_n^→̣λ_n as →̣0 for all n. We need the following stability result.


Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Assume that τ>1 and μ̅_0>0 are chosen to satisfy (<ref>).  Assume also that the combination parameters
{λ_n^}̣ are chosen to depend continuously on $̣ as→̣0and satisfy (<ref>),
(<ref>) and (<ref>). Then for alln≥ 0there hold

    ζ_n^→̣ζ_n,    z_n^→̣z_n,   ξ_n^→̣ξ_n     x_n^→̣x_n  as →̣0.

The result is trivial for n=0. We next assume that the result is true for all 0≤ n≤ m and show that
the result is also true for n = m+1. We consider two cases.

 Case 1: F(z_m)=y. In this case we have μ_m=0 and F(z_m^)̣ -y^→ 0
as →̣0 by the continuity of F and the induction hypothesis z_m^→̣z_m. Thus

    ξ_m+1^-̣ξ_m+1 =ζ_m^-̣ζ_m -μ_m^Ḷ(z_m^)̣^* J_s^ (F(z_m^)̣-y^)̣,

which together with the definition of μ_m^$̣ and the induction hypothesisζ_m^→̣ζ_mimplies that

    ξ_m+1^-̣ξ_m+1≤ζ_m^-̣ζ_m + C_0 μ̅_1 F(z_m^)̣-y^^p-1→ 0
      →̣0.

Consequently, by using the continuity of∇Θ^*we havex_m+1^=̣∇Θ^*(ξ_m+1^)̣→∇Θ^*(ξ_m+1) =x_m+1as→̣0. Recall thatζ_m+1^=̣ξ_m+1^+̣λ_m+1^(̣ξ_m+1^-̣ξ_m^)̣,     
z_m+1^=̣∇Θ^*(ζ_m+1^)̣.We may use the conditionλ_m+1^→̣λ_m+1to conclude thatζ_m+1^→̣ζ_m+1andz_m+1^→̣z_m+1as→̣0.

 Case 2: F(z_m) y. In this case we haveF(z_m^)̣-y^≥τ$̣ for small >̣0. Therefore

    μ_m^= min{μ̅_0 F(z_m^)̣-y^^p(s-1)/L(z_m^)̣^*J_s^ (F(z_m^)̣-y^)̣^p, μ̅_1}F(z_m^)̣-y^^p-s,
    μ_m    = min{μ̅_0 F(z_m)-y^p(s-1)/L(z_m)^*J_s^ (F(z_m)-y)^p, μ̅_1}F(z_m)-y^p-s.

If L(z_m)^*J_s^(F(z_m)-y) 0, then, by the induction hypothesis on z_m^$̣, it is easily seen thatμ_m^→̣μ_mas→̣0.
IfL(z_m)^*J_s^(F(z_m)-y) = 0, thenμ_m= μ̅_1 F(z_m)-y^p-sandμ_m^=̣μ̅_1 F(z_m^)̣-y^^p-sfor small>̣0. This again implies thatμ_m^→̣μ_mas→̣0.
Consequently, by utilizing the continuity ofF,L,J_s^and∇Θ^*and the induction hypotheses,
we can conclude thatξ_m+1^→̣ξ_m+1,x_m+1^→̣x_m+1,ζ_m+1^→̣ζ_m+1andz_m+1^→̣z_m+1as→̣0.


We are now in a position to give the main convergence result on our method (<ref>).


Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Assume that τ>1 and μ̅_0>0 are chosen to satisfy (<ref>).  Assume also that the combination parameters
{λ_n^}̣ are chosen to depend continuously on $̣ as→̣0and satisfy (<ref>),
(<ref>) and (<ref>). Letn_$̣ be chosen according to the discrepancy principle (<ref>).
Then there exists a solution x_* of (<ref>) in B_2ρ(x_0)∩(Θ) such that

    lim_→̣0x_n_^-̣x_*=0         lim_→̣0 D_ξ_n_^Θ(x_*, x_n_^)̣ =0.

If in addition (L(x^†))⊂(L(x)) for all x∈ B_3ρ(x_0), then x_*=x^†.


Let x_* be the solution of (<ref>) determined in Theorem <ref>.
We complete the proof by considering two cases.

 Case 1: Assume first that {y^_̣l} is a sequence satisfying y^_̣l-y≤_̣l
with _̣l→ 0 such that n_l:=n__̣l→n̂ as l→∞ for some finite
integer n̂. We may assume n_l = n̂ for all l. According to Lemma <ref>, we have
ζ_n_l^_̣l→ζ_n̂,    z_n_l^_̣l→ z_n̂,   ξ_n_l^_̣l→ξ_n̂     x_n_l^_̣l→ x_n̂   l→∞.

By the definition of n_l: = n__̣l, we have

    F(z_n_l^_̣l)-y^_̣l≤τ_̣l.

By using the similar argument for deriving (<ref>) we have
F(x_n_l^_̣l)-F(z_n_l^_̣l)≤ C F(z_n_l^_̣l)-y^_̣l
for some universal constant C. Thus
F(x_n_l^_̣l)-y^_̣l≤ (1+C) F(z_n_l^_̣l)-y^_̣l≤ (1+C) τ_̣l.

Taking l→∞ and using the continuity of F gives F(x_n̂) =y. Thus x_n̂ is a solution
of (<ref>) in B_2ρ(x_0)∩(Θ). By the monotonicity of {D_ξ_nΘ(x_n̂, x_n)}
with respect to n, we then obtain

D_ξ_nΘ(x_n̂, x_n) ≤ D_ξ_n̂Θ(x_n̂, x_n̂)=0,   ∀ n≥n̂.

Therefore x_n = x_n̂ for all n≥n̂. Since Theorem <ref> shows that x_n → x^* as n→∞,
we must have x_n̂ = x^* and thus x_n_l^_̣l→ x^* as l→∞.
This together with the lower semi-continuity of Θ shows that

    0    ≤lim inf_l→∞ D_ξ_n_l^_̣lΘ(x_*, x_n_l^_̣l)
    ≤lim sup_k→∞ D_ξ_n_l^_̣lΘ(x_*, x_n_l^_̣l)
       ≤Θ(x_*) -lim inf_l→∞Θ(x_n_l^_̣l)
    - lim_l→∞łξ_n_l^_̣l, x_*-x_n_l^_̣l   ≤Θ(x_*)-Θ(x_*)=0,.

Therefore  lim_l→∞ D_ξ_n_l^_̣lΘ(x_*, x_n_l^_̣l)=0.

 Case 2: Assume {y^_̣l} is a sequence satisfying y^_̣l-y≤_̣l
with _̣l→ 0 such that n_l:=n__̣l→∞ as l→∞.
Let n be any fixed integer, then n_l>n for large l. It then follows from Lemma <ref> that

    D_ξ_n_l^_̣lΘ(x_*, x_n_l^_̣l) ≤ D_ξ_n^_̣lΘ (x_*, x_n^_̣l)
    =Θ(x_*) -Θ(x_n^_̣l) -łξ_n^_̣l, x_*-x_n^_̣l.̊

By using Lemma <ref> and the lower semi-continuity of Θ we obtain

    0    ≤lim inf_l→∞ D_ξ_n_l^_̣lΘ (x_*, x_n_l^_̣l)
    ≤lim sup_l→∞ D_ξ_n_l^_̣lΘ (x_*, x_n_l^_̣l)
       ≤Θ(x_*) -lim inf_l→∞Θ(x_n^_̣l) -lim_l→∞łξ_n^_̣l, x_*-x_n^_̣l   ≤Θ(x_*) -Θ(x_n) -łξ_n, x_*-x_n   = D_ξ_nΘ(x_*, x_n) .

Since n can be arbitrary and Theorem <ref> implies that D_ξ_nΘ(x_*, x_n)→ 0
as n→∞, we therefore have lim_l→∞ D_ξ_n_l^_̣lΘ (x_*, x_n_l^_̣l)=0
and hence lim_l→∞x_n_l^_̣l - x_* =0.


 §.§ DBTS: the choice of λ_n^$̣

In this section we will discuss the choice of the combination parameterλ_n^$̣ which leads to a convergent
regularization method.

In Remark <ref> we have briefly discussed how to choose the combination parameter leading to the formulae
(<ref>) and (<ref>). However, these choices of λ_n^$̣ decrease to0as→̣0,
and consequently the acceleration effect will decrease as→̣0as well. Therefore, it is necessary to
find out other strategy for generatingλ_n^$̣ such that (<ref>) and (<ref>) hold.
We will adapt the discrete backtracking search (DBTS) algorithm introduced in <cit.> to our situation.
To this end, we take a function q:ℕ→ℕ that is non-increasing and

    ∑_i=0^∞ q(i)<∞.

The DBTS algorithm for choosing the combination parameter λ_n^$̣ in our method (<ref>)
is formulated in Algorithm <ref> below. Comparing with the one in <cit.>, there are two modifications:
The first modification is the definition ofβ_nin which we placeβ_n(i) ≤ n/(n+)instead ofβ_n(i) ≤ 1;
this modification gives the possibility to speed up convergence by making use of the Nesterov's acceleration strategy. The second modification
is in the “Else" part, where instead of settingλ_n^=̣0we calculateλ_n^$̣ by (<ref>); this
modification can provide additional acceleration to speed up convergence.


Givenξ_n^$̣,ξ_n-1^$̣, τ, δ, c_1, ν, q:ℕ→ℕ,
i_n-1^∈̣ℕ, j_max∈ℕSet γ_1 = c_1 p^* (2c_0)^p^*-1/νCalculateξ_n^-̣ξ_n-1^ and define, with ≥ 3,

    β_n(i) = min{q(i)/ξ_n^-̣ξ_n-1^, p^*(2 c_0)^p^*ρ^p/4 ξ_n^-̣ξ_n-1^^p^*, n/n+}.
Forj=1,…, j_max     Set λ_n^=̣β_n(i_n-1^+̣j);
     Calculate ζ_n^=̣ξ_n^+̣λ_n^(̣ξ_n^-̣ξ_n-1^)̣ and z_n^=̣∇Θ^*(ζ_n^)̣;
     Calculate μ_n^$̣ by (<ref>);
    Ify^-̣F(z_n^)̣≤τ$̣        λ_n^=̣ 0;
        i_n^=̣ i_n-1^+̣j;
        break;
    Else if(λ_n^+̣(λ_n^)^p^*) ξ_n^-̣ξ_n-1^^p^*≤γ_1 μ_n^F(z_n^)̣-y^^s        i_n^=̣ i_n-1^+̣ j;
        break;
    Else         calculate λ_n^$̣ by (<ref>);i_n^=̣ i_n-1^+̣j_max;
    End IfEnd ForOutput:λ_n^$̣, i_n^$̣Discrete backtracking search (DBTS-Θ) algorithm for λ_n^$̣,n≥ 1.

We need to show that the combination parameterλ_n^$̣ chosen by Algorithm <ref> satisfies
(<ref>), (<ref>) and (<ref>). From Algorithm <ref>
it is easily seen that 0≤λ_n^≤̣β(i_n^)̣. Therefore (<ref>) holds automatically.
When F(z_n^)̣-y^≤τ$̣, Algorithm <ref> givesλ_n^=̣0which ensures (<ref>)
hold. WhenF(z_n^)̣-y^ >τ$̣, Algorithm <ref>
either finds λ_n^$̣ of the formβ_n(i_n^)̣to satisfy (<ref>) or givesλ_n^$̣
by (<ref>) which again satisfies (<ref>). Thus Algorithm <ref> always produces
a λ_n^$̣ satisfying (<ref>).

For the noise-free case, Algorithm <ref> either produceλ_n=β_n(i_n)orλ_n = 0. Thus0≤λ_n ≤β_n(i_n). By the definition ofi_none can see thati_n≥ i_n-1+1and thusi_n≥ n.
Therefore, by usingβ_n(i) ≤ q(i)/ξ_n-ξ_n-1and the monotonicity ofq, we have

    ∑_n=0^∞λ_nξ_n-ξ_n-1   ≤∑_n=0^∞β_n(i_n)ξ_n-ξ_n-1≤∑_n=0^∞ q(i_n)
    ≤∑_n=0^∞ q(n) <∞.

Therefore (<ref>) is satisfied.

We can not use Theorem <ref> to conclude the regularization property of the two-point method (<ref>) when the
combination parameter is determined by Algorithm <ref> because the produced parameterλ_n^$̣ is not
necessarily continuously dependent on $̣ as→̣0. In fact,λ_n^$̣ may have many different cluster points
as →̣0. Using these different cluster points as the combination parameter in (<ref>) may lead to many
different iterative sequences for noise-free case. We need to consider all these possible iterative sequences altogether.
We will use Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) to denote the set consisting of all the iterative sequences
{(ξ_n, x_n)}⊂^*× defined by (<ref>), where the combination parameters {λ_n}
are chosen to satisfy

    1/p^*(2 c_0)^p^*-1(λ_n^p^*+λ_n) ξ_n-ξ_n-1^p^*≤c_1/νμ_n F(z_n)-y^s

and

    0≤λ_n ≤min{q(i_n)/ξ_n-ξ_n-1, p^*(2c_0)^p^*ρ^p/4 ξ_n-ξ_n-1^p^*,
    n/n+}

with a sequence {i_n} of integers satisfying i_0=0 and 1≤ i_n-i_n-1≤ j_max for all n.

Given a sequence {(ξ_n, x_n)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0), it is easy to check that the corresponding
combination parameters {λ_n} satisfy (<ref>), (<ref>) and (<ref>).
Therefore we may use Theorem <ref> to conclude the convergence of {x_n}. In the following we will show a uniform
convergence result.


Let all the conditions in Theorem <ref> hold. If (L(x^†)) ⊂(L(x)) for all x ∈ B_2ρ(x_0)∩(Θ),
then for any >0 there is an integer n() such that for any sequence {(ξ_n, x_n)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) there holds D_ξ_nΘ(x^†, x_n) < for all n≥ n().

Assume that the result is not true. Then there is an ε_0>0 such that for any ℓ≥ 1
there exist {(ξ_n^(ℓ), x_n^(ℓ))}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) and n_ℓ>ℓ such that

    D_ξ_n_ℓ^(ℓ)Θ (x^†, x_n_ℓ^(ℓ)) ≥ε_0.

We will construct, for each n=0, 1, ⋯, a strictly increasing subsequence {ℓ_n,k} of positive integers and
(ξ̂_n, x̂_n) ∈^* × such that


  (i){(ξ̂_n, x̂_n)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0);


  (ii) for each fixed n there hold x_n^(ℓ_n,k)→x̂_n, Θ(x_n^(ℓ_n,k)) →Θ(x̂_n) and ξ_n^(ℓ_n,k)→ξ̂_n as k→∞.

Assume that the above construction is available, we will derive a contradiction. According to (i), we may use Theorem <ref>
to conclude that D_ξ̂_nΘ(x^†, x̂_n) → 0 as n→∞. Thus we can pick a large
integer n̂ such that

D_ξ̂_n̂Θ(x^†, x̂_n̂)< ε_0/2.

Thus

    ε_0/2     > (D_ξ̂_n̂Θ(x^†, x̂_n̂)
    -D_ξ_n̂^(ℓ_n̂, k)Θ(x^†, x_n̂^(ℓ_n̂, k)))
    +D_ξ_n̂^(ℓ_n̂, k)Θ(x^†, x_n̂^(ℓ_n̂, k))   
        = -D_ξ_n̂^(ℓ_n̂, k)Θ(x̂_n̂, x_n̂^(ℓ_n̂, k))
    + łξ̂_n̂-ξ_n̂^(ℓ_n̂, k), x̂_n̂ -x^†
           + D_ξ_n̂^(ℓ_n̂, k)Θ(x^†, x_n̂^(ℓ_n̂, k)).

According to (ii), we have

-D_ξ_n̂^(ℓ_n̂, k)Θ(x̂_n̂, x_n̂^(ℓ_n̂, k))
+ łξ̂_n̂-ξ_n̂^(ℓ_n̂, k), x̂_n̂ -x^†→̊0    k→∞.

Hence we can find k̂ with ℓ̂:= ℓ_n̂, k̂≥n̂ such that

-D_ξ_n̂^(ℓ̂)Θ(x̂_n̂, x_n̂^(ℓ̂))
+ łξ̂_n̂-ξ_n̂^(ℓ̂), x̂_n̂ -x^†≥̊-_0/2.

Consequently, it follows from (<ref>) that

D_ξ_n̂^(ℓ̂)Θ(x^†, x_n̂^(ℓ̂))<_0.

Note that n_ℓ̂>ℓ̂= ℓ_n̂, k̂≥n̂. Thus, by the monotonicity
of {D_ξ_n^(ℓ̂)Θ(x^†, x_n^(ℓ̂))} with respect to n, we can obtain

D_ξ_n_ℓ̂^(ℓ̂)Θ(x^†, x_n_ℓ̂^(ℓ̂))
≤ D_ξ_n̂^(ℓ̂)Θ(x^†, x_n̂^(ℓ̂)) < ε_0

which is a contradiction to (<ref>) with ℓ=ℓ̂.

We turn to the construction of {ℓ_n,k} and (ξ̂_n, x̂_n), for each n=0, 1, ⋯, such that (i) and (ii) hold.
We use a diagonal argument. For n=0, we take (ξ̂_0, x̂_0)=(ξ_0, x_0) and ℓ_0,k=k for all k.
Since ξ_0^(k)=ξ_0 and x_0^(k)=x_0, (ii) holds automatically for n=0.

Next, assume that  we have constructed {ℓ_n,k} and (ξ̂_n, x̂_n) for all 0≤ n≤ m.
We will construct {ℓ_m+1, k} and (ξ̂_m+1, x̂_m+1). Recall the combination parameter λ_m^(ℓ_m,k)
and the integer i_m^(ℓ_m,k) involved in the definition of (ξ_m+1^(ℓ_m,k), x_m+1^(ℓ_m,k)). We have

0≤λ_m^(ℓ_m,k)≤m/m+     0 ≤ i_m^ (ℓ_m,k)≤ m j_max
for all k. Thus, we may pick a subsequence of {ℓ_m,k}, denoted by {ℓ_m+1, k}, such that

    lim_k→∞λ_m^(ℓ_m+1,k) =λ̂_m     
    i_m^(ℓ_m+1,k) = î_m   ∀ k

for some number 0≤λ̂_m ≤ m/(m+) and some integer î_m. We define
ζ̂_m = ξ̂_m + λ̂_m (ξ̂_m - ξ̂_m-1),     ẑ_m = ∇Θ^*(ζ̂_m).

By the induction hypothesis, (<ref>) and the continuity of ∇Θ^* we have

    ζ_m^(ℓ_m+1,k)→ζ̂_m     
    z_m^(ℓ_m+1,k)→ẑ_m     k →∞.

We define
μ̂_m = {[ min{μ̅_0 F(ẑ_m)-y^p(s-1)/L(ẑ_m)^* J_s^(F(ẑ_m)-y)^p, μ_1}F(ẑ_m)-y^p-s                                                           F(ẑ_m)  y,;                                                                    0                                                          F(ẑ_m) = y. ].

By using (<ref>), the continuity of F, L and J_s^, and the similar argument in the proof of Lemma <ref>,
we can show that

    μ_m^(ℓ_m+1,k)→μ̂_m     k→∞.

We next define

    ξ̂_m+1 = ζ̂_m - μ̂_m L(ẑ_m)^* J_s^ (F(ẑ_m)-y)     x̂_m+1 = ∇Θ^*(ξ̂_m+1).

By using (<ref>), (<ref>), and the continuity of F, L, J_s^ and ∇Θ^*,
it follows immediately that

    ξ_m+1^(ℓ_m+1,k)→ξ̂_m+1     x_m+1^(ℓ_m+1,k)→x̂_m+1   k→∞.

By the lower semi-continuity of Θ we then have

    Θ(x̂_m+1) ≤lim inf_k→∞Θ(x_m+1^(ℓ_m+1,k)).

On the other hand, by the definition of x_m+1^(ℓ_m+1,k) we have

x_m+1^(ℓ_m+1,k) = min_x∈{Θ(x) - łξ_m+1^(ℓ_m+1,k), x}.

Thus
Θ(x_m+1^(ℓ_m+1,k)) - łξ_m+1^(ℓ_m+1,k), x_m+1^(ℓ_m+1,k)≤̊Θ(x̂_m+1) - łξ_m+1^(ℓ_m+1,k), x̂_m+1.̊
Therefore, by virtue of (<ref>), we have

    lim sup_k→∞Θ(x_m+1^(ℓ_m+1,k))
       ≤Θ(x̂_m+1)+lim sup_k→∞łξ_m+1^(ℓ_m+1,k), x_m+1^(ℓ_m+1,k) -x̂_m+1    = Θ(x̂_m+1).

This together with (<ref>) implies that

    lim_k→∞Θ(x_m+1^(ℓ_m+1,k)) = Θ(x̂_m+1).

We thus complete the construction of {ℓ_m+1,k} and (ξ̂_m+1, x̂_m+1).

We still need to show that λ̂_m satisfies the requirements specified in the definition of
Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) in order to guarantee that the constructed sequence
is indeed in Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0).
Recall the definition of the sequences in Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) we have

    1/p^*(2 c_0)^p^*-1((λ_m^(ℓ_m+1, k))^p^*+λ_m^(ℓ_m+1, k))
    ξ_m^(ℓ_m+1, k)-ξ_m-1^(ℓ_m+1, k)^p^*
       ≤c_1/νμ_m^(ℓ_m+1, k)F(z_m^(ℓ_m+1, k))-y^s

and

0≤λ_m^(ℓ_m+1, k)≤min{q(i_m^(ℓ_m+1, k))/ξ_m^(ℓ_m+1, k)-ξ_m-1^(ℓ_m+1, k),
p^*(2c_0)^p^*ρ^p/4 ξ_m^(ℓ_m+1, k)-ξ_m-1^(ℓ_m+1, k)^p^*, m/m+}
with 1≤ i_m^(ℓ_m+1, k)-i_m-1^(ℓ_m+1, k)≤ j_max for all k. By using (<ref>), (<ref>)
and (<ref>), we may take k→∞ in the above two inequalities to conclude that

    1/p^*(2 c_0)^p^*-1(λ̂_m^p^*+ λ̂_m)
    ξ̂_m-ξ̂_m-1^p^*≤c_1/νμ̂_m F(ẑ_m)-y^s

and

0≤λ̂_m ≤min{q(î_m)/ξ̂_m-ξ̂_m-1,
p^*(2c_0)^p^*ρ^p/4 ξ̂_m-ξ̂_m-1^p^*, m/m+}
Furthermore, we can assume that {ℓ_m,k} was taken so that i_m-1^(ℓ_m,k) = î_m-1 for all k for some
integer î_m-1. Recall that {ℓ_m+1, k} is a subsequence of {ℓ_m,k}, we must have
1≤î_m - î_m-1 =  i_m^(ℓ_m+1, k)-i_m-1^(ℓ_m+1, k)≤ j_max. Note that i_0^(ℓ)=0
for all ℓ, we also have î_0 =0. The proof is therefore complete.

Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Let τ>1 and μ̅_0>0 be chosen to satisfy (<ref>). Let {y^_̣l} be a sequence of noisy data
satisfying y^_̣l-y≤_̣l with _̣l→ 0 as l→∞. Assume that the combination parameters
{λ_n^_̣l} are chosen by Algorithm <ref> with i_0^_̣l =0. Then, for any n≥ 0, by taking a
subsequence of {y^_̣l} if necessary, there is a sequence {(ξ_k, x_k)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0)
such that

    x_k^_̣l→ x_k,   ξ_k^_̣l→ξ_k
        Θ(x_k^_̣l)→Θ(x_k)   as  l→∞

for all 0≤ k≤ n.

We will use an induction argument on n. When n=0, nothing needs to be proved since x_0^_̣l=x_0 and
ξ_0^_̣l=ξ_0. Assume next that the result is true for n = m, i.e. there is a sequence
{(ξ_k, x_k)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) such that

x_k^_̣l→ x_k,   ξ_k^_̣l→ξ_k     Θ(x_k^_̣l) →Θ(x_k)    l →∞
for all 0≤ k ≤ m. In order to show that the result is also true for n = m+1, we will obtain a sequence from
Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) by retaining the first m+1 terms in {(ξ_k, x_k)}
and modifying the remaining terms. It suffices to redefine ξ_m+1 and x_m+1 because we can then apply
the method (<ref>) with λ_n =0 for n≥ m+1 to produce the remaining terms.

Note that the combination parameter λ_m^_̣l determined by Algorithm <ref> satisfies

    1/p^*(2c_0)^p^*-1((λ_m^_̣l)^p^* + λ_m^_̣l) ξ_m^_̣l-ξ_m-1^_̣l^p^*≤c_1/νμ_m^_̣lF(z_m^_̣l)-y^_̣l^s

and

    0≤λ_m^_̣l≤min{q(i_m^_̣l)/ξ_m^_̣l-ξ_m-1^_̣l,
    p^*(2c_0)^p^*ρ^p/4ξ_m^_̣l-ξ_m-1^_̣l^p^*, m/m+} 
       λ_m^_̣l = min{(2 c_0)^p^*-1 p^* Mτ^p _̣l^p/2 ξ_m^_̣l -ξ_m-1^_̣l^p^*, m/m+}

with 1 ≤ i_m^_̣l - i_m-1^_̣l≤ j_max. By taking a subsequence of {y^_̣l} if necessary,  we have

    lim_l→∞λ_m^_̣l = λ_m     
    i_m^_̣l = i_m    l

for some number 0≤λ_m ≤ m/(m+a) and some integer i_m. We now define ζ_m, z_m, μ_m, ξ_m+1 and
x_m+1 by (<ref>) with n=m. By using (<ref>) and the similar argument in the proof of Lemma <ref>,
we have

    ζ_m^_̣l→ζ_m,    z_m^_̣l→ z_m,   μ_m^_̣l→μ_m,   ξ_m+1^_̣l→ξ_m+1     x_m+1^_̣l→ x_m+1

as l→∞. Furthermore, we can use the similar argument for proving (<ref>) to obtain
Θ(x_m+1^_̣l) →Θ(x_m+1) as l→∞. Now by taking l→∞
in (<ref>) and (<ref>) and using the induction hypothesis, we can conclude that λ_m
and i_m satisfy the requirements specified in the definition of Γ_μ̅_0, μ̅_1, ν q(ξ_0, x_0).


We are now ready to show the regularization property of the method (<ref>) when the combination parameter λ_n^$̣
is chosen by Algorithm <ref>.


Let  be reflexive, let  be uniformly smooth, and let Assumption <ref> and Assumption <ref> hold.
Let τ>1 and μ̅_0>0 be chosen to satisfy (<ref>). Assume that the combination parameters {λ_n^_̣l}
are chosen by Algorithm <ref> with i_0^_̣l =0. Let n_$̣ be the integer determined by the
discrepancy principle (<ref>). If(L(x^†)) ⊂(L(x))for allx∈ B_3ρ(x_0), thenlim_→̣0x_n_^-̣ x^† =0     lim_→̣0 D_ξ_n_^Θ(x^†, x_n_^)̣ =0.We will show that

    D_ξ_n_^Θ(x^†, x_n_^)̣→ 0   →̣0

by a contradiction argument. If it is not true, then there is a number >0 and a subsequence {y^_̣l} of {y^}̣
satisfying y^_̣l-y≤_̣l with _̣l→ 0 as l→∞ such that

    D_ξ_n__̣l^_̣lΘ(x^†, x_n__̣l^_̣l) ≥   l.

We need to consider two cases.

 Case 1.{n__̣l} has a finite cluster point n̂. For this case, by using Lemma <ref> and by
taking a subsequence if necessary, we have n__̣l=n̂ for all l and there is a sequence
{(ξ_k, x_k)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) such that

x_n̂^_̣l→ x_n̂,   ξ_n̂^_̣l→ξ_n̂    Θ(x_n̂^_̣l) →Θ(x_n̂)    l →∞.

Using the similar argument in the proof of Theorem <ref> we can obtain
x_n̂ = x^†. Therefore

x_n__̣l^_̣l→ x^†,    ξ_n__̣l^_̣l→ξ_n̂    Θ(x_n__̣l^_̣l) →Θ(x^†)    l→∞.

Consequently D_ξ_n__̣l^_̣lΘ(x^†, x_n__̣l^_̣l) → 0 as l→∞ which contradicts (<ref>).

 Case 2.lim_l→∞ n__̣l=∞. According to Proposition <ref>, there is an integer
n() such that

    D_ξ_n(ε)Θ(x^†, x_n()) <  {(ξ_n, x_n)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0).

For this n(), by using Lemma <ref> and by taking a subsequence if necessary, we can find
{(ξ_n, x_n)}∈Γ_μ̅_0, μ̅_1, ν, q(ξ_0, x_0) such that

    x_n^_̣l→ x_n,   ξ_n^_̣l→ξ_n     Θ(x_n^_̣l) →Θ(x_n)
       l →∞

for 0≤ n≤ n(). Since lim_l→∞ n__̣l=∞, we have n__̣l>n() for large l.
Thus, it follows from (<ref>) in Proposition <ref> that

D_ξ_n__̣l^_̣lΘ(x^†, x_n__̣l^_̣l) ≤ D_ξ_n()^_̣lΘ(x^†, x_n()^_̣l).

Letting l→∞ and using (<ref>) and (<ref>) we can obtain

    lim sup_l→∞ D_ξ_n__̣l^_̣lΘ(x^†, x_n__̣l^_̣l)
    ≤lim sup_l→∞ D_ξ_n()^_̣lΘ(x^†, x_n()^_̣l)
    = D_ξ_n()Θ(x^†, x_n()) <

which is a contradiction to (<ref>).

Combining the above two cases, we therefore obtain (<ref>). By virtue of the p-convexity of Θ,
we then have x_n_^ -x^†→ 0 as →̣0.
 A two-point gradient method in Hilbert spaces was considered in <cit.> in which the combination parameter
is chosen by a discrete backtracking search (DBTS) algorithm. The regularization property was proved under the condition
that the noise-free counterpart of the method never terminates at a solution of (<ref>) in finite many steps if the
combination parameter is chosen by the DBTS algorithm. This technical condition seems difficult to be verified
because the exact data y is unavailable. In Theorem <ref> we removed this condition by using a uniform
convergence result established in Proposition <ref>.

§  NUMERICAL SIMULATIONS


In this section we will present numerical simulations on our TPG-DBTS method, i.e. the two point gradient method (<ref>)
with the combination parameterλ_n^$̣ chosen by the DBTS algorithm (Algorithm <ref>). In order to
illustrate the performance of TPG-DBTS algorithm, we will compare the computational results with the ones obtained by the
Landweber iteration  (<ref>) and the Nesterov acceleration of Landweber iteration, i.e. the method (<ref>) with λ_n^=̣ n/(n+) for some ≥ 3. In order to be fair, the step sizes μ_n^$̣ involved in all these methods are
computed by (<ref>) and all the iterations are terminated by the discrepancy principle withτ = 1.05.

A key ingredient for the numerical implementation is the determination ofx= ∇Θ^*(ξ)for any givenξ∈^*which is equivalent to solving the minimization problem

    x=min_z∈{Θ(z) -łξ, z}.

For some choices ofΘ, this minimization problem can be easily solved numerically. For instance, when=L^2(Ω)and the sought solution is piecewise constant, we may choose

    Θ(x) = 1/2βx_2^2 + |x|_TV

with a constantβ>0, where|x|_TVdenotes the total variation ofx. Then the minimization problem (<ref>)
becomes the total variation denoising problem

    x =  argmin_z∈ L^2(Ω){1/2βz-βξ_2^2 + |z|_TV}

which is nonsmooth and convex. Note that for thisΘ, Assumption <ref> holds withp=2andc_0 = 1/2β.
Many efficient algorithms have been developed for solving (<ref>), including the fast iterative shrinkage-thresholding
algorithm <cit.>, the alternating direction method of multipliers <cit.>, and the primal dual hybrid gradient (PDHG)
method <cit.>.

In the following numerical simulations we will only consider the situation that the sought solution is piecewise constant.
We will use the PDHG method to solve (<ref>) iteratively. Our simulations are performed via MATLAB R2012a
on a Lenovo laptop with Intel(R) Core(TM) i5 CPU 2.30GHz and 6GB memory.



 §.§ Computed tomography


Computed tomography (CT) consists in determining the density of cross sections of a human body by measuring the attenuation of X-rays as
they propagate through the biological tissues. Mathematically, it requires to determine a function supported on a bounded domain from
its line integrals <cit.>. In order to apply our method to solve the CT problems, we need a discrete model. In our numerical experiment,
we assume that the image is supported on a rectangular domain inℝ^2which is divided intoI× Jpixels so that the discrete
image has sizeI× Jand can be represented by a vectorx∈ℝ^NwithN=I× J. We further assume that
there aren_θprojection directions and in each direction there arepX-rays emitted. We want to reconstruct the image by using
the measurement data of attenuation along the rays which can be represented by a vectorb∈ℝ^MwithM=n_θ× p.
According to a standard discretization of the Radon  transform (<cit.>), we arrive at a linear algebraic systemFx=b,whereFis a sparse matrix of sizeM× Nwhose form depends on the scanner geometry.


 Computed tomography (β=1, τ=1.05).
_̣rel    method         n_$̣    CPU time (s)x_n^-̣x^†_2/x^†_20.05           Landweber          109         22.33          0.19431     
      Nesterov           34          7.45           0.17573     
      TPG-DBTS           34          11.99          0.17573     

0.01           Landweber          489         106.97         0.06410      
      Nesterov           79          16.25          0.05923      
      TPG-DBTS           79          29.51          0.05923      

0.005          Landweber          879         197.65         0.03578      
      Nesterov           109         22.25          0.03165      
      TPG-DBTS           109         43.02          0.03165      

0.001          Landweber          3299        775.10         0.00694      
      Nesterov           247         51.05          0.00521      
      TPG-DBTS           247         105.29         0.00521      

0.0005         Landweber          5703        1328.31        0.00326     
      Nesterov           351         82.93          0.00199     
      TPG-DBTS           351         153.79         0.00199     



In the numerical simulations we consider only test problems that model the standard 2D parallel-beam tomography. The true image is
taken to be the modified Shepp-Logan phantom of size256× 256generated by MATLAB. This phantom is widely used in evaluating
tomographic reconstruction algorithms. We use the full angle model with 45 projection angles evenly distributed between 1 and 180 degrees,
with 367 lines per projection. The function  in MATLAB package AIR TOOLS <cit.> is used to generate the
sparse matrixF, which has the sizeM=16515andN=66536. Letx^†denote the vector formed by stacking all the columns of
the true image and letb = F x^†be the true data. We add Gaussian noise onbto generate a noisy datab^$̣ with relative noise level
δ_rel=b^-̣b_2/b_2 so that the noise level is =̣_̣relb_2. We will use b^$̣ to reconstructx^†.
In order to capture the feature of the sought image, we takeΘto be the form (<ref>) withβ =1.
In our numerical simulations we will useξ_0 =0as an initial guess. According to (<ref>) we needμ̅_0 <2 (1-1/τ)/β.
Therefore we take the parametersμ̅_0andμ̅_1in the definition ofμ_n^$̣ to be μ̅_0 = 1.8(1-1/τ)/β
and μ̅_1 = 20000. For implementing TPG-DBTS method with λ_n^$̣ chosen by Algorithm
<ref>, we takej_max = 1,=5,γ_0 = 0.1in (<ref>),γ_1 = 0.4; we also choose the functionq: ℕ→ℕbyq(m) = m^-1.1. For implementing the Nesterov acceleration of Landweber iteration, we takeλ_n^=̣ n/(n+)with = 5. During the computation, the total variation denoising problem (<ref>) involved in
each iteration step is solved approximately by the PDHG method after 100 iterations.



    < g r a p h i c s >

The computed tomography using noisy data with relative noise level _̣rel = 0.01.


The computational results by TPG-DBTS, Landweber, and Nesterov acceleration of Landweber are reported in Table <ref>, including
the number of iterationsn_$̣, the CPU running time and the relative errors x_n_^-̣x^†_2/x^†_2, using noisy data with
various relative noise level _̣rel>0. Table <ref> shows that both TPG-DBTS and Nesterov acceleration, terminated by the discrepancy principle,
lead to a considerable decrease in the number of iterations and the amount of computational time, which demonstrates that these two methods
have the striking acceleration effect. Moreover, both TPG-DBTS and Nesterov acceleration produce more accurate results than Landweber iteration.
With the above setup, our computation shows that TPG-DBTS produces the combination parameter λ_n^$̣ which is exactly
same as the combination parametern/(n+)in Nesterov acceleration in each iteration step. Therefore, TPG-DBTS and Nesterov acceleration
require the same number of iterations and produce the same reconstruction result. Because TPG-DBTS spends more time on determiningλ_n^$̣,
the Nesterov acceleration requires less amount of computational time than TPG-DBTS. However, unlike TPG-DBTS, there exists no convergence result
concerning Nesterov acceleration for ill-posed inverse problems.

In order to visualize the reconstruction accuracy of the TPG-DBTS method, we plot in Figure <ref> the true image, the reconstruction
result by TPG-DBTS using noisy data with relative noise level _̣rel = 0.01, the curve of λ_n^$̣ versusn, and the relative errorx_n^-̣x^†_2/x^†_2versusnfor TPG-DBTS and Landweber iteration.




 §.§ Elliptic parameter identification


We consider the identification of the parametercin the elliptic boundary value problem

    {[ - u + cu = f        in Ω,;          u=g        on ∂Ω ].

from anL^2(Ω)-measurement of the stateu, whereΩ⊂ℝ^dwithd≤ 3is
a bounded domain with Lipschitz boundaryΩ,f∈ H^-1(Ω)andg∈ H^1/2(Ω). We assume that
the sought parameterc^†is inL^2(Ω). This problem reduces to solvingF(c) =uif we define the
nonlinear operatorF: L^2(Ω)→  L^2(Ω)by

    F(c): = u(c) ,

whereu(c)∈ H^1(Ω)⊂  L^2(Ω)is the unique solution of (<ref>). This operatorFis well defined on

    : = { c∈  L^2(Ω) : c-ĉ_ L^2(Ω)≤_0 ĉ≥0, a.e.}

for some positive constant_0>0. It is well-known (<cit.>) that the operatorFis weakly closed
and Fréchet differentiable withF'(c) h = v         F'(c)^* σ = - u(c) wforc∈andh, σ∈ L^2(Ω), wherev, w∈ H^1(Ω)are the unique
solutions of the problems

    {[ - v + cv = -h u(c)   Ω,;                 v=0   Ω ].     {[ - w + cw = σ  Ω,;          w=0   Ω ].

respectively. Moreover,Fsatisfies Assumption <ref> (c).

 2-dimensional elliptic parameter estimation (β=10, τ=1.05).
$̣      methodn_$̣    CPU time (s)   c_n^-̣c^†_L^2

0.005          Landweber          34          9.41          0.26884      
      Nesterov           20          5.74          0.24815      
      TPG-DBTS           21          6.22          0.24815      

0.001          Landweber          184         54.49         0.12831      
      Nesterov           72          20.58         0.11466      
      TPG-DBTS           72          24.33         0.11466      

0.0005         Landweber          263         74.15         0.11096      
      Nesterov           111         29.96         0.09485      
      TPG-DBTS           111         41.19         0.09485      

0.0001         Landweber         1042         324.49        0.08107     
      Nesterov           222         65.26         0.07688     
      TPG-DBTS           222         103.64        0.07688     

0.00005        Landweber          2346        785.07        0.07221     
      Nesterov           452         147.03        0.06609     
      TPG-DBTS           434         242.24        0.06342     




In our numerical simulation, we consider the two-dimensional problem withΩ = [0, 1]× [0,1]and
the sought parameter is assumed to be

    c^†(x,y) = {[                                    1, if  (x-0.65)^2 + (y-0.36)^2 ≤ 0.18^2,;                                  0.5, if  (x-0.35)^2 + 4(y-0.75)^2 ≤ 0.2^2,;                                    0,                            elsewhere. ].

Assumingu(c^†) = x+y, we add random Gaussian noise to produce noisy datau^$̣ satisfying u^-̣u(c^†)_L^2(Ω)≤$̣
with various noise level>̣0. We will useu^$̣ to reconstruct c^†. In order to capture the feature of the sought parameter,
we take Θ to be the form (<ref>) with β = 10. We will use the initial guess ξ_0=0 to carry out the iterations.
The parameters μ̅_0 and μ̅_1 in the definition of μ_n^$̣ are taken to beμ̅_0 = (1-1/τ)/βandμ̅_1 = 20000. For implementing TPG-DBTS method withλ_n^$̣ chosen by Algorithm
<ref>, we take j_max = 1, =5, γ_0 = 0.1 in (<ref>), γ_1 = 0.3; we also choose the function
q: ℕ→ℕ by q(m) = m^-1.2. For implementing the Nesterov acceleration of Landweber iteration, we take
λ_n^=̣ n/(n+) with = 5. In order to carry out the computation, we divide Ω into 128× 128 small squares
of equal size and solve all partial differential equations involved approximately by a multigrid method (<cit.>) via finite difference
discretization. The total variation denoising problem (<ref>) involved in each iteration step is solved by the PDHG method
after 200 iterations.



    < g r a p h i c s >

The 2-dimensional elliptic parameter identification using noisy data with noise level δ = 0.0001.

In Table <ref> we report the computational results by TPG-DBTS, Landweber, and Nesterov acceleration of Landweber, including
the number of iterations n_$̣, the CPU running time and the absolute errorsc_n_^-̣c^†_L^2(Ω), for various
noise level>̣0. Table <ref> shows that both TPG-DBTS and Nesterov acceleration, terminated by the discrepancy principle,
reduce the number of iterations and the amount of computational time significantly, and produce more accurate results
than Landweber iteration. This demonstrates that these two methods have a remarkable acceleration effect.
In order to visualize the reconstruction accuracy of the TPG-DBTS method, we plot in Figure <ref> the true solution, the reconstruction
result by TPG-DBTS with noise level=̣ 0.0001, the curve ofλ_n^$̣ versus n, and the error c_n^-̣c^†_L^2(Ω) versus n
for TPG-DBTS and Landweber iteration.





 §.§ Robin coefficient reconstruction


We consider the heat conduction process in a homogeneous solid rod located on the interval [0,π]. If the endpoints
of the rod contacts with liquid media, then the convective heat transfer occurs. The temperature field of u(x,t) during
a time interval [0,T] with a fixed time of interest T>0 can be modeled by

    {[                                u_t-a^2 u_xx = 0,                                x∈(0,π), t∈(0,T);; u_x(0,t) = f(t),   u_x(π,t) + σ(t)u(π,t) = φ(t),                                         t∈[0,T];;                                 u(x,0) = u_0(x),                                        x∈ [0,π]. ].

The function σ(t)≥ 0 represents the corrosion damage, which is interpreted as a Robin coefficient of energy exchange.
We will assume f, φ and u_0 are all continuous. Notice that if σ(t) is given, the problem (<ref>) is a
well-posed direct problem. The inverse problem of identifying the Robin coefficient σ(t) requires additional data to be
specified. We consider the reconstruction of σ(t) from the temperature information measured at the boundary

    u(0,t) = g(t),   t∈[0,T] .

Define

    : = {σ∈ L^2[0,T], 0<σ_-≤σ≤σ_+, a.e. in [0,T]},

and define the nonlinear operator F:σ∈→ u[σ](0,t)∈ L^2[0,T], where u[σ] denotes the unique solution of (<ref>).
Then the above Robin coefficient inversion problem reduces to solving the equation F(σ) = g. We refer to <cit.> for the
well-posedness of F and the uniqueness of the inverse problem in the L^2 sense.
By the standard theory of parabolic equation, one can show that F is Fréchet differentiable in the sense that

    F(σ+h) - F(σ) - F'(σ) h _L^2(0,T) = o(h_L^2(0,T))

for all σ,σ+h∈, where [F'(σ)h](t) = w(0,t) and w is the unique solution of

    {[                                             w_t-a^2 w_xx = 0,                                             x∈(0,π), t∈(0,T);; w_x(0,t) = 0,     w_x(π,t) + σ(t) w(π,t) = - h(t) u[σ](π, t),                                                      t∈[0,T];;                                                   w(x,0) = 0,                                                     x∈ [0,π]. ].

In addition, the adjoint of the Fréchet derivative is given by

    [F'(σ)^*ζ](t) =  u[σ](π,t)v(π,τ),

where v(x,t) solves the adjoint system

    {[                             -v_t-a^2 v_xx = 0,                              x∈(0,π), t∈(0,T);; v_x(0,t) = ζ(t),    v_x(π,t) + σ(t) v(π,t) =0,                                       t∈[0,T];;                                    v(x,T) = 0,                                      x∈ [0,π]. ].
Robin coefficient reconstruction (β=1, τ=1.05).
δ    Method        n_$̣   CPU time (s)σ_n^-̣σ^†_L^20.1           Landweber       242          1.18        0.120974     
      Nesterov       76           0.56        0.119851     
      TPG-DBTS       76           1.03        0.119851       

0.01          Landweber       1431         5.33        0.042751     
      Nesterov       240          1.14        0.038788     
      TPG-DBTS       235          2.63        0.038193       

0.001          Landweber      12033        42.14       0.007127      
      Nesterov       706          2.76        0.002754       
      TPG-DBTS       732          7.44        0.002569        

0.0001         Landweber      31021        109.85      0.000702     
      Nesterov       959          3.65        0.000236     
      TPG-DBTS       2113         18.89       0.000102      




In our numerical simulations, we takea=5,T=1, and assume the sought Robin coefficient is

    σ^†(t) = {[               1.5,     0≤ t ≤ 0.1563,;                 2, 0.1563 <t≤ 0.3125,;               1.2, 0.3125 <t≤ 0.5469,;               2.5, 0.5469 <t≤ 0.6250,;               1.8, 0.6250 <t≤ 0.7813,;                 1,      0.7813 <t≤ 1. ].

We also assume that the exact solution of the forward problem (<ref>) withσ= σ^†is

    u(x,t) = e^-a^2tsin x+x^2 + 2a^2t

through which we can obtain the expression of(f(t),u_0(x),φ(t))and the inversion inputg(t):=u(0, t).
We add random Gaussian noise ongto produce noisy datag^$̣ satisfying g^-̣g_L^2(0, T)≤$̣ with various noise
level>̣0. We will useg^$̣ to reconstruct σ^†. In order to capture the feature of the sought Robin coefficient,
we take Θ to be the form (<ref>) with β = 1. We will use the initial guess ξ_0=0 to carry out the computation.
The parameters μ̅_0 and μ̅_1 in the definition of μ_n^$̣ are taken to beμ̅_0 = 1.8(1-1/τ)/βandμ̅_1 = 20000. For implementing TPG-DBTS method withλ_n^$̣ chosen by Algorithm
<ref>, we take j_max = 2, =5, γ_0 = 0.1 in (<ref>), γ_1 = 0.4; we also choose the function
q: ℕ→ℕ by q(m) = m^-1.1. For implementing the Nesterov acceleration of Landweber iteration, we take
λ_n^=̣ n/(n+) with = 5. During the computation, the initial-boundary value problems for parabolic equation are
transformed into integral equations by the potential theory (<cit.>) and then solved by a boundary element method by dividing
[0, T] into N = 64 subintervals of equal length. The total variation denoising problem (<ref>) involved in each iteration
step is solved approximately by the PDHG method after 200 iterations.

In Table <ref> we report the computational results by TPG-DBTS, Landweber, and Nesterov acceleration of Landweber, using noisy data
for various noise level >̣0, which clearly demonstrates the acceleration effect of TPG-DBTS and Nesterov acceleration and shows that
these two methods have superior performance over Landweber iteration. In Figure <ref> we also plot the computational results by TPG-DBTS
using noisy data with noise level =̣ 0.001. We note that the combination parameter λ_n^$̣ produced by TPG-DBTS may be different fromn/(n+)for somen, but eventuallyλ_n^$̣ becomes the same as the combination parameter n/(n+) in Nesterov acceleration.





    < g r a p h i c s >

Robin coefficient reconstruction with δ = 0.001.

§  ACKNOWLEDGEMENT


The work of M Zhong is partially supported by the National Natural Science Foundation of China (No. 11501102) and Natural Science Foundation of Jiangsu Province (No. BK20150594).
The work of W Wang is partially supported by the National Natural Science Foundation of China (No. 11871180) and Natural Science Foundation of Zhejiang Province (No. LY19A010009).
The work of Q Jin is partially supported by the Future Fellowship of the Australian Research Council (FT170100231).


99AP2016
H Attouch and J. Peypouquet,  The rate of convergence of Nesterov’s accelerated forward-backward method is actually
faster than O(1/k^2),  SIAM J. Optim.,  26 (2016), 1824–1834.

bt09
A. Beck and M. Teboulle,   A fast iterative shrinkage-thresholding algorithm for linear inverse
problems, SIAM J. Imag. Sci.,  2 (2009), 183–202.

bt09(2)
A. Beck and M. Teboulle,  Fast gradient-based algorithms for constrained total variation image
denoising and deblurring problems, IEEE Trans. Image Process, 18 (2009),  2419–2434.

bh12
R. Bot and T. Hein,  Iterative regularization with a geeral penalty term|theory and applications
to L^1 and TV regularization, Inverse Problems, 28 (2012), 104010.

BPCPE2011
S. Boyd, N. Parikh, E. Chu, B. Peleato and J. Eckstein,
 Distributed optimization and statistical learning via the alternating direction method of multipliers,
 Foundations and Trends^® in Machine Learning, 3 (2011), 1–122.

ehn96
H. W. Engl, M. Hanke and A. Neubauer,  Regularization of Inverse Problems, Dordrecht: Kluwer,
1996.

H2016
W. Hackbusch,  Iterative Solution of Large Sparse Systems of Equations, Second edition, Applied
Mathematical Sciences, 95. Springer, 2016. xxiii+509 pp.

H1991
M. Hanke,  Accelerated Landweber iterations for the solution of ill-posed equations,
Numer. Math., 60 (1991), 341–373.

H1997
M. Hanke,  A regularizing Levenberg-Marquardt scheme with applications to inverse groundwater
filtration problems, Inverse Problems, 13 (1997), 79–95.

hns95
M. Hanke, A. Neubauer and O. Scherzer,  A convergence analysis of the Landweber iteration
for nonlinear ill-posed problems, Numer. Math., 72 (1995), 21–37.

Hansen2012
P. C. Hansen and M. Saxild-Hansen,  AIR Tools — a MATLAB package of algebraic iterative
reconstruction methods, J. Comput. Appl. Math., 236 (2012), 2167–2178.

HJW2015
M. Hegland, Q. Jin and W. Wang,  Accelerated Landweber iteration with convex penalty for
linear inverse problems in Banach spaces, Appl. Anal., 94 (2015), 524–547.

hk10
T. Hein and K. S. Kazimierski,   Accelerated Landweber iteration in Banach spaces, Inverse Problems,
26 (2010), 1037–1050.

hr17
S. Hubmer and R. Ramlau,  Convergence Analysis of a Two-Point Gradient Method for Nonlinear Ill-Posed Problems,
Inverse Problems, 33 (2017), 095004.

hr18
S. Hubmer and R. Ramlau,  Nesterov's accelerated gradient method for nonlinear ill-posed problems with a locally
convex residual functional, Inverse Problems, 34 (2018), 095003.

j12
Q. Jin,  Inexact Newton-Landweber iteration for solving nonlinear inverse problems in Banach
spaces, Inverse Problems, 28 (2012), 065002.

j16
Q. Jin,  Landweber-Kaczmarz method in Banach spaces with inexact inner solvers,
Inverse Problems, 32 (2016), 104005.

jw13
Q. Jin and W. Wang,  Landweber iteration of Kaczmarz type with general non-smooth convex penalty functionals,
Inverse Problems, 29 (2013), 085011.

jy16
Q. Jin and H. Yang,  Levenberg-Marquardt method in Banach spaces with general convex regularization
terms, Numer. Math. 133 (2016),  655–684.

jz13
Q. Jin and M. Zhong,  On the iteratively regularized Gauss-Newton method in Banach spaces with applications to parameter
identification problems, Numer. Math., 124 (2013), 647–683.

jz14
Q. Jin and M. Zhong,  Nonstationary iterated Tikhonov regularization in Banach spaces with
uniformly convex penalty terms,  Numer. Math.,  127 (2014) 485–513.

k15
B. Kaltenbacher,  A convergence rates result for an iteratively regularized Gauss-Newton
Halley method in Banach space, Inverse Problems, 31 (2015), 015007.

kss09
B. Kaltenbacher, F. Schöpfer and T. Schuster,  Iterative methods for nonlinear ill-posed problems in
Banach spaces: convergence and applications to parameter identification problems, Inverse Problems,
25 (2009), 065003.

K1989
R. Kress,  Linear Integral Equations, Springer-Verlag, Berlin, 1989.

N2001
F. Natterer,  The Mathematics of Computerized Tomography, SIAM, Philadelphia, 2001.

n83
Y. Nesterov,  A method of solving a convex programming problem with convergence
rate O(1/k^2),  Soviet Mathematics Doklady, 27 (1983), 372–376.

sls06
F. Schöpfer, A. K. Louis and T. Schuster,  Nonlinear iterative methods for linear ill-posed
problems in Banach spaces, Inverse problems, 22 (2006), 311–329.

sls09
F. Schöpfer and T. Schuster,  Fast regularizing sequential subspace optimization in Banach spaces,
Inverse Problems, 25 (2009), 015013.

wwh18
J. Wang, W. Wang and B. Han,  Accelerated Landweber iteration with non-smooth
convex penalty for nonlinear inverse problems in Banach spaces, preprint 2018.

wl17
Y. C. Wang and J. J. Liu,  Identification of non-smooth boundary heat dissipation by partial boundary data,
Applied Mathematics Letters,  69 (2017),  42–48.

ZA02
C. Zălinscu,  Convex Analysis in General Vector Spaces, World Scientific Publishing Co., Inc., River Edge, New Jersey, 2002.



zw19
M. Zhong and W. Wang,  A regularizing multilevel approach for nonlinear inverse problems, Applied Numerical Mathematics, 32:(2019), 297–331 


ZC2008 M. Zhu and T. F. Chan,  An efficient primaldual hybrid gradient algorithm for Total Variation
image restoration, CAM Report 08-34, UCLA (2008).