Tangential cone condition and Lipschitz stability for the full waveform forward operator in the acoustic regime

Time-domain full waveform inversion (FWI) in the acoustic regime comprises a parameter identification problem for the acoustic wave equation: pressure waves are initiated by sources, get scattered by the earth’s inner structure, and their reflected parts are picked up by receivers located on the surface. From these reflected wave fields the two parameters, density and sound speed, have to be reconstructed. Mathematically, FWI reduces to the solution of a nonlinear and ill-posed operator equation involving the parameter-to-solution map of the wave equation. Newton-like iterative regularization schemes are well suited and well analyzed to tackle this inverse problem. Their convergence results are often based on an assumption about the nonlinear map known as tangential cone condition. In this paper we verify this assumption for a semi-discrete version of FWI where the searched-for parameters are restricted to a finite dimensional space. As a byproduct we establish that the semi-discrete seismic inverse problem is locally Lipschitz stable, in particular, it is conditionally well-posed.


Introduction
Time-domain full waveform inversion (FWI) is the up-to-date geophysical imaging technique being capable to exploit the full information content of recorded seismic waves which have been excited locally by controlled sources or globally by earthquakes, see, e.g. [9,28] for a presentation from a geophysical point of view. Mathematically, FWI entails a parameter identification task for the underlying wave propagation model. In the acoustic regime, where the medium does not support shear stress, wave propagation is modeled by the acoustic/compressional wave equation and the searched-for parameters are pressure wave speed v p and bulk density . In this paper we rigorously explore a mathematical aspect of FWI in the acoustic setting.
Let S be the corresponding parameter-to-solution (parameter-to-state) map that maps (v p , ) to (p, v), the pressure and velocity fields (solution of the acoustic wave equation with respect to v p and ). Then, FWI in the acoustic regime means the solution of the (nonlinear) equation ΨS(v p , ) = Ψ(p, v), (1.1) where the observation operator Ψ models the measurement process (in the geophysical language, Ψ(p, v) collects the seismograms). Since the inverse problem (1.1) is locally ill-posed (in the infinite-dimensional setting, see [15,16]), it requires to be regularized. Newton-like iterative regularization schemes are well-established workhorses for stably solving nonlinear ill-posed problems like (1.1), see, e.g. [10,12,19,21,24]. To simplify the notation let us consider a generic setting where F : D(F) ⊂ V → W denotes the underlying nonlinear operator between Banach spaces. Assume we want to solve F(x) = y for a given y ∈ W by a Newton-like scheme. The generic iteration rule reads: choose a starting guess x 0 ∈ D(F) and iterate according to Here, F is the Fréchet-derivative (F-derivative) of F. The various methods differ in how they determine the Newton step s k from the locally linearized equation. Their convergence analysis, however, often relies on the same structural assumption which is known under the name tangential cone condition (TCC). It can be traced back to [22] and reads: F satisfies the TCC at x + ∈ int(D(F)) if for an η < 1 (sufficiently small) where B r (x + ) is the open ball in V of radius r > 0 about x + (sometimes the TCC is formulated in a ball with respect to a Bregman distance).
In the fully continuous (infinite dimensional) setting only a few academic examples of nonlinear ill-posed problems are known for which TCC holds, see, e.g. [10, section 4] and [12, section 2.4], but consult [13] for recent developments. However, in a semi-discrete setting the situation is more relaxed. For instance, a semi-discrete TCC has been derived for the inverse problem of the complete electrode model in 2D-electrical impedance tomography [18]. It turns out that injectivity of F (x + ) is essentially sufficient not only to yield a semi-discrete TCC but also a Lipschitz stability like where c > 0 is a constant. We will demonstrate this implication under rather general assumptions. Semi-discrete Lipschitz estimates and conditional well-posedness for various inverse problems have already been derived, e.g. in [1][2][3][4][5] and we add time-domain FWI in the acoustic regime to this list. In fact, if we confine v p and to suitable finite dimensional spaces, the F-derivative of ΨS is one-to-one. We need to emphasize that, in general, neither (1.3) nor (1.2) carry over to the continuous setting when the finite dimension of the semi-discrete setting is increased. Typically, the constant in (1.3) blows up while the radius r decreases. This will occur for our FWI application, since otherwise (1.1) would be locally well-posed in the infinite-dimensional formulation.
The presentation of our findings is organized as follows. In the next section we set the stage by introducing the acoustic wave equation as a first order system and by recalling existence and uniqueness results. Then, in section 3 we define our semi-discrete model where sources are fired in one part Σ of the propagation medium D and the resulting wave fields are recorded at a different part Ω of D. Here, sound speed and bulk density are expressed as linear combinations of smooth basis functions which are locally independent in D: if a linear combination vanishes on an open subset of D it vanishes globally in D. For this model we formulate the seismic inverse problem and characterize the F-derivative of the forward map by a different but akin acoustic wave equation. For this wave equation we show a fundamental property in proposition 3.1: there is a source supported in Σ such that the wave field of the F-derivative does not vanish identically on Ω. The proof is based on Holmgren's uniqueness theorem and the propagation of singularities along bicharacteristics of the wave operator. Finally, section 4 presents the main result (theorem 4.1) which originates from a local uniqueness statement that even holds when only one single source is fired (remark 4.4).
So as not to distract the reader from the overall picture we kept the main part of the paper rather short by moving technical and auxiliary material to three appendices: appendix A contains the proof of proposition 3.1. In appendix B we prove Lipschitz continuity of the F-derivative of the forward map within an abstract setting. Therefore, theorem B.2 covers other first order systems as well. The final appendix C includes likewise an auxiliary statement which is nevertheless interesting in its own right: a semi-discrete mapping whose F-derivative is injective and continuous, satisfies the TCC and is locally Lipschitz stable (lemma C.1).

The setting
We consider the acoustic wave equation as a first order system. Let p : , 3}, be the pressure and the velocity field, respectively, where D ⊂ R d is a bounded, connected domain with a piecewise C 1 -boundary. Then, is the P-wave or bulk modulus. We work with c rather than with v p , only to simplify the presentation. The wave equations (2.1) and (2.2) can be written as initial value problem Let us define the space X = L 2 (D) × L 2 (D, R d ) and its subset with ∂D = ∂D D∪ ∂D N . The operator A : D(A) ⊂ X → X is maximal monotone, see the beginning of appendix B for a definition. 3 and see, e.g. [16].
In the remainder of this work we assume the environment to be at rest before we fire the source: (2.8)

The semi-discrete full waveform forward map
is an open set where the sources can be initiated. As we can recover only finitely many degrees of freedom we restrict the parameters of (2.1) and (2.2) to a finite dimensional space. To this end we set where the functions {ϕ j : j = 1, . . . , M} are locally independent over D, that is, if a linear combination vanishes on a nonempty open subset Ω of D then the linear combination must be trivial: Concrete examples for V include: (b) Real-analytic radial basis functions: let ϕ : R d → R be a positive definite and radially symmetric function, see, e.g. [29, chapter 6]. For pairwise different knots ξ j ∈ D, j = 1, . . . , M, the translates {ϕ(· − ξ j ) : j = 1, . . . , M} are linear independent over D. If ϕ is additionally real-analytic then these translates are also locally independent. In fact, any linear combination of these translates is itself an analytic function and as such zero everywhere in D if it vanishes on a nonempty open subset. For instance, the Gaussian ϕ(x) = exp(−γ|x| 2 ), γ > 0, and the multiquadrics ϕ(x) = 1/(1 + |x| 2 ) β , β > 0, have the required properties and are, moreover, positive.
In seismic exploration only part of the wave field can be measured. To model this fact, we introduce the observation (restriction) operator Ψ : The following property of the wave system (3.2) and (3.3) is fundamental for our main result in theorem 4.1 below. Its technical and somewhat lengthy proof is content of appendix A. For its formulation we introduce the space in D, where the infimum is taken over all smooth curves γ in D satisfying γ(α) = x and γ(β) = y. For E ⊂ D open and nonempty, let for all x, y ∈ D.
Hence, in this case, we can replace the condition on T in the proposition above by the stronger condition which is easier to check than (3.5) when there is a lower bound for v p . Basically, T must be large enough so that the wave triggered in Σ reaches the whole of D in the observation period.
We set Φ = Ψ • F. Then, FWI in a semi-discrete acoustic regime consists in finding (c, )

Remark 3.3.
From a numerical point of view, our choice of global basis functions for discretizing c and seems a bit far-fetched. The straightforward approach would be, for instance, to take indicator functions subordinate to the mesh of the used finite element discretization of (2.1) and (2.2). In remark A.3 of appendix A we will address this issue in greater detail.

Injectivity, Lipschitz stability, and tangential cone condition
In a first step towards the TCC we verify injectivity of Φ (c, ). We cast the injectivity problem into an operator framework: in view of (2.8), the map f → (p, v) is linear and we redefine Φ by We right away state our main result of this work, namely that a Lipschitz estimate like (1.3) and a TCC like (1.2) hold for Φ. and Please observe that the TCC (4.4) is slightly stronger than the generic formulation (1.2) insofar as the involved constant decreases with the distance of the arguments (c 1 , 1 ) and (c 2 , 2 ).
The proof of theorem 4.1 is based upon two preparatory results which allow us to apply lemma C.1 of appendix C to Φ. First, we show a local uniqueness result as an immediate consequence of proposition 3.1. W) is an injective mapping and we have that Proof. Assume the minimum to be zero. As Second, we present a continuity result for Φ whose proof is given in appendix B.

Theorem 4.3. We have Lipschitz continuity of the
The involved constant only depends on T, λ − , and λ + .
for one fixed source f in (2.1) and (2.2). In the language of the geophysical community, Φ models a one-shot experiment whereas Φ describes a multi-shot experiment. The local uniqueness result of corollary 4.2 holds accordingly for Φ provided the fired single source f coincides with one of those whose existence for (c, ) is guaranteed by proposition 3.1. Further, theorem 4.1 carries over to the single-source experiment as well (provided the 'right' source is fired).

Conclusion
In this paper we have verified that the full waveform forward operator Φ, see (4.1), satisfies the TCC. Thus, the local convergence of many Newton-like iteration schemes is guaranteed for recovering discrete approximations to the wave speed and density from seismic recordings.
Our achievement should be improved by future research in several ways: (a) As already indicated by remark 3.3, local basis functions have some advantages over global functions for representing wave speed and density. Therefore, it would be desirable to extend the TCC to this setting. In remark A.3 we discuss the resulting challenge in more detail. (b) A further goal is to incorporate more realistic wave propagation models like the elastic or viscoelastic wave equations where the latter equation is the most accurate model right now accounting for attenuation and dispersion. Both models fit in the abstract framework of appendix B, see [17]. Hence, a version of theorem 4.3 holds for them as well. It remains to verify local injectivity via an adapted version of proposition 3.1, which we are confident is valid for the linear elastic model. Moreover, a local Holmgren theorem for viscoelasticity was recently reported in [6]. (c) The Lipschitz estimate (4.3) is local, that is, it holds in a neighborhood of (c + , + ) only. We would like to extend it to an estimate which holds globally in V 2 + since then we would have a global uniqueness result for the semi-discrete seismic inverse problem in the acoustic regime.

Acknowledgments
Funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)-Project-ID 258734477-SB 1173. We thank Roland Schnaubelt for initiating our collaboration and Holger Wendland for discussing some aspects of radial basis functions with us.

Data availability statement
No new data were created or analysed in this study.

Appendix A. Proof of proposition 3.1
The proof of proposition 3.1 is based on the unique continuation and on the propagation of singularities for our system. It will be advantageous to reduce the system (2.1) and (2.2) to two second-order systems with scalar principal part: The principal part of the equation of p is the second-order hyperbolic operator which shows that the second-order system for v has also the second-order hyperbolic operator ∂ 2 t − v 2 p (x)Δ as its principal part. After this preparation we can formulate a precise result on unique continuation for our system, see [8, theorem 1.1]. For that we recall also the definition of the Riemannian metric (3.4).
then the function (p, v) vanishes at 'half time', that is, Moreover, for any T 1 > 0 with T > T 1 + 2 dist(D, E), we have that The proof of this result consist of two distinct parts. At first one proves a local uniqueness theorem (Holmgren uniqueness theorem), i.e. the unique continuation across non-characteristic surfaces. In the second step, one uses geometry to produce a global result, see [20].
Remark A.2. Our approach using the second-order system has the advantage that we could invoke Tataru's uniqueness theorem which proves local uniqueness even for coefficients which are analytic in time and C 1 in space or vice versa [25,26]. However, Tataru's result applies only to scalar second-order operators.
Even though (A.1) is a second-order system, its principal part is a scalar second-order hyperbolic operator. With a small adjustment, Tataru's approach can be applied [7]. Now we turn to the actual proof of proposition 3.1 which is divided into two steps. First, we will show that there exists a forcing term f such that for T > dist(D, Σ), the wave field u = (p, v) does not vanish in (0, T) × Ω. We argue by contradiction. Suppose that u ≡ 0 in (0, T) × Ω. Then, by lemma A.1 there exists ε > 0 such that where λ ∈ C ∞ 0 (0, ε) and g ∈ H 1 (D) with support in Σ, and consider the initial-boundary value problem with initial data p(0, x) = g(x)/c(x), ∂ t p(0, x) = 0, and boundary data The boundary data are inferred from (2.1) and (2.2). This problem has a unique solution p ∈ C([0, ε], H 1 (D)). Furthermore, we define v(t, x) : Then u = ( p, v) satisfies system (2.1) and (2.2) with f = 0. We will use propagation of singularities to establish that u is not zero in (0, ε) × (D\Σ). Indeed, the singularities of the initial data g will travel along the null bicharacteristics of the hyperbolic operator q(x, τ , ξ) = ∂ 2 t − v 2 p (x)Δ which is the principal part of (A.3), see, e.g. [27] or [11, chapter 23]. We point out that this result is true as long as v p ∈ C ∞ (D).
The null bicharacteristics γ : From q(x; τ , ξ) = 0 we infer that τ = ±|ξ|v p (x). Hence, over each point (x, ξ) at t = 0 there are two bicharacteristics. Furthermore, from the ODE we infer that τ (s) = τ for all s and thus, t = 2τ s = ±2|ξ|v p (x)s. So, in both bicharacteristics one can introduce t as a parameter, that is, γ ± (t) = (t, x ± (t), ±v p (x)|ξ|, ξ ± (t)). By the chain rule If (x, ξ) is in the wave front set of g, then the segments of the two null bicharacteristics γ ± in ((0, ε) × D) × (R × R d ), with initial (x, ξ) at t = 0, will be in the wave front set of u. The x component of the bicharacteristic is a geodesic of the metric (3.4). There exist points (x, ξ) ∈ Σ × R d such that at least one of the two bicharacteristics with the initial data (x, ξ) will satisfy x + (t) ∈ D\ Σ or x − (t) ∈ D\ Σ for some 0 < t < ε. This can be seen as follows. Let δ > 0 and let F ⊂ Σ with a smooth boundary such that there exists x ∈ ∂F such that inf z∈D\Σ |x − z| < δ. Choose ξ = ν x . Then the bicharacteristic starting at (x, ξ) satisfies x − (t) ∈ D\ Σ for some t ∈ (0, ε). Note that δ > 0 depends on v p C 1 (D) in view of (A.5) and also on ε > 0.
Suppose now that g is supported in Σ and that its wave front set contains such the point (x, ξ). Then the projection of the wave front set of the solution u = ( p, v) into (0, ε) × D ) must contain the points (t, x − (t)) and thus, the solution cannot vanish in (0, ε) × (D\ Σ).
The solution of (2.1) and (2.2) can now be expressed by Duhamel's principle via Indeed, one computes where we used that λ(0) = 0. Moreover, in view of (A.3) and (A.4) we have that which yield Since u(t, x) is not identically zero for all x ∈ D\ Σ and t ∈ (0, ε), there exists a function λ ∈ C ∞ 0 (0, ε) such that u will not vanish in (0, ε) × (D\ Σ). After reversing the shift in time, this contradicts (A.2) and we have proved that u cannot vanish in (0, ε) × Ω.
In the final step of the proof we validate that the solution (p, v) of (3.2) and (3.3) does not vanish identically on (0, T) × Ω. Assume the contrary. Then, it follows from (3.2) and (3.3) (or, more precisely, from its integrated version (2.7)) that Thus, we must have that Let V loc := span{p j χ D j : j = 1, . . . M} where χ D j denotes the indicator function of D j and p j is a polynomial.
If we now represent c and in V loc , we can prove, by a slight modification of our arguments from above, that for any h ∈ V 2 loc there is a forcing term f and a time T > 0 such that the solution (p, v) of (2.1) and (2.2) does not vanish in (0, T) × (supp h 1 ∪ supp h 2 ). Thus, for each h we can guarantee that at least one of the forcing terms in (3.2) and (3.3) is active. This result is, however, not sufficient to carry over proposition 3.1 (and hence theorem 4.1) to V loc . It remains to show that for each h there is one forcing term f for (2.1) and (2.2) such that the induced forcing terms in (3.2) and (3.3) guarantee (p, v) not to vanish in (0, T) × Ω. We strongly conjecture this to be a fact, unfortunately, we are unable to give rigorous arguments at present.
Even if we succeed, theorem 4.1 might hold only for the multi-shot operator Φ since the applied source f depends on h and one source might not serve all h ∈ V 2 loc . This is then in contrast to the global ansatz functions where we could verify the TCC also for the one-shot operator Φ, see remark 4.4.

Appendix B. A continuity result
In this appendix we verify that Φ : V 2 + ⊂ V 2 → L(V 2 , W) defined in (4.2) is a Lipschitz continuous mapping. We first provide a result for the abstract evolution equation in the spirit of [16].
with u ∈ C ([0, T], X) being the mild (in fact the classical) solution of where u is the classical solution of (B.1) with respect to f. Moreover, F is Lipschitz continuous, that is, The involved constant only depends on T, β − , and β + .
Proof. We can be brief in proving F-differentiability as we will rely on results from [16]. A close inspection of the proofs of lemma 3.3 and theorem 3.6 of [16] yields, for H sufficiently small, that which is the claimed differentiability. which is the claimed Lipschitz continuity.
Finally, we want to emphasize that we can replace the injectivity assumption in the above lemma by Lipschitz stability. Indeed, continuity and (C.1) imply (C.2). Note that Lipschitz stability is known for a variety of semi-discrete inverse problems. We refer, e.g. to [1][2][3][4][5].