Stable reconstruction of simple Riemannian manifolds from unknown interior sources

Consider the geometric inverse problem: There is a set of delta-sources in spacetime that emit waves travelling at unit speed. If we know all the arrival times at the boundary cylinder of the spacetime, can we reconstruct the space, a Riemannian manifold with boundary? With a finite set of sources we can only hope to get an approximate reconstruction, and we indeed provide a discrete metric approximation to the manifold with explicit data-driven error bounds when the manifold is simple. This is the geometrization of a seismological inverse problem where we measure the arrival times on the surface of waves from an unknown number of unknown interior microseismic events at unknown times. The closeness of two metric spaces with a marked boundary is measured by a labeled Gromov--Hausdorff distance. If measurements are done for infinite time and spatially dense sources, our construction produces the true Riemannian manifold and the finite-time approximations converge to it in the metric sense.


Introduction
We study a geometric inverse problem arising from global seismology.Multiple sources go off at unknown times and we measure the arrival times, not knowing which arrivals belong to the same event.The aim is to extract as much information as possible about the planet, modeled here by a simple Riemannian manifold M of dimension n = dim(M ) ≥ 2, and the sources from this information.
Mathematically, the problem boils down to this: Suppose we know the boundary of a Riemannian manifold and for some unknown collection of interior points we know the union of their boundary distance graphs, each shifted with an unknown offset.What can we say about the manifold and these special points from this data?A more detailed description of the data can be found in Section 2 below.
The distances between the source points turn out to be determined exactly by the data (Theorem 2), but there is no hope of reconstructing the manifold if the source set is finite.However, the source set can be regarded as a discrete metric space, and this space approximates the Riemannian manifold.It is not only close as a metric space, but we can also assign approximate boundary points with good accuracy (Theorem 9).To describe the similarity of two metric spaces "with the same boundary", we define a labeled Gromov-Hausdorff distance.This extends the classical Gromov-Hausdorff distance and compares both the similarity of the metric spaces and the sameness of the boundaries -with a fixed model space for the boundary.It is in the sense of this distance that our approximation is good, and the quality of the approximation can be estimated explicitly directly from data.
The sources can be point sources in space time produced by a suitable stochastic process (Proposition 12), for example.We present a method to look at the data on the boundary without knowing the manifold and reconstruct an approximate manifold -and get an explicit error bound on our reconstruction.This increases the applicability of geometric inverse problems to the messy cases with real data, but we do not attempt to optimize the bounds.When it comes to explicit geometric estimates, this paper should be seen as a proof of concept.
We assume the source events to be discrete in the spacetime, but their projections to the space can be dense.In this case the approximate reconstruction is actually perfect, and even the smooth and Riemannian structure is determined (Theorem 1).If measurements are made for finite but increasing time, the approximate reconstructions converge to the true manifold in the labeled Gromov-Hausdorff distance if the sources are spatially dense (Theorem 11).
Estimation of the quality of the approximation from data alone is lossy.To prove that it is not hopelessly so, we prove that by making the source points dense enough the observed density can be made arbitrarily small.There is an explicit estimate between the true and observed density in both directions (Theorems 1 and 10).
We solve a problem stemming from physics by developing geometry rather than by applying ad hoc tricks.These developments include the introduction of the labeled Gromov-Hausdorff distance and various quantitative descriptions of simplicity of a Riemannian manifold.
From a more applied point of view, the key advances are that we introduce of a new and more physically relevant problem and provide quantitative and fully data-driven error estimates.
The precise definitions and statements of results are given in Section 2. We will discuss the results and their context in Section 3. The overall strategy of proofs is given in Section 4 and the actual proofs follow in subsequent sections.
Acknowledgements.MVdH was supported by the Simons Foundation under the MATH + X program, the National Science Foundation under grant DMS-1815143, and the corporate members of the Geo-Mathematical Imaging Group at Rice University.JI was supported by the Academy of Finland (projects 332890 and 336254).ML was supported by Academy of Finland (projects 284715 and 303754).TS was supported by the Simons Foundation under the MATH + X program and the corporate members of the Geo-Mathematical Imaging Group at Rice University.The authors want to thank Peter Caday and Vitaly Katsnelson for useful discussions.

Definitions and theorems
2.1.Arrival time data.Let M be a Riemannian manifold with boundary.We denote by d : M × M → R its Riemannian distance.The spacetime M × R comes with two natural projections, π onto M and τ onto R.
Let S ⊂ M × R be a set of sources.The source s ∈ S goes off at the point π(s) at the time τ (s).The data measured is the arrival times of signals at the boundary, given as the set Q(S) = {(x, τ (s) + d(x, π(s)); x ∈ ∂M, s ∈ S}. ( We emphasize that we know the set Q(S) but not the set S. We only observe the data as this point set without any labels to tell which arrivals correspond to the same sources.The arrival time function of a source s ∈ S is the function a s : ∂M → R given by a s (x) = τ (s) + d(π(s), x).Let us denote the graph of the function a s by G(a s ) ⊂ ∂M × R. The data set can be rewritten as Q(S) = s∈S G(a s ).
Let φ : ∂M 1 → ∂M 2 be a diffeomorphism between two boundaries.If A ⊂ ∂M 2 × R, we define The set A would typically be the data set Q(S 2 ) or a subset thereof.When there are two manifolds, we decorate all the objects with the subscripts 1 and 2 when clarity requires so.
A subset of the spacetime M ×R is said to be discrete when it has no accumulation points.If M is compact, it follows that any time slice M ×[a, b] contains only finitely many points of a discrete set.
Our results apply to so-called simple manifolds.A Riemannian manifold is called simple if it is compact with strictly convex boundary (in the sense of definiteness of the second fundamental form) and the exponential map at every point is a diffeomorphism on its maximal domain of definition.Consequently there are no conjugate points and any two points on the manifold, including the boundary, are joined by a unique geodesic depending smoothly on the endpoints.
A key concept in our analysis is density.We say that a subset A ⊂ X of a metric space X is ε-dense if the balls B(a, ε) with a ∈ A cover X. Equivalently, every point in X has a point of A at a distance below ε.
We set out to study whether and how a Riemannian manifold M can be determined by the set Q(S) and prior knowledge of the boundary ∂M .
2.2.Precise determination.We begin by studying what can be determined from the data with no error.
The source sets S i ⊂ M i × R can be finite or infinite.However, the discreteness assumption implies that every S i ∩ (M i × [−T, T ]) is finite and so the source sets are at most countably infinite.Although the sources are not dense in the spacetime, the source points in space can be dense.In that case the same data determines the whole manifold up a diffeomorphism.Theorem 1.Let M 1 and M 2 be simple Riemannian manifolds so that the diameter of both manifolds is at most C diam and the sectional curvature is bounded from above by C sec+ > 0 so that C diam C sec+ < π.
(3) Let φ : ∂M 1 → ∂M 2 be a Riemannian isometry.Suppose the two sets S i ⊂ int(M i )× R are discrete in M i × R. Let Q(S i ) ⊂ ∂M i × R be defined as above.
The condition (3) in most of our theorems always holds when the manifold is a priori known to have non-positive sectional curvature.See Section 3.2.4 for details.
The main novelty of Theorem 1 is in allowing the data set Q(S) to be a union rather than an indexed collection of graphs of the arrival time functions.The rest of our results do not have a similar precedent, and our method of proof for Theorem 1 is new and constructive.See the discussion in Section 3 below for details.
The next result concerns a situation with a finite amount of sources.Informally, the next theorem states that knowledge of ∂M as a smooth manifold and the set Q(S) determines the distances between the sources, the time differences, and two kinds of distance differences at the boundary.More formally: Theorem 2. Let M 1 and M 2 be simple Riemannian manifolds and φ : there is a bijection ξ : ) for all x ∈ ∂M 1 and s, r ∈ S 1 , and x, y ∈ ∂M 1 and s ∈ S 1 .
While measuring for infinite time gives full uniqueness, it is not realistic.Therefore we consider measurements for a finite but increasing time and show that there is an approximate reconstruction and in the limit of infinite time it tends to the correct one in a suitable sense.To be able to state the results, we need to set up a way to compare the true manifold to a discrete approximation.Before pursuing this direction, we record a result that helps clean up the data by disentangling the union of graphs into separate graphs.The arrival time functions a s are easily verified to be smooth when π(s) / ∈ ∂M and the manifold is simple.
Proposition 3. Let M be a simple Riemannian manifold.Suppose S ⊂ int(M )×R is discrete in M × R. Then the set Q(S) and the Riemannian structure of ∂M determine the set determines the set of connected components of the sets G(a s ) ∩ Ω, s ∈ S.

2.3.
Labeled Gromov-Hausdorff distance.This subsection is devoted to the metric geometry we will use to compare our approximate reconstructions to the true manifold.Definition 4. Let X and Y be compact metric spaces and L any set (understood as a set of labels).Let α : L → X and β : L → Y be any functions.We define the labeled Gromov-Hausdorff distance between (X, α) and (Y, β) to be Here d Z H is the Hausdorff distance on the metric space Z.If the set of labels is empty, then this distance reduces to the usual Gromov-Hausdorff distance; in this context we understand the supremum of an empty set to be zero.
We will be comparing our discrete reconstruction to the manifold M with labels given by the known set L = ∂M , labeled by the inclusion ι : ∂M → M .To make it meaningful to state that the labeled Gromov-Hausdorff distance is small, we must ensure that this is a well-behaved concept of distance.
The point of this definition is that in addition to getting a good metric approximation of the manifold in the sense of Gromov-Hausdorff distance, we want to get knowledge of how well the known boundary ∂M sits inside the approximation.This is what the concept is designed to do.Proposition 5.The labeled Gromov-Hausdorff distance is symmetric and satisfies the triangle inequality.Moreover, if and only if there is an isometry h : X → Y so that h • α = β.
2.4.Quantitative simplicity.We will make approximate reconstructions of the whole manifold, and we can quantify the error of the reconstruction precisely once we have a priori bounds on curvature, diameter, and similar geometric properties.Definition 6.We say that a manifold M with boundary has bounded geometry with constants The diameter of M is at most C diam .
(2) The sectional curvature is always in the range [−C sec− , C sec+ ].
(3) For any x ∈ M and η 1 , η 2 ∈ T x M we have as long as the exponentials are defined.
(4) Every Jacobi field J with J(0) = 0 along any unit speed geodesic satisfies for all t > 0 for which the geodesic is defined.(5) If h 1 and h 2 are, respectively, the first and second fundamental forms of the boundary ∂M , then The distances on M and along ∂M satisfy for all x, y ∈ ∂M .(7) Take any unit speed geodesic γ on M and let ρ be the distance function ρ(x) = d(x, γ(0)).For any t > 0, let w ∈ T γ(t) M be any vector orthogonal to γ(t).Then the Hessian H ρ of ρ satisfies We will derive a number of other estimates from these with constants depending on the ones listed here.See Section 9 for details.The constants appearing in the derived estimates are all given in Appendix B. We will refer to the appendix whenever a new constant is introduced in a claim.All claims are given with simple constants, so there will be a number of relations between the various constants as is clear from the appendix.
It should be noted in part 7 of the definition that w is a vector and H ρ w is a covector.The inner product is their duality pairing, which we find most convenient to write as w, H ρ w instead of (H ρ w)(w) or H ρ (w, w).
Having these estimates and constants is not an added assumptions, but merely a quantification of the simplicity of the manifold M .Our next result justifies that the constants of Definition 6 can well be called constants of quantitative simplicity.Proposition 7. Every simple manifold has bounded geometry in the sense of Definition 6 for some constants.
2.5.Approximate determination.Theorems 2 and 1 state what can be determined exactly from our data.We now turn to studying how well the whole manifold can be reconstructed with finite data.
The source set is as before.The measurements at the boundary start at time t = 0 and continue up to time t = T .From this finite amount of data we can then draw an approximate conclusion about the geometry of our manifold.The error can be quantified very concretely in terms of the geometric bounds.
The only points we know on the manifold are the finite number of sources that have produced a boundary signal on the time interval [0, T ] and the a priori known boundary itself.Therefore what we have is a discrete approximation of the smooth Riemannian manifold.In addition to getting a metric approximation of the manifold itself, we need to approximately embed the known boundary ∂M into the approximation, and so we will use the Gromov-Hausdorff distance labeled by ∂M .
Let us denote the set of spatial source points by P = π(S).To each p ∈ P is associated a boundary distance function r p : ∂M → R given by r p (x) = d(x, p).Let c(p) ⊂ ∂M be the set of critical points of this function.By compactness #c(p) ≥ 2.
Consider a simple manifold of bounded geometry in the sense of Definition 6.At x ∈ c(p) let λ(p, x) be the smallest eigenvalue of the Hessian of r p at x.This requires that ∂M is known as a Riemannian manifold, not just as a smooth one.For y ∈ c(p) we define and for all x ∈ ∂M E Finally, we let These three types of Es will be all determined by the data and the constants of bounded geometry.Observe that if the Hessians of the distance functions are never large enough, all these Es may well be infinite.For any p, q ∈ P and r, s > 0 we define the lentil The lentil is said to have thickness One can check that on a simple manifold the thickness is the length of the segment of the geodesic γ p,q between the two points that lies within the lentil, provided that r, s < d(p, q).To see this, observe that a point γ p,q (t) is in the lentil if and only if t < r and d(p, q) − t < s.
Our inversion procedure starts with preprocessing the data.Proposition 3 gives on aspect of it, and the next one says that the derived quantities introduced above are determined by the data.Proposition 8. Let M 1 and M 2 be simple Riemannian manifolds.Let φ : ∂M 1 → ∂M 2 be a smooth isometry.Suppose that the two sets R be defined as above.Let ξ : S 1 → S 2 be the bijection of Theorem 2 and let Ξ : P 1 → P 2 be the map for which Ξ • π 1 = π 2 • ξ.Then: (1) The critical point sets correspond via c(Ξ(p)) = φ(c(p)) for all p ∈ P 1 .The boundary distance quantities correspond via E 1 (p, y) = E 2 (Ξ(p), φ(y)) for all p ∈ P and y ∈ c(p).Moreover and E 1 (x) = E 2 (φ(x)) for all x ∈ ∂M 1 and E 1 = E 2 .
(2) If we fix any ε 1 > 0 and define then Γ 2 = Ξ(Γ 1 ).(3) Fix any ε 1 > 0 and δ > 0. Define Γ i for i = 1, 2 as above.Take any x 1 , y 1 ∈ Γ 1 so that d(x 1 , y 1 ) > δ, and pick any r ∈ (δ, d 1 (x 1 , y 1 )).Set s = d 1 (x 1 , y 1 ) − r + δ and let where x 2 = Ξ(x 1 ) and y 2 = Ξ(y 1 ).Then and so Now that we know the data to determine the relevant auxiliary quantities uniquely, we move to working on a single manifold.The next theorem states that under known bounds on the geometry and metrically known boundary one can estimate directly, from the data, how dense the source points are.For our purposes this is to assign approximate boundary points, and estimate the metric distance of an approximate reconstruction from the true manifold.
In practice, the set P in the next theorem stands for the spatial source set π(S) as above, but the claim is valid for any subset of int(M ).Theorem 9. Let M be simple a Riemannian manifold with the constants of Definition 6 (by Proposition 7) satisfying the condition (3), and let P ⊂ int(M ) be any set.This theorem uses the constants given in equations (154), (156), (158), and (160).
Fix any ε 1 > 0 and define the set Γ ⊂ P as in (17).Set ε 2 = ε 1 + E and δ = C 1 ε 2 .Suppose that for every x, y ∈ Γ with d(x, y) > δ and r ∈ (δ, d(x, y)) we have L x,y r,d(x,y)−r+δ ∩ P = ∅. Then: (1) (2) Given any α : ∂M → P so that for all x ∈ ∂M , the pair (P, α) is close to (M, ι), with the inclusion ι : ∂M → M , in the sense of the labelled Gromov-Hausdorff distance: A map α satisfying (23) will be constructed explicitly from the data.This is explained in the proof of the theorem in Section 9.4.
The first part of this theorem says that under an additional curvature assumption the quantities defined by the data (due to Proposition 8) give an explicit estimate for the density of the source points.The second part is an estimate of the labeled Gromov-Hausdorff distance between the approximating discrete space P with the approximate boundary described by α and the true manifold M with boundary ∂M .The key assumption is that all lentils of thickness δ, which is proportional to ε 2 , contain a source point, and by Proposition 8 this is determined uniquely by the data.
The quality of the estimates given in claim 1 of Theorem 9 depend on the choice of ε 1 .When faced with a practical reconstruction task, one ought to vary ε 1 and see which choice gives the best error bound ε 2 while still having non-empty lentils.The only parameter we choose is ε 1 ; the second parameter ε 2 which plays a more central role is based on the observed quantity E of (14).Heuristically, ε 1 describes how close the chosen source points are to the boundary and ε 2 how far from any given boundary point these near-boundary sources are.These constants are then used to establish estimates for interior density of the sources.
The estimates of claim 1 of Theorem 9 are vacuous if ε 2 = ∞.For the theorem to be meaningful, we need to ensure that for some ε 1 we indeed have ε 2 < ∞.To be able to get an arbitrarily good estimate for density and metric distance, we also need to ensure that ε 2 will become as small as we like if the source set P ⊂ M is dense enough.This is a reverse density estimate: If the source set is very dense, it will appear dense by our estimates.
and ε 1 = C 5 ε, then: (1) The quantity ε 2 = ε 1 + E as used in Theorem 9 satisfies (2) Denote δ = C 1 ε 2 .For any x, y ∈ P that are more than distance δ apart and any r ∈ (δ, d(x, y)) the lentil L x,y r,d(x,y)−r+δ contains an element of P .(3) The set P satisfies the assumptions of Theorem 9 for the given choice of ε 1 ; especially the lentil intersection property (21) holds.Therefore one can conclude from the data set Q(S) and the geometric constants alone that P is ε-dense and the discrete approximation is ε-good.
Choosing ε 1 correctly is important for Theorem 10.If ε 1 is too large, the resulting estimates are not tight enough.If it is too small, then too few boundary sources are identified and not all lentils contain source points.
2.6.Convergence of approximations.The next theorem treats the case of gradually acquired and interpreted data.Assuming the sources are eventually dense, a finite amount of data determines an approximate manifold.The data itself gives an estimate for how good the approximation is and the approximations converge to the true manifold with boundary in the labeled Gromov-Hausdorff sense.
Theorem 11.Let M be a simple Riemannian manifold with geometry bounded in the sense of Definition 6 satisfying the condition (3).Let S ⊂ int(M ) × R be a countable discrete set of sources so that the projection π(S ∩ (M × [0, ∞))) ⊂ M is dense.
Then for each T > 0 the boundary ∂M as a Riemannian manifold and the set determine a finite metric space M T and a map α T : ∂M → M T so that lim Thus (M T , α T ) → (M, ι) as T → ∞ in the labeled Gromov-Hausdorff sense.
Finally, we point out that if the geometric assumptions of Theorem 11 are satisfied, then the assumptions on the sources are almost certainly satisfied in a simple stochastic model.That is, for almost any source pattern a finite amount of data determines an approximation and the approximation converges to the true solution with increasing data in a manner quantifiable from data itself.
We say that two measures are uniformly comparable when they are absolutely continuous with respect to each other and the Radon-Nikodym derivatives both ways are bounded.Let ν be the product of the volume measure of M and the Lebesgue measure of R, and let µ be a measure on M × R uniformly comparable with ν.A homogeneous Poisson point process on the spacetime M ×R with intensity λ > 0 is such a point process that the number of points in the set A ⊂ M × R follows the Poisson distribution with parameter λμ(A).Proposition 12. Let M be a simple Riemannian manifold with bounded geometry satisfying (3).Let µ be a measure on M ×R which is uniformly comparable with the natural product measure ν.Suppose the set S ⊂ M × R is given by a homogeneous Poisson point process with intensity λ > 0 and the measure µ.
That is, if the events S ⊂ M × R are given by a Poisson process, then almost surely the reconstruction from data for finite but increasing time converges to the true solution as time increases.

Discussion
3.1.Geometric inverse problems and seismology.An elastic body -e.g. a planet -can be modeled as a manifold, where distance is measured in travel time: The distance between two points is the shortest time it takes for a wave to go from one point to the other.If the material is isotropic or elliptically anisotropic, then this elastic geometry is Riemannian.However, this sets a very stringent assumption on the stiffness tensor describing the elastic system, and Riemannian geometry is therefore insufficient to describe the propagation of seismic waves in the Earth.If no structural assumptions on the stiffness tensor apart from the physically necessary symmetry and positivity properties are made, this leads necessarily modeling the planet by a Finsler manifold as was explained in [13].
An isotropically elastic medium carries pressure (P ) and shear (S ) wave speeds that are conformally Euclidean metrics.Of these two the P -waves are faster [8].In order to be true to the isotropic elasticity we should measure both P -and S -wave arrivals.As we have many sources going off in different locations and we do not know which arrival is related to which source we cannot a priori catalog arrivals of different wave types.We simplify this aspect of the problem by disregarding polarizations and considering only one type of wave.
In the two following subsections we review some seismologically relevant geometric inverse problems on Riemannian and Finsler manifolds.We also draw relations between the mathematical data, considered in this paper, and some real world seismic measurements.In the final subsection we discuss our contributions to field of inverse problems and possible future research directions suggested by the results of this paper.
3.1.1.Geometric inverse problems on Riemannian manifolds.The problem of determination of the isometry type of a compact Riemannian manifold from its boundary distance data was introduced for the first time in [30].The reconstruction of the smooth atlas on the manifold and the metric tensor in these coordinates was originally considered in [27].In contrast to the paper at hand these earlier results do not need any extra assumption for the geometry, but have a complete data in the sense that the distance function r p to the boundary is known for any interior source point p ∈ int(M ).
The problem of boundary distance data is related to many other geometric inverse problems.For instance, it is a crucial step in proving uniqueness for Gel'fand's inverse boundary spectral problem [27].Gel'fand's problem concerns the question whether the data (∂M, (λ determine (M, g) up to isometry, when (λ j , φ j ) are the Dirichlet eigenvalues and the corresponding L 2 -orthonormal eigenfunctions of the Laplace-Beltrami operator.
Belishev and Kurylev provide an affirmative answer to this problem in [4].
In [28] the authors studied a question of approximating a Riemannian manifold under the assumption: For a finite set of receivers R ⊂ ∂M one can measure the travel times {r p (z); z ∈ R} for any point p ∈ int(M ).In this problem arrival times (with known initial times) are measured from all interior sources at a finite number of boundary points, whereas in the present paper we measure arrival times (with unknown initial times) from a finite number of interior sources at all boundary points.The authors construct an approximate finite metric space M ε , with the same cardinality as the receiver set R, and show that if the Hausdorff distance of R and ∂M is less than ε > 0 then the Gromov-Hausdorff distance of M and M ε is proportional to some positive power of ε.In contrast to the present paper an independent travel time measurement was made for any interior source point.
We recall that for a source s ∈ S we do not know the initial time τ (s), but due to Proposition 3 we can recover the arrival times a s (z) = d(z, π(s)) + τ (s) for each z ∈ ∂M .Taking the difference of the arrival times one obtains a boundary distance difference function for all z 1 , z 2 ∈ ∂M , and this is independent of the initial time.In [36] it is shown that if U ⊂ N is a compact subset of a closed Riemannian manifold (N, g) and int(U ) = ∅, then distance difference data determine (N, g) up to an isometry.This result was recently generalized to complete Riemannian manifolds [25] and for Riemannian manifolds with boundary [15,24].These results require sources at all interior points, unlike the ones in the present paper.
If the sign in the definition of the distance difference functions is changed, we arrive at the distance sum functions, for all x ∈ M and z 1 , z 2 ∈ ∂M .These functions give the lengths of the broken geodesics, that is, the union of the shortest geodesics connecting z 1 to x and the shortest geodesics connecting x to z 2 .Also, the gradients of D + x (z 1 , z 2 ) with respect to z 1 and z 2 give the velocity vectors of these geodesics.The inverse problem of determining the manifold (M, g) from the broken geodesic data, consisting of the initial and the final points and directions, and the total length, of the broken geodesics, has been considered in [31].In [31] the authors show that broken geodesic data determine the boundary distance data and use then the results of [27,30] to prove that the broken geodesic data determine the Riemannian manifold up to an isometry.A different variant of broken geodesic data was recently considered in [40].
The Riemannian wave operator is a globally hyperbolic linear partial differential operator of real principal type.Therefore, the Riemannian distance function and the propagation of a singularity initiated by a point source in space time are related to one another.We let u be the solution of the Riemannian wave equation with a point source s ∈ S. In [17,19] it is shown that the image, Λ, of the wavefront set of u, under the musical isomorphism with the image of the unit sphere S π(s) M at π(s) under the geodesic flow of g.Thus Λ ∩ ∂(SM ), where SM is the unit sphere bundle of (M, g), coincides with the exit directions of geodesics emitted from π(s).In [37] the authors show that if (M, g) is a compact smooth non-trapping Riemannian manifold with smooth strictly convex boundary, then generically the scattering data of point sources (∂M, R ∂M (M )) determine (M, g) up to isometry.Here, R ∂M (x) ∈ R ∂M (M ) for x ∈ M stands for the collection of tangential components to boundary of exit directions of geodesics from x to ∂M .
If M is an open subset of a complete Riemannian manifold (N, g), then one more important data set, given by the wavefront set of u, is related to a generalized spheres of radius r > 0, that is given by the formula and exp p is not singular at v}.
In [11] the authors show that the spherical surface data determine the universal cover space of N .If a generalized sphere S(p, r) is given the authors show that there exists a specific coordinate structure in a neighborhood of any maximal normal geodesic to S(p, r) such that in these coordinates metric tensor g can can be determined.However this does not determine g globally.The authors provide an example of two different metric tensors which produce the same spherical surface data.
A classical geometric inverse problem, that is closely related to the distance functions, asks: Does the Dirichlet-to-Neumann mapping of a Riemannian wave operator determine a Riemannian manifold up to isometry?For the full boundary data this problem was solved originally in [4] using the Boundary control method.Partial boundary data questions have been studied for instance in [35,43].Recently [33] extended these results for connection Laplacians.Lately also inverse problems related to non-linear hyperbolic equations have been studied extensively [32,38,57].For a review of inverse boundary value problems for partial differential equations see [34,55].
Maybe the most studied geometric inverse problem formulated with the distance functions is the Boundary rigidity problem.This problem asks: Does the boundary distance function, that gives a distance between any two boundary points, determine (M, g) up to an isometry?In an affirmative case (M, g) is said to be boundary rigid.For a general Riemannian manifold the problem is false: Suppose the manifold contains a domain with very slow wave speed, such that all the geodesics starting and ending at the boundary avoid this domain.Then in this domain one can perturb the metric in such a way that the boundary distance function does not change.It was conjectured in [42] that for all compact simple Riemannian manifolds the answer is affirmative.In two dimensions it was solved in [47].For higher dimensional case the problem is still open, but different variations of it has been considered for instance in [7,9,52,53].
A general feature of geometric inverse problems related to seismology is that each point source is analyzed in isolation and there is a source in every point on the manifold.We drop both of these assumptions to step towards a physically more accurate model.The source set can be finite and the data does not a priori make a distinction between different events.This distinction is achieved by Proposition 3. Another novelty is that our discrete approximation is quantitatively stable in the sense of the labeled Gromov-Hausdorff distance.All our estimates are given by the data and some geometric constants and do not depend on other a priori knowledge of the manifold M .
3.1.2.Related geometric inverse problems on Finsler manifolds.In [13] we studied the recovery of a compact Finsler manifold from its boundary distance data.In contrast to earlier Riemannian results [27,30] the data only determines the topological and smooth structures, but not the global geometry.The Finsler function F : T M → [0, ∞) can be however recovered in a closure of the set G(M, F ) ⊂ T M , that consists of points (p, v) ∈ T M such that the corresponding geodesic γ p,v is distance minimizing to the terminal boundary point.We also showed that if the set T M \ G(M, F ) is non-empty then any small perturbation of F in this set leads to a Finsler metric whose boundary distance data agrees with the one of F .If G(M, F ) = T M , then the boundary distance data determines (M, F ) up to a Finsler isometry.For instance the isometry class of any simple Finsler manifold is determined by this data.The same is not true if only the boundary distance function is known [23].Thus a simple Finsler manifold is never boundary rigid.In [14] we utilized the main result of [13] and generalized the result of [31], about the broken geodesic data, on reversible Finsler manifolds, satisfying a convex foliation condition.
Although simple Finsler manifolds are not boundary rigid there are results considering their rigidity questions for some special Finsler metrics.For instance it was shown in [44] that Randers metrics F i = G i + β i indexed with i ∈ {1, 2} with simple and boundary rigid Riemannian norm G i (x, v) = g ij (x)v i v j and closed one-form β i , have the same boundary distance function if and only if G 1 = Ψ * G 2 for some boundary fixing diffeomorphism Ψ : M → M and β 1 − β 2 = dφ for some smooth function φ vanishing on ∂M .It is worth of mentioning that analogous results have been presented earlier on a Riemannian manifold in the presence of a magnetic field [10,2].
The sphere data described in connection to the Riemannian results of [11] above has also been studied on Finsler manifolds [12].Knowledge of the spheres uniquely determine the fundamental tensor and the curvature operator along any geodesic passing through the known domain, but in contrast to Riemannian geometry this information is insufficient for a full reconstruction of the universal cover.
3.1.3.Geophysical literature.In fact, the inverse problem considered here can be directly related to seismology.The sources correspond with microseismic events with unknown locations and origin times while the metric corresponds with the wave speed.Here, we consider only one wave speed associated with elastic Pwaves, but incorporating a second wave speed associated with elastic S -waves is quite straightforward.In the past decade there has been extensive research on the joint recovery of wave speed and event locations, with mining, geothermal and hydraulic-fracturing applications of induced seismicity but also in studies of the crustal structure.Here, we present a comprehensive analysis providing fundamental insight in the feasibility of succeeding in this, while focusing on recovering the wave speed.We use geometric data rather than wave fields; in applications, the extraction of these have been routine, for example, with template matching [18].
Initial empirical studies [26,5] assumed a simple wave speed model varying in one coordinate only and one or two strings of receivers (in boreholes) penetrating the manifold.Strategies have been broadly based on optimization, such as variations of gradient descent [49,21], with intertwining updating [41] and a neighborhood algorithm [54] employing different optimization criteria [58,51] or via intermediate (approximate) interior wave-field recovery [50].A basic statistical, Bayesian framework has been developed in parallel [60].

Discussion of technical assumptions.
In this subsection we present a discussion of many technical assumptions we made in our theorems.Our focus here is on the key phenomenology of having multiple sources with unknown times and we have chosen to simplify other aspects of the problem.

Density of source points.
We assume in all our results that the sources are discrete in the spacetime M × R.This is physically reasonable and also practical.If the sources were dense or even accumulated somewhere, Proposition 3 would become far more complicated.
However, there is little hope of full uniqueness with finitely many source points.The sources can be discrete in spacetime but dense in space, and this is indeed the setting of Theorems 1 and 11 in which the Riemannian manifold is determined uniquely either from full time data or asymptotically from increasing time data.
3.2.2.Full data and convexity.All our results concern full data in the sense that arrivals are recorded on the whole boundary ∂M .The full boundary is crucial for Proposition 3 which disentangles the data into a collection of graphs, and we make so heavy use of differential tools on the boundary throughout the paper that a discrete subset of the boundary is beyond current reach.
Let us construct an explicit example of a surface M and a subset Γ ⊂ ∂M so that our results fail with data recorded only on Γ.Every pair of points on a smooth compact Riemannian manifold with boundary is connected by a C 1 -smooth distance minimizing curve [1].We choose our a manifold to be the horseshoe-shaped domain of Figure 1.
A domain where partial data is insufficient.We split the domain M into two pieces M 1 and M 2 with respect to the line (red dotted line) that is normal to ∂M at x 0 ∈ ∂M (blue dot).Then we choose a domain Γ ⊂ ∂M 1 (red arch) so that any minimizing curve joining a point on Γ and a point in M 2 touches the boundary near x 0 .The curve P ⊂ M is any involute of the boundary, meaning that the distance from all points on P to x 0 is the same.
Because P is an involute of the boundary as in the figure, d(x, p) = d(x, q) for any x ∈ Γ and p, q ∈ P .Therefore from the point of view of our data, the set P appears to collapse to a point.
Worse, if Γ happens to be an involute as well, then all the distance functions to points p ∈ M 2 are constants in Γ.As the unknown origin times produce unknown constant offsets to the boundary distances, all points in M 2 look alike when seen from an involutive Γ in this sense.
The problems of Figure 1 also illustrate the problems lack of convexity of the boundary may cause.
Similar problems may arise in higher dimensions if data is recorded on too small sets.Let the manifold M be the closed unit ball in R n and consider partial data on Γ = ∂M ∩ H, where H ⊂ R n is a hyperplane through the origin.Let p 1 and p 2 be two points in int(M ) \ H situated symmetrically about H. Then the boundary distance functions of p 1 and p 2 agree on Γ.The sets M and Γ ⊂ ∂M have a reflection symmetry which leaves the data invariant, making it impossible to distinguish the two sides.
3.2.3.Infinitely many lentils.The crucial condition for the estimates in Theorem 9 was that all of the lentils meet the spatial source set P .There is an infinite family of lentils, so we have a large number of conditions to check from our data.The lentils are open and cover a certain compact subset of M (excluding a layer near the boundary), so in fact using a finite number will suffice.However, it is not easy to identify a covering collection of lentils -for which one would then easily check whether they contain source points -from data.
3.2.4.Conjugate points and boundary sources.If there is a source point on the boundary, the graph of the distance function is singular at that point.The separation of the set Q(S) into graphs in Proposition 3 relies on smoothness.If several corners happen to coincide, especially if n = 2, it may be difficult to choose how to continue the graphs correctly.
This problem becomes far worse if there are conjugate points.They also cause the boundary distance function to be non-smooth but can do so for several different points.This makes both disentangling the data into graphs and the analysis of those graphs substantially more complicated.
The condition (3) on the constants of Definition 6 is assumed in Theorems 1, 9, 10, and 11.This condition is used in Lemma 24 and Proposition 30 when the manifold is compared to model manifold with constant sectional curvature.For the comparison to work, we need the model manifold of the same diameter to not have conjugate points, and this is exactly what the condition ensures.
The condition (3) always holds in negative curvature.If the manifold is known to have non-positive sectional curvature and an explicit bound on the diameter, one can simply choose the sectional curvature upper bound C sec+ > 0 to be small enough to satisfy (3).The condition also holds on all simple manifolds of constant sectional curvature.
The constants C a of (151), C b of (152), and C c of (153), as well as all constants derived from them, become worse when C diam C sec+ gets close to π.

3.2.5.
Constants of quantitative simplicity.The constants of quantitative simplicity have to be the same for both manifolds M i in Theorem 9.The set Γ i and the approximate boundary inclusion map α depend on these constants.The distances between the source points are determined irrespective of the constants by Theorem 2, but the constants have an effect on the approximate boundary structure and the estimates on the quality of the reconstruction.
3.2.6.Scaling of small quantities.Consider the case when all the various epsilons are very small.In the setting of Theorem 9 we have δ ≈ ε 2 ≈ ε 2 and in Theorem 10 we have δ ≈ ε 2 ε 2 ≈ ε.Ideally, we would expect to see no second powers so that all concepts of density -the size of lentils, the true density, the observed quantity ε 2 , the density estimated from data -to be bi-Lipschitz equivalent to each other.The difference in scaling is all due to the lentils being shaped so that This scaling is easy to verify in Euclidean geometry, and the relevant Riemannian aspects are covered in Lemmas 24 and 25, and Proposition 30.
3.2.7.Generalization to simple Berwald manifolds.As discussed above, our methods depend on the smoothness of the distance function.Thus these techniques are not applicable for general Riemannian metrics and relaxing the simplicity assumption would likely require different methods.There is however hope to apply these techniques to study analogous questions in the case of certain simple Finsler metrics.
Any constant speed geodesic γ of a smooth connected and complete Finsler manifold (N, F ) is given by the solution of the geodesic equations, in local coordinates, γk (t) + Γ k ij (γ, γ) γi (t) γj (t) = 0.Here Γ k ij (x, v) are the coefficients of the Chern connection on the tangent bundle [3].Due to the dependence of the directional variable v ∈ T x N a reverse curve of a Finsler geodesic does not need to be a geodesic.If the Chern connection coefficients are directionally independent the Finsler metric F is called Berwald.Thus Berwald manifolds have a well defined canonical covariant differential operator for tensor fields of any order.Therefore Berwald metrics carry a Gauss type of formula to connect boundary and interior Hessians for the distance function.This is not true for a general Finsler metrics and the lack of Gauss formula would cause issues in several parts of our proofs.It is worth of mentioning that Berwald metrics are not just a theoretical curiosity but have a connection to linear elasticity.It has been observed that elastic Finsler metrics, arising from transverse isotropic medium in weak anisotropy, are actually Berwaldian [59].

Plan of proofs
Our data is given in the form of a set, a union of graphs, and we begin by disentangling it into the separate graphs.We prove Proposition 3 to this effect in Section 5.
In Section 6 we turn to finding distances between source points and other information that the data gives exactly.With the help of Proposition 3, we may start with the knowledge of the arrival time functions a s .Each graph corresponds to a unique source point, giving us the first claim.Computing suitable differences of arrival times, especially along the unique geodesic between two given source points, we cancel out unwanted contributions and find the distance between the two points.Once we have these first two claims, the rest of Theorem 2 follows straightforwardly.
In Section 7 we depart from differential geometry to the geometry of metric spaces and study the labeled Gromov-Hausdorff distance in more detail, including the proof of Proposition 5.
The most substantial part of our proofs is in estimating the density of the sources, and we will do this in Section 9 using the explicit bounds on geometric quantities found in Section 8 -including the proof of Proposition 7.
The first task is to estimate the density of sources near the boundary, which practically amounts to estimating how near the boundary a given source point is.This is done in terms of the curvature of the graph of the arrival time function; when the point is very close to the boundary, the Hessian of the said function blows up.We will prove that d(p, x) ≤ E(p, x) for any source point p ∈ P and a critical point x ∈ ∂M of its arrival time function.This gives an estimate in terms of data for how close to any given boundary point must there be a source point.The initial data is first processed into the form of E(p, x) and quantities derived from it, and we verify in Proposition 8 that the necessary auxiliary quantities are indeed determined by our data.
To get a density estimate deep inside the manifold, a different set of tools must be employed.A key tool is the lentil, an intersection of two metric balls.Whether a source point belongs in a lentil defined by two other source points is entirely decided by the distances between the three points, and by Theorem 2 this information is indeed determined by the data.We will show that with a suitable choice of δ > 0 (as given in Theorem 9) the lentils cover the deep interior of the manifold but are small.This requires estimates in two directions for the lentils so that they are large enough to cover the relevant subset of M without gaps but are small enough so that they provide a useful density estimate.If every point is contained in a lentil and that lentil contains a source point, then the density of sources is bounded by the diameters of the lentils.The more complicated side of lentils is to ensure that they cover enough of the manifold, and this relies on a number of estimates based on bounded geometry.These density estimates constitute a proof of Theorem 9.
Density estimates turn out to be simpler in the reverse direction, showing that the methods used in the proof of Theorem 9 are fine enough to give a decent estimate of density.These reverse estimates prove Theorem 10.
In Section 10 we turn our attention to matters of convergence.The proof of Theorem 11 is mostly based on Theorem 9 and Theorem 10; as T increases, the source points become denser and denser both in reality and in terms of data-driven estimates, and thus the quality of the discrete approximation improves.
Theorem 1 follows a similar path and can be regarded essentially as a corollary of Theorem 11, and thus it has a relatively simple proof now that all the tools are available.The metric convergence results show that the two manifolds of Theorem 1 are isometric as metric spaces, and by the Myers-Steenrod theorem the isometry is in fact smooth.Proposition 12 gets its proof, using similar ideas and basic properties of Poisson point processes, at the very end.

Separation of graphs on a bundle
Let us denote the scalar second fundamental form by h 2 and the boundary metric (the first fundamental form) by h 1 .If p ∈ int(M ), ρ : M → R is the distance function ρ(x) = d(p, x) at p then the corresponding boundary distance function r : M → R, is just the restriction of ρ on the boundary ∂M .By the Gauss formula for hypersurfaces (see e.g.[39, Theorem 8.13 a]) and the definition of Hessians we have that We will make use of this throughout the paper.
Our method of proving Proposition 3 is by means of lifting the data from the boundary ∂M × R of the spacetime to a suitable bundle where the graphs do not intersect.We start by showing that when s = s , then the graphs of a s and a s can only have up to first order tangency when they intersect.Therefore 2-jets will separate them.Lemma 13.Let M be a simple Riemannian manifold and s, s ∈ int(M ) × R be two distinct points in the spacetime.Then if a s (x) = a s (x) at some point x ∈ ∂M , then the gradient or the Hessian of the two functions will differ at x.
Proof.Take any two s, s ∈ S so that a s and a s , their differentials, and their Hessians agree at a point x ∈ ∂M .We aim to show that s = s .
First, observe that the differential da s (x) is the covector corresponding to the tangential component of the velocity of the unit speed geodesic from π(s) to x. Therefore da s (x) = da s (x) implies that p := π(s) and p := π(s ) are on the same geodesic starting from x.
Let γ be this unit speed geodesic starting at x.We have γ(t 0 ) = p and γ(t 0 ) = p for some t 0 , t 0 > 0.
Consider the Hessian H := H ρp of the interior distance function ρ p : M → R defined by ρ s (x) = d(x, p) and its counterpart H corresponding to p .If a normal Jacobi field J along γ satisfies HJ(t) = −D t J(t) for some t < t 0 , then J(t 0 ) = 0, (see for instance [39,Proposition 11.2]).Similarly, H J(t) = −D t J(t) for any t < t 0 implies J(t 0 ) = 0.
Both Hessians satisfy H ( ) γ(t) = − γ(t) for t < min(t 0 , t 0 ).As the geodesic meets ∂M transversally at x and the boundary Hessians of a s and a s agree there, we have in fact the equality H(x) = H (x) for the Hessian operators on T x M .This follows from (37).Now take any nonzero w ∈ T x M normal to γ(0) and let J be the Jacobi field along γ with J(0) = w and D t J(0) = Hw = H w. As observed above, this implies J(t 0 ) = 0 and J(t 0 ) = 0. Due to the lack of conjugate points p = γ(t 0 ) = γ(t 0 ) = p .
As the two source points are equal distance t 0 = t 0 from x and the two arrivals are at the same time, we also have τ (s) = τ (s ).This concludes the proof of s = s .
Proof of Proposition 3. Let us abbreviate Q := Q(S) ⊂ ∂M × R. As S ⊂ M × R is discrete and M is compact, there are only finitely many points of S in any interval M × [a, b].By simplicity there is no geodesic longer than the diameter of the manifold.Thus Q written as is a locally finite union on ∂M × R.
Let σ(Q) be the set of "smooth points of Q", where Q is given locally as a single graph.By definition σ(Q) ⊂ Q is open.It is also dense, for otherwise there would be two graphs that coincide in an open set, contradicting Lemma 13.
Let E be the bundle over ∂M × R with fibers This E is the product of the bundles T x M and T x M ⊗ T x M pulled to ∂M × R over the projection ∂M × R → ∂M .For any smooth function h : U → R defined in an open set U ⊂ ∂M we the lift its graph as LG(h) ⊂ E so that This defines a smooth submanifold of E.
In this way we lift all of σ(Q There are no second order intersections by lemma 13, so the smooth submanifolds L(G(a s )) ⊂ E are pairwise disjoint.Therefore Q determines the set Projecting from E down to ∂M × R gives the graphs G(a s ) and thus also the functions a s , proving the first part of the claim.
For the second claim we simply localize (41) to be over Ω ⊂ ∂M × R, and we get a disjoint union of the sets L[G(a s ) ∩ Ω].A single graph G(a s ) may be cut into several pieces by this procedure and we do not necessarily know which pieces correspond to the same source.Therefore what we obtain is the collection of connected components of graphs restricted to Ω.
It was proven in [30] that the boundary distance function r x : ∂M → R of x ∈ int(M ) determines the point x uniquely.The next lemma improves this slightly: the boundary distance function modulo constants is enough for this uniqueness.Lemma 14.Let M be a simple Riemannian manifold and x, y ∈ int(M ) any two distinct points.Denote by r x and r y the boundary distance functions ∂M → R from these points.The difference function r x − r y cannot be a constant function.
Proof.Suppose that the difference function r x − r y takes the constant value c ∈ R despite x = y.If c = 0, the boundary distance functions coincide and thus by [30] x = y, yields a contradiction.Therefore we may suppose that c = 0.
Let γ be the maximal geodesic through the points x and y with endpoints x , y ∈ ∂M ordered so that x is closer to x and y to y.By assumption we have Due to simplicity and the order of the points, we have d(x , y) = d(x , x) + d(x, y).Therefore 0 = d(x, y) + c.The same calculation with the roles of x and y reversed shows that 0 = d(x, y) − c.These two together imply that d(x, y) = c = 0, which is a contradiction.

Exact observables
Proof of Theorem 2. Due to Proposition 3 the data determines the arrival time functions a s (x) = d(π(s), x) + τ (s) at all points x ∈ ∂M for every s ∈ S. It also follows from the same proposition that distinct sources have distinct arrival time functions.Therefore now that the two manifolds M 1 and M 2 and their source sets S 1 and S 2 have the same data modulo identification by φ, we have a bijection between the source sets given by identifying the corresponding arrival time functions which coincide.This proves claim 1.
We will prove all subsequent claims by describing how the arrival time functions determine the quantities in question.As these functions coincide on the two manifolds, the reconstructed quantities coincide.We therefore drop the subscripts and work with a single manifold.We can think of the source points being indexed by their arrival time graphs.To simplify notation, we denote the source points as p s = π(s) and source times as t s = τ (s).
The arrival time functions now have the form a s (x) = d(p s , x) + t s .Two source points p s and p r coincide if and only if a s − a r is a constant function -this is a straightforward corollary of Lemma 13 and Lemma 14.We pass to a subset of S so that each source point p s is only present once; the claims extend easily to the duplicated source points.
For any s, r ∈ S we define the functions f rs : ∂M × ∂M → R by These functions are determined by the data.Taking s = r proves claim 5.
From now on we suppose that s = r.The differentials da s and da r agree at x ∈ ∂M if and only if the geodesics joining x to p s and p r start in the same direction at x. Therefore the differentials agree at exactly two points, the endpoints of the unique maximal geodesic through the points p s and p r .Let these boundary points be x and y.
Depending on how the four points (p s , p r , x, y) are ordered on the geodesic, we have Thus claim 2 is given by This determination of distances between source points will play a central role in the proofs of the other theorems.By switching the points x and y if needed (depending on the sign in ( 45)), we can assume to have d(p r , y) = d(p r , p s ) + d(p s , y).Then for any z ∈ ∂M , we have As everything on the left-hand side is determined by the data, so is then the righthand side, and we have obtained the function f rs : ∂M → R given by for all s, r ∈ S.This proves claim 4. Finally, we note that for any z ∈ ∂M we have and claim 3 follows.
As mentioned, the proof of Theorem 1 will be postponed until we have approximation tools available.

The labeled Gromov-Hausdorff distance
Our proof of the basic properties of the labeled Gromov-Hausdorff distance is similar to well known proofs of the usual Gromov-Hausdorff distance (see for instance [20,Chapter 3] or [6, Section 7.3]) but with the labels interwoven into it.For the sake of completeness we record our proof in its entirety, although the aspects unrelated to labels are indeed well known.
Proof of Proposition 5.The symmetry of the labelled Gromov-Hausdorff distance is evident.The triangle inequality can be verified as in [6,Proposition 7.3.16].Also, if there is h as described, then choosing Z = Y , f = h, and g = id shows that the distance is zero.The only nontrivial claim is that if the distance is zero, then such an h exists.Let us prove that.
Take any compact metric space Z and isometric embeddings f : X → Z and g : Y → Z.We can define a semimetric (all properties of a metric but d(x, y) = 0 need not imply x = y) d on the disjoint union W = X Y by letting The embeddings are isometric, so d Z (f (a), f (b)) = d X (a, b) and similarly for Y and g.
Take any k ∈ N. By assumption we have a metric space Z k and thus a semimetric k for all ∈ L. These semimetrics agree with the metrics on X and Y when restricted to either set.
By compactness there is a finite set A k ⊂ X so that the semiballs cover X.Similarly, there is a finite A k ⊂ X so that the semiballs cover Y .Similarly, there are finite sets B k and B k in Y so that the semiballs of radius 1 k cover Y and X.We define We can make all of these choices so that A k ⊂ A k+1 and B k ⊂ B k+1 .
For any k ∈ N on the finite set j=k has a converging subsequence and the pointwise limit of semimetrics is a semimetric on A k B k .Constructing a diagonal sequence gives us a subsequence of (d j ) which converges pointwise on W = k∈N (A k B k ) to a semimetric δ.Because each d k agrees with the original metrics on X and Y , the only thing to inspect are the "cross-distances" between X and Y .
We can extend δ to all of W as follows.When x and y are both in X or both in Y , we use the metrics on these spaces.When x ∈ X and y ∈ Y , we pick for each A simple argument shows that this limit is independent of the choice of the approximating sequences.It is also straightforward to check that the extended δ is indeed a semimetric on W .We want to show that for every x ∈ X there is a unique y ∈ Y so that δ(x, y) = 0.The triangle inequality and δ agreeing with the metric on Y shows that the point is unique.For existence, there is a point in x k ∈ A k that is 1 k -close to x, and there is y k ∈ B k that is 1 k -close to x k .We can then set y = lim k→∞ y k .Similarly, each y ∈ Y has a unique x ∈ X so that δ(x, y) = 0.
This gives rise to a bijection h : X → Y that satisfies δ(x, h(x)) = 0.This can be checked to be an isometry.
We record two additional propositions.The proofs are immediate and we omit them.
Proposition 15.Let X and Y be compact metric spaces and L any set.Let α : L → X and β : L → Y be any two functions.If L ⊂ L, then As d ∅ GH (X, ∅; Y, ∅) is the usual Gromov-Hausdorff distance between X and Y , the labeled kind of convergence implies the usual kind of convergence.
Proposition 16.Let X be a compact metric space and Y ⊂ X a subset.Let α : L → X and β : L → Y be any two functions on a set L. (56)

Bounded geometry
We begin our study of bounded geometry by proving Proposition 7 with the help of two lemmas, and then we move on to finding further estimates based on these basic bounds.Lemma 17.Let M be a simple manifold, let U x ⊂ T x M be the maximal domain of definition of the exponential map, and let U ⊂ T M be the subset with fibers U x .The map exp × π : U → M × M that maps is a diffeomorphism and so has a smooth inverse θ : Proof.As each exp x : U x → M is a diffeomorphism on a simple manifold, the map exp × π is clearly smooth and bijective.What remains to check is the invertibility of the differential at every point (x, v) ∈ U .The differential has a convenient block structure due to π(x, v) = x being independent of v, and so d(exp × π) is invertible at (x, v) if and only if d exp x is invertible at v. The invertibility of d exp x is true by assumption.
Lemma 18.Let W 1 , . . ., W n be vector fields on a simple manifold M constituting a global orthonormal frame.Let Γ i jk (x, y) be the Christoffel symbol at x ∈ M with the normal coordinates centered at y ∈ M with W 1 (y), . . ., W n (y) as the coordinate basis for T y M .Then sup x,y∈M for all indices i, j, k.
There is a global orthonormal frame on every simple manifold given by vector fields W i as stated.Such a frame can be produced by the Gram-Schmidt method from a general frame coming from the trivializability of the tangent bundle.
Proof of Lemma 18.The vectors W 1 (y), . . ., W n (y) constitute an orthonormal basis for T y M .With the help of Lemma 17 we see that the basis vectors for T x M in the normal coordinates of y are given by Here and later in this proof we indicate the variable of differentiation by superscript in [square brackets] when needed for clarity.The invariant coordinate map for the normal coordinates of y is θ y : M → T y M given by θ y (x) = θ(x, y).Using the basis given by the frame, we get a proper coordinate map θ y : M → R n given by θ y i = W i (y) • θ y , which means Now we have a concrete description of the normal coordinates in terms of the frame.Let us denote . This is the differential of a coordinate map and thus a basis covector for T * x M induced by the normal coordinates of y.The corresponding basis vectors on T x M are w i (x, y).In terms of the coordinate map z = θ y : M → R n the Christoffel symbols are given by All the items in this formula depend smoothly on x and y, up to the boundary.The claim then follows from compactness of M × M .
Proof of Proposition 7. Parts 1, 2, 5, and 6 of Definition 6 follow easily from compactness.We will prove the other parts more carefully.
Part 3: Let us denote y i = exp x (η i ) ∈ M .In terms of Lemma 18 we have θ(y i , x) = η i .Let us equip T M with the Sasaki metric and M × M with the product metric.An important property of the Sasaki metric on T M is that it respects the natural distance on the fibers: the Sasaki distance between v, w ∈ T x M is simply |v − w| because the fibers are totally geodesic [48].Because θ is a diffeomorphism between compact manifolds with boundary, it is Lipschitz-continuous.The Lipschitz constant is denoted by C exp , that is the maximum value of the norm of the Jacobian of θ.
Part 4: Let us decompose an arbitrary Jacobi field J starting with J(0) = 0 into parallel and normal components relative to γ: As 2 and both components are Jacobi fields, it suffices to prove the estimate separately for parallel and normal Jacobi fields.
If J is tangential to γ then J(t) = tc γ(t) for some c ∈ R and |J(t)| 2 = c 2 t 2 .Thus for the tangential Jacobi fields the claim holds for any C JF ≥ 2. Now we assume that J is normal to γ.Let ρ(x) = d(x, γ(0)).Then where λ max (t) is the eigenvalue of the Hessian H ρ of ρ at γ(t) with the largest absolute value.Let us choose normal coordinates at γ(0) as in Lemma 18.In these coordinates we have ρ(x) = |x| e for all x ∈ M , where the subscript e stands for the Euclidean norm.Geodesics with initial unit speed v ∈ S γ(0) M are given by γ(t) = tv.
Therefore the Hessian of the distance function ρ is given in these coordinates by where The numbers A ij are all bounded uniformly for all γ(0) ∈ M and x ∈ M \{γ(0)} due to Lemma 18.Thus the Hessians H ρ , as quadratic forms, are uniformly bounded by some constant 1 2 C JF , and |λ max (t)| ≤ 1 2 C JF .The claim follows.Part 7: This is similar to part 4 and we continue with the same proof setup.The point where everything is evaluated is x = tv, so by equations ( 64) and (65) we have for all w ∈ γ(t The bound on the Christoffel symbols implies that Γ k ij (tv) ≤ C for some absolute constant C ≥ 0 independent of t, v, and γ(0).In the Euclidean norm v e = 1, so for all w ∈ γ(t) ⊥ .By the diffeomorphic nature of exp × π of the proof of Lemma 18 the Euclidean norm of w is bi-Lipschitz equivalent with the Riemannian norm of w, with bi-Lipschitz constants independent of x, and γ(0).This ensures that the desired two constants exist.
We will next establish a number of estimates based on the results of Proposition 7. Most of them we will rewrite in Section 9.1 in the language of lentils.Before getting started, we recall three standard comparison results, stated in a form suitable for our use.
We will only use Lemma 19 to compare to constant sectional curvature, and we state the next two lemmas only in this setting.We denote by M κ the complete Riemannian manifold with constant sectional curvature κ ∈ R, either the hyperbolic space, the Euclidean space, or the sphere, depending on the sign of κ.
Consider a simple Riemannian manifold M with sectional curvature bounded from above by C sec+ .For any three points x, y, z ∈ M we want to pick three points x, ỹ, z ∈ M Csec+ so that We denote distances on both manifolds by d, as there should be no confusion as to which space each point belongs.As long as M satisfies the estimate (3), such points on M Csec+ exist.
The following two lemmas compare triangles on M with corresponding triangles on M κ .They both follow from Rauch's comparison theorem, although not exactly the same version as Lemma 19 above.We recall from Definition 6 that the sectional curvature is bounded from above by the constant C sec+ > 0.
Lemma 20 ([16, Chapter 10, Proposition 2.5]).Let M be a simple manifold whose constants (of Definition 6 and Proposition 7) satisfy the constraint (3).Take any three points x, y, z ∈ M .Then there are points x, ỹ, z ∈ M Csec+ so that d(x, y) = d(x, ỹ), d(x, z) = d(x, z), and the angle is the same at x and x.The remaining distances satisfy d(y, z) ≤ d(ỹ, z).
Finally we give a comparison result for the angles of geodesic triangles.
Lemma 21 ([29, Theorem 2.7.6]).Let M be a simple manifold whose constants (of Definition 6 and Proposition 7) satisfy the constraint (3).Take any three points x, y, z ∈ M .Then there are points x, ỹ, z ∈ M Csec+ so that d(x, y) = d(x, ỹ), d(x, z) = d(x, z), and d(y, z) = d(ỹ, z), and the angle at x is at most that at x.We will then proceed to estimates specific to our setting.Lemma 22.Let M be a simple Riemannian manifold of bounded geometry.Let J be a Jacobi field along any constant (not necessarily unit) speed geodesic with for all t > 0, where the constant is given by (152).
Proof.Let us first prove the claim for unit speed geodesics.
If M has sectional curvature ≥ k, then by Lemma 19 every Jacobi field J with J(0) = 0 and |D t J(0)| = 1 satisfies |J(t)| ≤ (−k) −1/2 sinh(t √ −k) before a conjugate point.As M is simple, there are no conjugate points.The function x → x −1 sinh(x) is increasing for x ∈ [0, ∞) and diam(M ) ≤ C diam .Thus we get the desired estimate with the constant sup as claimed.
Let Γ(t, s) be a family of geodesics so that Γ(t, 0) = γ(t) and J(t) = ∂ s Γ(t, s)| s=0 for a Jacobi field J along γ satisfying the assumptions of the claim.Consider then the new Jacobi field Ĵ(t ) = ∂ s Γ(λt , λ −1 s )| s =0 along the geodesic γ.The scaling of the two parameters was chosen so that Ĵ is a vector field along γ and D t Ĵ(0) = D t J(0), so that both covariant derivatives have unit norm at γ(0) = γ(0).As we have Ĵ(t ) = λ −1 J(λt ), the already proven claim for J along γ gives (73) This is the claimed estimate for Ĵ.
Then we may estimate the distance between the two required points by taking the s-curve with t constant.First, the Jacobi field growth estimate of Lemma 22 gives By the uniform invertibility of the exponential maps due to bounded geometry we have as claimed.
9. Density estimates in bounded geometry 9.1.Lentils.We remind the reader that a lentil is defined simply as We will always assume that r, s < d(x, y).
We say that the thickness δ x,y r,s of this lentil is the length of the geodesic segment γ x,y ∩L x,y r,s , where γ x,y is the unit speed geodesic from x to y.Under the assumptions r, s < d(x, y) this leads to δ x,y r,s = r + s − d(x, y).
The assumptions r, s < d(x, y) imply that δ x,y r,s < r, s.The midpoint of the lentil is The transversal radius of the lentil L x,y r,s with midpoint m = m x,y r,s is R x,y r,s = sup{R > 0; exp m (η) ∈ L whenever η ∈ T m M , η ⊥ γ and |η| < R}. (82) That is, (the closure of) the lentil contains a geodesic disc -a ball of dimension n − 1 -of this radius normal to γ x,y .
Lemma 24.Consider a simple manifold of bounded geometry with the constants satisfying the condition (3).If a lentil L x,y r,s has thickness δ > 0, then its diameter satisfies diam(L x,y r,s where the constant is given by (170).
The condition (3) ensures that sin(C diam C sec+ ) > 0 and thus C e is a well defined and positive constant.
Proof of Lemma 24.The proof is composed of two parts.In the first part we show that if δ is small, then the opening angle of the lentil as seen from x is small.In the second part we use this to show that the lentil itself is small.
Let us denote L := L x,y r,s .Let γ = γ x,y be the geodesic between x and y and let := d(x, y).Take any point z ∈ L. The three points x, y, z form a triangle with sides , d(x, z) =: a, and d(y, z) =: b.Let us call the angle at x by the name β ∈ [0, π].
We will estimate this angle by comparing to a constant curvature reference manifold.As earlier, let M Csec+ be the sphere of the same dimension as M with constant sectional curvature C sec+ > 0. In the corresponding geodesic triangle with the same side lengths , a, b the corresponding angle is β.By Lemma 21 we have β ≤ β.
By the spherical law of cosines we have where we denoted ˜ = C sec+ , ã = a C sec+ , and b = b C sec+ .We will also use scaled versions of the two radii and the thickness: r = r C sec+ , s = s C sec+ , and δ = δ C sec+ .
By the assumption that z ∈ L we have a < r and b < s, and the reverse triangle inequality gives a > − s.These same inequalities hold when scaled with C sec+ , whence We will simplify this further using sin(r − δ) ≥ sin(r) − δ.Finally, we make use of the fact that sin(t) ≥ T −1 sin(T )t for all t ∈ [0, T ] when T ∈ (0, π).All the distances on the manifold are bounded by the diameter bound C diam , so in our setting T = C diam C sec+ , which was indeed assumed to be below π.Combining these, we find with the constant given by (171).(In the limit C sec+ → 0 we have C f = 1, and one can indeed verify the estimate with this constant using Euclidean comparison geometry when the sectional curvature is non-positive.)Using cos(β) ≤ 1 − 1 5 β 2 , which is valid for all β ∈ [0, π], we get which leads to This explicit estimate for β in terms of δ and r concludes the first part of the proof.Now consider any two points z 1 , z 2 ∈ L. Estimating their distance will amount to estimating the diameter of the lentil.Let β i be the angle at x between γ x,y and γ x,zi .
If sectional curvature and length of geodesics is bounded uniformly, then By (90) we have β i ≤ 5C f δ/r, and clearly d(ẑ 1 , ẑ2 ) < δ when ẑ2 is defined analogously to ẑ1 above.The triangle inequality gives Thus the diameter of the lentil can be estimated by as claimed.
Lemma 25.Consider a simple manifold M of bounded geometry with the constants satisfying the condition (3).Take any x, y ∈ M and consider the lentil L x,y r,s with r, s ∈ (0, d(x, y)).If m is the midpoint of L x,y r,s , then B(m, 1 2 δ) ⊂ L x,y r,s and thus R x,y r,s ≥ 1 2 δ.Proof.Recall that L x,y r,s = B(x, r) ∩ B(y, s).As the distances from m to x and y are r − 1 2 δ and s − 1 2 δ, respectively, it follows from the triangle inequality that B(m, 1 2 δ) is contained in both balls defining the lentil.The estimate for the transversal radius follows immediately.
Lemma 25 is sufficient for our needs, but we point out that the exponent of δ can be improved.We provide this observation as Proposition 30 in Appendix A.
The endpoints x and y of the relevant geodesics will be source points in Γ ⊂ B(∂M, ε 1 ).As the source set is countable in practice and can well be finite, the geodesics will not cover the deep interior of the manifold.But if the geodesics are fattened suitably, with the suitable level given by lemma 27, they do cover the desired set without gaps.Such a thickening is most naturally realized in our setting as a union of lentils with x, y, and δ > 0 fixed.Lemma 25 and Proposition 30 gives a lower bound on the fatness of this "sausage".9.2.Boundary density estimates.We will now describe how the quantities defined in equations ( 12), (13), and ( 14) estimate the distance of source points to the boundary.Consider a source s ∈ S located at the point p = π(s).The boundary distance function r p : ∂M → R is given by r p (x) = d(p, x) and the arrival time function a s is this function shifted by a constant.Proximity of a source point p ∈ int(M ) to the boundary causes the Hessians of both these functions to blow up.However, some of the curvature of the graph can be caused by the curvature of the boundary itself rather than the proximity of the source point.Bounded geometry gives a concrete bound on this effect and allows us to ensure with concrete estimates that the source point really is close to the boundary.
Recall that c(p) ⊂ ∂M is the set of critical points of the boundary distance function r p .Lemma 26.Let M be a simple Riemannian manifold with bounded geometry and ∅ = P ⊂ int(M ).Let ε 1 > 0. The quantities defined in equations ( 12), (13), and (14) satisfy the following estimates: (1) Take any p ∈ P .If x ∈ c(p), then d(p, x) ≤ E(p, x).
(3) For any x ∈ ∂M there is p ∈ P so that d(p, x) ≤ E + ε 1 .
Proof.Part 1: Take any p ∈ int(M ).Let ρ : M → R be the distance function ρ(x) = d(p, x) and r = ρ| ∂M .The Hessians H ρ and H r of these functions are symmetric quadratic forms on T M and T ∂M , respectively.Let γ be any unit speed geodesic with γ(0) = p and J a Jacobi field along γ so that J(0) = 0. Suppose γ meets the boundary orthogonally, which is equivalent with the exit point being a critical point of r.Then the Hessian of the distance function has the property (cf.proof of Lemma 13) for all t > 0 and thus by part 4 of Definition 6 we also have As there are no conjugate points, this amounts to in the sense of quadratic forms on γ⊥ ⊂ T γ(t) M .If λ(x) denotes the smallest eigenvalue of H r (x) at a boundary point x ∈ ∂M , then (37) yields for all p ∈ P and y ∈ c(p).Thus by the definition of E(x) in ( 13) there is a point p ∈ P with a distance to x less than E(x) + ε 1 to x.Part 3: Follows immediately from the previous one.
With the aid of Lemma 26 we can now prove Proposition 8.
Proof of Proposition 8. Proposition 3 implies that the data determines the arrival time functions a s for all s ∈ S, although we do not have a description of the index set S yet.As these functions differ from the boundary distance functions r s : ∂M → R, r s (x) = d(x, π(s)), only by a constant, we thus know the differential and the Hessian of each r s on all of ∂M .Therefore the data determines the critical points of these functions and the function E(p, y), for p = π(s) and y ∈ c(p), of (12).From this one can easily compute E(x) and E from ( 13) and ( 14).The set Γ is defined in terms of these quantities, so it is uniquely determined by the data as well.
Theorem 2 indicates that the data determines the pointwise spatial distances of the source points, and the last part of the claim follows.9.3.Interior density estimates.For density of sources deep in the manifold M , we use two estimates.The first one concerns the density of the set of geodesics connecting near-boundary source points in Γ.The second one ensures that the lentils cover enough of the manifold.Lemma 27.Let M be a simple manifold of bounded geometry and let ε 2 > 0. Let Γ ⊂ M be such that ∂M ⊂ B(Γ, ε 2 ).Then for every z ∈ int(M ) there are points x, y ∈ Γ so that Proof.Fix any x ∈ Γ and take any z ∈ int(M ).Let γ = γ x,z be a constant speed geodesic with γ(0) = x and γ(t) = z for some t > 0. We extend this geodesic beyond z so that it meets ∂M at some time t > t.We scale the constant speed so that t = 1 and we denote ŷ := γ(1) ∈ ∂M .There is y ∈ Γ so that d(y, ŷ) < ε 2 .Let γ : [0, 1] → M be the constant speed geodesic for which γ(0) = x and γ(1) = y.
By Lemma 23 we have As t < 1, we have thus as claimed.
Proof.We showed in Lemma 27 above that for any z ∈ int(M ) there are x, y ∈ Γ so that the distance from z to the trace of the geodesic γ x,y is at most C d ε 2 .We will choose h > 0 later (see (108) for an explicit expression) and require that d(x, z) ≥ h and d(y, z) ≥ h.These additional requirements may not hold for the pair of points x, y ∈ Γ without further assumptions for z.As both x and y are ε 1 -close to the boundary, the additional requirements d(x, z) ≥ h and d(y, z) ≥ h are certainly satisfied if d(z, ∂M ) ≥ ε 1 + h.This is why we placed a boundary distance assumption on z in the statement of this lemma.
Let z ∈ γ x,y be the nearest point on this geodesic to the point z.After choosing h large enough we obtain by the reverse triangle inequality This implies that z is an interior point of γ x,y , and the geodesics γ x,y and γ z,z meet orthogonally at z .To make sure that z is contained in a lentil L x,y r,s , with thickness δ and midpoint z , we require that the transversal radius of that lentil is more than d(z, z ).By Lemma 25 the aforementioned condition of the transversal radius is satisfied when The condition for z to be a midpoint of a lentil of thickness δ is that d(x, z ) > 1 2 δ and d(y, z ) > 1 2 δ, which by equation ( 103) are satisfied when Combining the two conditions (104) and (105) gives and therefore The conditions (104) and (105) then become so we may choose This is why we chose C 1 in (154) as we did.
Finally we note that by the choice of z ∈ γ x,y and due to (103), (108), and (110) we have d(x, y) > δ.If we set r = d(x, z )+ 1 2 δ and s = δ +d(x, y)−r = d(z , y)+ 1 2 δ, then the lentil L x,y r,s is of width δ and has z as a midpoint.Moreover, by the previous argument this lentil contains the point z, if d(z, ∂M ) ≥ ε 1 + C g ε 2 , as required in the claim of this lemma.
The choices of h and δ in the proof above are optimal up to a constant and the explicit choice simplifies our estimates.9.4.Global density estimates.With all the estimates we have collected, we are ready to prove both parts of Theorem 9. Claim 1 states that if all lentils of thickness δ meet a source point (which is verifiable from data), then the source set P 1 ⊂ M 1 has a concrete density estimate.The thickness δ is chosen carefully depending on ε 2 .Claim 2 states that under the assumptions of the previous claim the measurements define a discrete metric space which is a quantitatively good approximation for the true space.
Proof of Theorem 9.For claim 1, take any z ∈ M .We wish to show that there is p ∈ P so that d(z, p) < ε with ε chosen as in (22).We split the proof in two cases: near the boundary and deep within the manifold, and use different tools in either of these cases to verify the validity of the density claim.
The boundary case: If d(z, ∂M ) < ε 1 + C g ε 2 , then there is x ∈ ∂M so that d(z, x) < ε 1 +C g ε 2 .By claim 3 of Lemma 26, there is p ∈ P so that d(x, p) < E +ε 1 .Recall that ε 2 = E + ε 1 > ε 1 .By the triangle inequality we have This is a sufficient estimate near the boundary.The deep interior case: By definition of the set Γ in (17) and Lemma 26 we have Γ ⊂ B(∂M, ε 1 ) and ∂M ⊂ B(Γ, ε 2 ).If d(z, ∂M ) ≥ ε 1 + C g ε 2 , then Lemma 28 shows that there is a lentil L = L x,y r,s of the correct thickness δ = C 1 ε 2 < r, s so that z ∈ L. By the crucial assumption of the theorem, there is p ∈ P ∩ L. Lemma 24 then gives us This is a sufficient estimate for deep interior points.Now it remains to combine the estimates for points in the deep interior and near the boundary.Combining (111) and (112) gives for all z ∈ M .This proves that P ⊂ M is ε-dense with the choice of ( 22) and concludes the proof of claim 1.
The labeled Gromov-Hausdorff distance is now straightforward to estimate for claim 2. By the definition of E(x) in ( 13) we can choose a function α : ∂M → P that satisfies (23).This function can be actually identified from the data Q(S) as we know the first fundamental form of ∂M and the value of E(p, y) for all p ∈ π(S) and y ∈ c(p), given by ( 12), is obtained by the computing the gradient and the Hessian of the corresponding arrival time function.Then by Lemma 26 it follows that for each x ∈ ∂M there exists y ∈ c(α(x)) that satisfies By the choice of ε in (22) this proves estimate (24) and completes the proof.9.5.Reverse density estimates.We do not need any more preparations before the next proof.
Proof of Theorem 10.Consider a unit speed geodesic γ starting at a source point p = γ(0) ∈ P .Take any t > 0 and any w ∈ T γ(t) M orthogonal to γ(t To utilize the assumption that P ⊂ M is ε-dense with ε chosen by (25), suppose that d(∂M, p) < ε and let z p ∈ ∂M be the nearest boundary point to it.Estimate (118) then implies that when ε < C H1 /(C H2 + C SFF ).Due to our choice of ε, we have and so Consider then an arbitrary point x ∈ ∂M .There is a source point q ∈ P with d(x, q) < ε.Let z q ∈ ∂M be the closest boundary point to q.We have d(z q , q) < ε, and by bounded geometry (part 6 of Definition 6) By (121) we get As this bound it is independent of x ∈ ∂M , we have the same bound for E from ( 14): under the assumption (120).
Combining our choice ε 1 = C 5 ε with (124), we get All other estimates for ε 2 follow from this one.
To satisfy (26), we require that Thus (120) and (128) are satisfies due to (25), and by the previous remark the estimate (26) follows.
It remains to show that each lentil L = L x,y r,d(x,y)−r+δ with x, y ∈ P contain source points as claimed.The thickness of the lentil L is δ = C 1 ε 2 , so by Lemma 25 B(m, 1  2 δ) ⊂ L, where m is the midpoint of the lentil L. Therefore each lentil of thickness δ contains a ball of radius 1  2 δ.It follows from ε-density of P ⊂ M that that each such lentil contains a source point because This is why we chose C 5 so that 10. Convergence of discrete approximations for all x ∈ ∂M .Therefore the set Q(S, T ) of ( 27) contains the full graph of a s for all the sources s ∈ S with 0 < τ (s) < T − C diam .Let us denote the set of source points with their complete graphs contained in ∂M × (0, T ) by P T ⊂ π(S).
If T is too small, we may have P T = ∅.In this case we choose the metric approximation M T to be a set of one point and α T : ∂M → M T the constant map.When there are sources, we set M T = P T .
Take any ε > 0 and let ε be given by (25).As T >0 P T is dense in M by assumption, by compactness there is T (ε) > 0 so that P T is ε-dense for all T ≥ T (ε).
With the choices ε 1 = C 5 ε and ε 2 = ε 1 + E(T ), where E(T ) is defined as in (14), and δ = C 1 ε 2 we have by Theorem 10.By the last claim of Theorem 10 the assumption of Theorem 9 is satisfied.Thus the set P T is also ε-dense as required in claim (2) of Theorem 9. We choose a map α T : ∂M → P T so that ( 23) is satisfied when α = α T and P = P T .Finally the estimate ( 24) with (131) yields This concludes the the proof.
10.2.Dense sources.Now we are finally able to prove Theorem 1.We will make use of a version of the Myers-Steenrod Theorem [45,46] which states that a metric isometry between smooth Riemannian manifolds is necessarily smooth.In the case of simple manifolds it is straightforward to prove and the boundary causes no technical trouble.We record the proof here for the sake of completeness.
Lemma 29 (Myers-Steenrod Theorem).If Φ : M 1 → M 2 is a metric isometry between simple Riemannian manifolds, it is smooth up to the boundary and an isometry also in the sense that Φ * g 2 = g 1 .
Proof.It suffices to show that Φ is smooth.Then preserving distances quickly implies Φ * g 2 = g 1 via differentiating distance functions.Take any x 1 ∈ int(M 1 ) and let x 2 = Φ(x 1 ).If γ is any constant speed geodesic with γ(0) = x 1 , then Φ • γ is a constant speed geodesic starting at x 2 with the same speed.Let Φ * : T x1 M 1 → T x2 M 2 be the map that maps the velocity vectors of these geodesics to each other.That is, if ∂ t γ(t)| t=0 = v, then ∂ The map Φ * preserves norms and inner products, and therefore it preserves the (squared) distance between any pair of vectors in T x1 M 1 .Thus Φ * is an isometry between finite-dimensional inner product spaces, and as a Euclidean isometry fixing the origin it is linear and thus smooth.As the exponential maps are diffeomorphisms on their maximal domain of definition, the map Φ = exp x2 • Φ * • exp −1 x1 is smooth.
Proof of Theorem 1.Let ι 1 : ∂M 1 → M 1 be the usual inclusion map and let ι 2 = φ : ∂M 1 → M 2 (with extended codomain here for convenience) so that we have boundary inclusions ι i : ∂M 1 → M i with the same domain.
For any T > 0 consider data on the bounded set ∂M × (−T, T ).This time the time interval extends into both the past and the future.As in the proof of Theorem 11 above, we get a finite metric space M T and a map α T : ∂M 1 → M T so that d ∂M1 GH (M T , α T ; M i , ι i ) → 0 (135) as T → ∞ for both i = 1, 2. By Proposition 8 the approximating metric space M T constructed from data is the same for the two manifolds M 1 and M 2 .We then apply Proposition 5.By the triangle inequality as T → ∞ and so d ∂M1 GH (M 1 , ι 1 ; M 2 , ι 2 ) = 0. Now Proposition 5 gives an isometry Φ : M 1 → M 2 so that Φ • ι 1 = ι 2 .The condition Φ • ι 1 = ι 2 simply means that Φ| ∂M1 = φ.By Lemma 29 the isometry Φ is actually smooth.
It remains to check that the sources s = (p, t) ∈ S 1 ⊂ int(M 1 ) × R are mapped correctly.By Theorem 2 we know that there is a bijection, but we have to verify that it corresponds to the isometry Φ as claimed.In light of Proposition 3, the data can be seen as a collection of graphs.Let us denote A i = {a s ; s ∈ S i }.The two manifolds having equivalent data means that the map φ * : A 2 → A 1 that takes a → φ * a = a•φ is bijective.We will show that the maps b i : S i → A i with b i (s) = a s are bijective, and therefore the natural bijection between the sources S 1 → S 2 is ξ = b −1 2 • (φ * ) −1 • b 1 .We need this map to satisfy (4), which now amounts to a (p,t) = a (Φ(p),t) • φ (137) as functions ∂M 1 → R for all (p, t) ∈ S 1 .This follows straightforwardly from the definitions and Φ being an isometry.Let us then show that b i is bijective for both i.To this end, take any two distinct sources s, ŝ ∈ S i .The two arrival time functions are a s = r s + τ (s) and a ŝ = r ŝ + τ (ŝ).
If the spatial source points are different, π(s) = π(ŝ), then r s − r ŝ cannot be a constant function on ∂M due to Lemma 14. Therefore a s and a ŝ do not coincide.If the spatial points are the same, then r s − r ŝ is a constant function 0. As s = ŝ but π(s) = π(ŝ), we must have τ (s) = τ (ŝ) and thus a s − a ŝ is a non-zero constant function.
We have thus proven that s = ŝ implies b i (s) = b i (ŝ).This concludes the proof that b i is injective and thus the theorem is proven.10.3.Stochastic convergence.The proof of the stochastic result only requires checking that the point process almost surely produces the correct kind of source set.
Proof of Proposition 12. Countability of the set S ⊂ M × R given by the homogeneous Poisson point process is certain.We need to prove that the following properties hold almost surely: (1) π(S) ∩ ∂M = ∅.
(3) π(S ∩ (M × [0, ∞))) ⊂ M is dense.The probability that a measurable set A ⊂ M × R contains k ∈ N points of the source set S is For this and other basic properties of homogeneous Poisson point processes, see e.g.[56,22].The first property follows simply from µ(∂M × R) = 0.The second property follows from A ∩ S being almost certainly finite when µ(A) < ∞.
For the third property, let (x i ) ∞ i=1 be a dense sequence in M and define for any i, j ≥ 1 the set A i,j = B M (x i , j −1 ) × [0, ∞) ⊂ M × R.
The projection π(S ∩ (M × [0, ∞))) is dense in M if #(A i,j ∩ S) > 0 for all i and j.As this event for individual indices has probability 1, the event for density is a countable intersection of events of full probability.Therefore density is almost certain.
The conditions of Theorem 11 are thus met almost surely.

by lemma 22 .
Estimating the distance between the two points by the length of a curve joining them, we thus have d(z 1 , ẑ1 ) ≤ C b rβ 1 .