Inverse Spectral Problems for Collapsing Manifolds II: Quantitative Stability of Reconstruction for Orbifolds
Matti Lassas Jinpeng Lu Takao Yamaguchi
 















 






theorem
lemmaLemma[section]
proposition[lemma]Proposition
coro[lemma]Corollary
theorem[lemma]Theorem



definition
definition[lemma]Definition
question[lemma]Question
remark[lemma]Remark
problem[lemma]Problem

equationsection


ℂ
ℝ
ℚ
ℤ
ℕ





vol




    supp diamdistdiag det⟨⟨|⟩ε∂ℳℱf̂ĥδαγβ̱Hδ̣ yΛ_3Matti Lassas: Department of Mathematics and Statistics, University of Helsinki, FI-00014 Helsinki, FinlandEmail address: Jinpeng Lu: Department of Mathematics and Statistics, University of Helsinki, FI-00014 Helsinki, FinlandEmail address: Takao Yamaguchi: Institute of Mathematics, University of Tsukuba, Tsukuba 305-8571, JapanEmail address: 
                  Inverse Spectral Problems for Collapsing Manifolds II: Quantitative Stability of Reconstruction for Orbifolds
                    Matti Lassas,  Jinpeng Lu,  Takao Yamaguchi 
                 
                
                
                
                Dedicated to the memory of Yaroslav Kurylev 
                    
                =============================================================================================================
                
                
                We consider the inverse problem of determining the metric-measure structure of collapsing manifolds from local measurements of spectral data.
                In the part I of the paper, we proved the uniqueness of the inverse problem and a continuity result for the stability in the closure of Riemannian manifolds with bounded diameter and sectional curvature in the measured Gromov-Hausdorff topology.
                In this paper we show that when the collapse of dimension is 1-dimensional, it is possible to obtain quantitative stability of the inverse problem for Riemannian orbifolds.
                The proof is based on an improved version of the quantitative unique continuation for the wave operator on Riemannian manifolds by removing assumptions on the covariant derivatives of the curvature tensor.
                
                
                § INTRODUCTION
                
                
                We consider the class of connected, closed, smooth Riemannian manifolds (M,g) defined by
                
                    dim(M)=n,    |R(M)|≤Λ^2,    diam(M)≤ D,
                
                where R(M),  diam(M) are the sectional curvature and diameter of M.
                Denote by μ_M the normalized Riemannian measure,
                
                    dμ_M=1/Vol(M)dV_g,
                
                where Vol(M) is the Riemannian volume of (M,g), and dV_g is the Riemannian volume element of (M,g).
                We denote this class of manifolds (M,g) equipped with their normalized measure μ_M by M M (n,Λ,D). Let M M(n,Λ,D) be its closure with respect to the measured Gromov-Hausdorff topology (<cit.>).
                
                It is well-known that a sequence of n-dimensional manifolds in the class M M(n,Λ, D) can collapse to a lower dimensional space when the injectivity radius of the sequence of manifolds goes to zero.
                A metric-measure space (X,μ)∈ M M(n,Λ, D) has the stratification 
                
                     X = S_0(X) ⊃ S_1(X) ⊃⋯⊃ S_d(X),
                       
                with the following property: if S_j(X)∖ S_j+1(X) is non-empty, then it is a 
                (d-j)-dimensional Riemannian manifold, where d=dim(X)≤ n. The regular part of X is the set
                X^reg:=S_0(X) ∖ S_1(X) which is an open d-dimensional Riemannian manifold of Zygmund class C^2_∗, that is, the transition maps between coordinate charts are of class C^3_* and the metric tensor in these charts is of class C^2_* (<cit.>).
                Recall that the Zygmund space C^2_*(X^reg) has the relation C^1,1(X^reg)⊂ C^2_*(X^reg)⊂ C^1,α(X^reg) for any 0<α<1.
                The singular set of X is the complement of the regular part, denoted by X^sing:=X∖ X^reg or S. In particular, the singular set X^sing has dimension at most d-1, and μ(X^sing)=0.
                The measure μ is absolutely continuous with respect to the Riemannian measure dV_g on X^reg, with the density function ρ_X∈ C^2_*(X^reg).
                In particular, if (X)=n-1 then X is a Riemannian orbifold <cit.>. 
                Recall that a metric space X is called a Riemannian orbifold if, roughly speaking, for each p∈ X there exists a neighborhood U of p such that U is isometric to the quotient of a Riemannian manifold by an action of a finite group of isometries (see e.g. <cit.> or Appendix <ref> for a formal definition).
                As a basic example, the product manifold 𝕊^2× [0,ϵ] with points ( x, 0) ∈𝕊^2 ×{0} and (e^2π i/m· x, ϵ) ∈𝕊^2 ×{ϵ} identified collapses to a Riemannian orbifold 𝕊^2/ℤ_m as ϵ→ 0, where e^2π i/m· x stands for the rotation on the unit sphere 𝕊^2⊂ℝ^3 by the angle 2π/m along the z-axis.
                The orbifold 𝕊^2/ℤ_m, having the shape of a rugby ball, has conic singular points only at the north and south poles. 
                Relevant figures and an example of collapsing manifolds in physics can be found in Section <ref>.
                
                
                
                As shown in <cit.>, the Dirichlet's quadratic form A[u,v]=⟨d|u,dv _L^2(X^reg,dμ)^2,
                u,v∈ C^0,1(X),  defines
                a self-adjoint operator Δ_X on (X, μ), which we call the weighted Laplacian, denoted by Δ_X.
                In a local coordinate (x^1,⋯,x^d) of (X^reg,g_X), it has the form
                
                    Δ_X u= 1/ρ_X |g_X|^1/2∑_j,k=1^d / x^j(ρ_X |g_X|^1/2 g^jk_X / x^k  u  ),    |g_X|=((g_X)_jk).
                
                Note that in the part I of the paper <cit.> we used the nonnegative definite Laplacian, while in this paper we use analyst's nonpostive definite Laplacian. 
                Denote by λ_j the j-th eigenvalue of the weighted Laplacian Δ_X and by ϕ_j the corresponding orthonormalized eigenfunction in L^2(X,μ),
                
                    -Δ_X ϕ_j=λ_j ϕ_j,    j=0,1,⋯.
                
                Fukaya proved in <cit.> that the j-th eigenvalue, for any j, is a continuous function on  M M(n,Λ,D) with respect to the measured Gromov-Hausdorff topology.
                
                
                We consider the following generalization of Gel'fand's inverse problem <cit.> in the class M M(n,Λ,D).
                
                
                Let (X,μ)∈ M M(n,Λ,D) and p∈ X. Suppose that one can measure the spectral data {λ_j,ϕ_j|_B(p,r)}_j=0^∞ in a ball B(p,r) of radius r>0. 
                
                Do these data uniquely determine the metric structure of X and the measure μ?
                
                
                In the special case of X being a priori a Riemannian manifold of dimension n, the problem is reduced to the classical case that was solved in <cit.>, with quantitative stability of log-log type proved in <cit.>. Stronger Hölder type of stability estimates can be obtained in <cit.> if additional geometric assumptions are assumed, e.g., if the metric is close to being simple.
                In general, the positive answer to this problem and an abstract continuity for the stability are proved in the part I of this paper <cit.>. In this part II, we investigate if a quantitative stability estimate for the problem can be obtained in the general case when X is not a manifold.  
                We also discuss applications of Gel'fand's inverse problem for collapsing manifolds in manifold learning and physics in Section <ref>.
                
                
                For the purpose of obtaining a stability estimate, we need to impose bounds for the covariant derivatives of the curvature tensor.
                
                Following similar notations in <cit.>, we say a Riemannian manifold M is Λ_k-regular if
                
                
                    ∇^i R(M) ≤Λ_k,  for  i=1,⋯,k.
                
                We denote the class of connected closed smooth Riemannian manifolds satisfying (<ref>) and (<ref>) by M M(n,Λ,Λ_k,D), and denote its closure by M M(n,Λ,Λ_k,D) with respect to the measured Gromov-Hausdorff topology.
                Our main result is to show that when dim(X)=n-1, i.e., X is an orbifold, it is possible to obtain a quantitative stability estimate for the reconstruction of the regular part from incomplete spectral data.
                
                
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose (X)=n-1 and Vol_n-1(X) ≥ v_0.
                Let σ∈ (0,1), and U⊂ X^reg be an open subset containing a ball B_X(p,r_0).
                Then there exists δ=δ(X,r_0,σ)>0, such that
                the finite interior spectral data {λ_j,ϕ_j|_U}_j=0^δ^-1 for δ<δ determine a finite metric space X such that
                
                    d_GH(X∖ S_σ,δ,X) < C_1(X,σ) (loglog|logδ| )^-C_2,
                
                where S_σ,δ is a subset (depending on σ,δ) of the C_3σ^1/4-neighborhood of the singular set of X, and X∖ S_σ,δ is equipped with the restriction of the metric of X.
                
                The constant C_1(X,σ) depends on σ,r_0,n,Λ,,D,v_0 as well as the space X, and C_1(X,σ)→∞ as σ→ 0. The constant C_2 depends on n, and C_3 depends on n,Λ,D.
                
                Moreover, suppose that the complete spectral data of (X,μ), (X',μ')∈ M M(n,Λ,D) on open sets U⊂ X, U'⊂ X' are equivalent by a homeomorphism Ψ_U: U→ U', namely λ_j=λ'_j and Ψ_U^*(ϕ'_j|_U')=ϕ_j|_U for all j∈, where λ'_j,ϕ'_j are the eigenvalues and orthonormalized eigenfunctions of the weighted Laplacian on (X',μ').
                Then there exists a measure-preserving isometry F: X→ X' such that F|_U=Ψ_U and μ=F^* (μ').
                
                
                
                
                
                Recall that the singular set of X is a union of submanifolds having dimension at most n-2.
                Note that the incomplete data case of the theorem above is of an asymptotic nature as the reconstruction is done to a fixed orbifold.
                Our method can also be extended to the case where the finite interior spectral data are given up to small noise.
                The case for complete spectral data in the theorem above, valid for any dimension of collapse, is due to the part I of this paper <cit.>.
                
                
                
                
                
                The basic idea for proving Theorem <ref> is as follows. Let M_i∈ M M(n,Λ,D) be a sequence of Riemannian manifolds converging to X. We consider the orthonormal frame bundles FM_i of M_i endowed with the natural Riemannian metric. Note that O(n) isometrically acts on FM_i and M_i=FM_i/O(n). Passing to a subsequence, we may assume that FM_i converge to a space Y, and X=Y/O(n). 
                The purpose of this construction is that unlike X, the limit space Y is always a Riemannian manifold, see <cit.>, and the Riemannian metric on Y is of class C^2 under the Λ_k-regular condition (<ref>), see Section <ref>.
                Thus, we can apply the quantitative unique continuation for the wave equation and a quantitative version of the boundary control method on Y to reconstruct the interior distance functions on U⊂ X in Section <ref> and <ref>. Next, using the interior distance functions, we approximately locate the singular set of X and reconstruct the smooth metric structure away from the singular set in Section <ref>. 
                However, due to the lack of higher order curvature bounds on the limit space Y, we need to improve the earlier results on the quantitative unique continuation <cit.> by removing the assumptions on the covariant derivatives of the curvature tensor.
                This result is stated as follows.
                
                
                Let (Y,g) be a connected, closed, smooth manifold with C^2-smooth Riemannian metric g, and let ρ_Y∈ C^1(Y) be a density function on Y.
                Consider the wave equation with the weighted Laplacian Δ_Y on Y,
                
                    (∂_t^2-Δ_Y) u= f,
                
                where Δ_Y is given by (<ref>) with ρ_X replaced by the density function ρ_Y.
                
                Let (Y,g) be a connected, closed, smooth manifold with C^2-smooth Riemannian metric g, satisfying 
                
                    dim(Y)=n,    |R(Y)|≤Λ^2,    diam(Y)≤ D,    Vol(Y)≥ c_0.
                
                Let ρ_Y∈ C^1(Y) be a density function on Y and consider the weighted Laplacian operator Δ_Y.
                Suppose u∈ H^1(Y×[-T,T]) is a solution of the non-homogeneous wave equation (<ref>) with f∈ L^2(Y× [-T,T]). Let V⊂ Y be an open subset with smooth boundary.
                If the norms satisfy
                
                    u_H^1(Y×[-T,T])≤ E,  u_H^1(V× [-T,T])≤ε_0,
                
                then for 0<h<h_0, we have
                
                    u_L^2(Ω(h))≤ C_4 exp(h^-C_5)E/(log(1+h^-1E/f_L^2(Y×[-T,T])+h^-2ε_0)) ^1/2  .
                
                The domain Ω(h) is defined by
                
                    Ω(h)={(x,t)∈ (Y ∖ V)× [-T,T]: T-|t|-d_Y(x,V) > √(h)}.
                
                The constants C_4,h_0 explicitly depend on n,T,Λ,D,c_0 and ρ_Y_C^0,1, and C_5 explicitly depends on n. The norms are with respect to the Riemannian volume element on (Y,g).
                
                
                
                
                The proof of Theorem <ref> is based on the smoothening of Riemannian metric <cit.> and is done in Section <ref>.
                
                
                Acknowledgements.
                The authors thank Atsushi Kasue for helpful discussions.
                 M.L. and J.L. were supported by PDE-Inverse project of the European Research Council of the European Union, 
                 the FAME-flagship and the grant 336786 of the Research Council of Finland. Views and opinions expressed are
                those of the authors only and do not necessarily reflect those of the European
                Union or the other funding organizations.  
                Neither the European Union nor the other funding organizations can be held responsible for them. 
                T.Y. was supported by JSPS KAKENHI Grant Number 21H00977.
                
                
                
                
                
                § PRELIMINARIES
                
                
                 §.§ Geometric structure of collapsing manifolds
                
                
                Given n∈_+, Λ,D>0, let M(n,Λ,D) denote the class of connected, closed, smooth Riemannian manifolds (M,g) satisfying
                
                    dim(M)=n,    |R(M)|≤Λ^2,    diam(M)≤ D,
                
                and let M(n,Λ,D) denote its closure in the Gromov-Hausdorff topology.
                Suppose that X ∈ M(n,Λ,D) is a compact metric space which is the Gromov-Hausdorff limit of closed Riemannian manifolds M_i ∈ M(n,Λ,D).
                
                
                Let FM_i be the orthonormal frame bundle of M_i. 
                Fixing a Riemannian metric on O(n), FM_i has a natural Riemannian structure with uniformly bounded sectional curvature by O'Neill's formula,
                
                    |R(FM)| ≤Λ_F^2=C(n,Λ).
                
                The dimension and diameter of FM satisfy
                
                    (FM)=n+(O(n)),  diam(FM)≤ D_F=C(n,D).
                
                
                
                The limit space X has the following geometric structure due to Fukaya.
                
                
                Let X∈ M(n,Λ,D). Then there exists a smooth manifold Y with C^1,α-Riemannian metric, for any 0<α<1, 
                on which O(n) acts as isometries in such a way that X is isometric to Y/O(n), and (Y)=(X)+(O(n)).
                
                
                
                
                The manifold Y is constructed as the Gromov-Hausdorff limit of the orthonormal frame bundles FM_i of M_i.
                Here we briefly explain why Y is a manifold. 
                Fix any point p∈ X and put  p_i:= ψ_i(p), where ψ_i:X→ M_i is
                an _i-Gromov-Hausdorff approximation with lim_i→∞_i = 0.
                Let B ⊂^n be the open ball around the origin O in ^n of radius π/Λ, and
                let exp_i: B→ M_i be the composition of the exponential map exp_p_i:
                 T_p_i(M_i)→ M_i and a linear isometric embedding 
                B→ B( O,π/Λ)⊂ T_p_i(M_i).
                Since |R(M_i)| ≤Λ^2, the exponential map on B(O,π/Λ) has maximal rank, and 
                we have the pull-back metric g̃_i := exp_i^*(g_i) on B⊂^n.
                Moreover, the injectivity radius of (B,g̃_i) is uniformly bounded from below, e.g. <cit.>.
                Therefore, we may assume that (B, g̃_i) converges to a
                C^1,α-metric (B, g̃_0) in the C^1-topology, 
                for any 0<α < 1 (<cit.>).
                
                Let G_i denote the set of all  isometric 
                embeddings γ:(B', g̃_i) → (B, g̃_i) such that
                exp_i ∘γ = exp_i on B', where B':= B( O, π/2Λ)⊂^n. 
                Then G_i forms a local pseudogroup, see e.g. <cit.>.
                Passing to a subsequence, we may assume that
                G_i converges to a local pseudogroup   G consisting of isometric 
                embeddings γ:(B', g̃_0) → (B, g̃_0), and that
                (B,g̃_i, G_i) converges to (B, g̃_0,G) in the 
                equivariant Gromov-Hausdorff topology. The quotient space (B',g̃_i)/G_i=B(p_i, π/2Λ)⊂ M_i 
                converges to (B',g̃_0)/G,
                which implies that 
                (B', g̃_0)/G=B(p, π/2Λ) ⊂ X.
                (See <cit.> for further details on basic properties of the equivariant Gromov-Hausdorff convergence.)
                
                The fact that Y is a Riemannian manifold can be seen as follows.
                The pseudogroup action of G_i on (B', g̃_i) induces an isometric pseudogroup
                action, denoted by Ĝ_i, on the frame bundle F(B', g̃_i) of
                (B', g̃_i) defined by differential. Therefore, 
                F(B', g̃_i)/Ĝ_i = FB(p_i, π/2Λ). Passing to a subsequence,
                we may assume that  (F(B', g̃_i), Ĝ_i) 
                converges to (F(B', g̃_0), Ĝ) in the equivariant
                GH-topology, where Ĝ denotes the
                isometric pseudogroup action on F(B',g̃_0) induced from that of 
                G on (B', g̃_0). 
                Since the action of G on (B', g̃_0) is isometric, the action of Ĝ on  F(B', g̃_0)
                is free. 
                
                
                
                
                
                
                Therefore, F(B',g̃_0)/Ĝ is a Riemannian manifold, and
                so is Y. The regularity of the metric follows from <cit.>.
                
                
                
                
                
                From the construction above, the diameter of Y is clearly uniformly bounded:
                
                    diam(Y)≤ C(n,D).
                
                Furthermore, under the additional -regular condition
                
                    ∇^i R(M) ≤,   for   i=1,2,3,
                
                
                
                
                the lifted metrics g̃_i are of C^4,α in harmonic coordinates in the construction above (<cit.>), which, passing to a subsequence, converge to a metric g̃_0 ∈ C^4,α in the C^4-topology. Since the induced pseudogroup action Ĝ on the orthonormal frame bundle F(B',g̃_0) is isometric and free, the Riemannian manifold Y has C^4,α-metric tensor and in particular, the sectional curvature of Y is defined.
                
                
                Let F M(n,Λ,,D) denote the set of orthonormal frame bundles FM of closed Riemannian manifolds M in the class M(n,Λ,D) satisfying the additional -regular condition (<ref>), and let F M(n,Λ,,D) denote its closure in the Gromov-Hausdorff topology. 
                
                The sectional curvature of any Y∈ F M(n,Λ,,D) is uniformly bounded:
                
                    |R(Y)| ≤ C(n, Λ,,D).
                
                
                
                The lower curvature bound is immediate: Y has curvature bounded from below in the sense of Alexandrov (e.g. <cit.>).
                The argument for the upper sectional curvature bound is due to <cit.>.
                Suppose that the sectional curvature is unbounded at y_k∈ Y_k, then <cit.> implies that the isotropy group G_y_0 at lim_k→∞ y_k=y_0∈ Y_0=lim_k→∞ Y_k has positive dimension. However, the (induced) pseudo-group action Ĝ on the orthonormal frame bundle is free, so the isotropy group is trivial, which is a contradiction.
                Note that <cit.> relies on the sectional curvature estimate in <cit.>, and the latter requires that Y has C^4-smooth Riemannian metric.
                
                
                
                
                The next lemma concerns the relation between the volumes of X and Y.
                
                
                Let X ∈ M(n,Λ,D) and Y ∈ F M(n,Λ,D) be defined as in Theorem <ref>. Then we have the following identity:
                
                
                
                
                    Vol_d+(O(n)) (Y) =  Vol(O(n))   Vol_d(X).
                
                
                
                
                
                
                Let M_i∈ M(n,Λ,D) be a sequence of manifolds converging to X, and, passing to a subsequence, the orthonormal frame bundles FM_i∈ F M(n,Λ,D) converges to Y.
                We have fibrations
                f_i:M_i→ X and f̃_i:FM_i→ Y such that 
                f_i∘π_i=π∘f̃_i, where π_i:FM_i→ M_i is the projection. On the other hand, for every x∈ X, as i→∞,
                (f_i∘π_i)^-1(x) and (π∘f̃_i)^-1(x)
                converge to O(n) and π^-1(x) respectively.
                Thus, π^-1(x) is isometric to O(n). π:Y→ X=Y/O(n) is the projection along 
                O(n)-orbits, and the distance of X is the orbit 
                distance, d_X(x,x')=d_Y(π^-1(x),π^-1(x')) for all
                x,x'∈ X.
                Thus, |dπ(v)|=1 if v is normal to the O(n)-orbit. 
                
                Then the claim follows from the co-area formula.
                
                
                
                
                In addition, we assume a lower bound for the volume of X,
                
                    Vol_d(X)≥ v_0,    d=(X).
                
                Thus, Lemma <ref> gives us a lower bound on the volume of Y.
                Then by Cheeger's injectivity radius estimate (<cit.>) and Lemma <ref>, we have
                
                    inj(Y)≥ i_0=i_0(n,Λ,,D,dim(X),v_0).
                
                
                
                 §.§ Density functions and the weighted Laplacian
                
                
                In <cit.>, Fukaya introduced a topology on the set of compact metric spaces equipped with Borel measures, called the measured Gromov-Hausdorff topology. A sequence of compact metric-measure spaces (X_i,μ_i) with μ_i(X_i)=1 is said to converge to (X,μ) in the measured Gromov-Hausdorff topology, if X_i converges to X in the Gromov-Hausdorff topology and, in addition, the measure (ψ_i)_∗(μ_i) converges to μ in the weak^* topology, where ψ_i:X_i→ X is an ε-Gromov-Hausdorff approximation.
                When X_i is a Riemannian manifold, the measure μ_i is taken to be the normalized Riemannian measure (<ref>).
                
                
                We denote by M M (n,Λ,D) the set of closed Riemannian manifolds M in the class M (n,Λ,D) equipped with their normalized Riemannian measure μ_M, and by M M(n,Λ,D) its closure with respect to the measured Gromov-Hausdorff topology.
                Similarly for the orthonormal frame bundles, we denote by F M M(n,Λ,D) the set of orthonormal frame bundles FM of closed Riemannian manifolds M∈ M (n,Λ,D) equipped with the normalized Riemannian measure on FM, and denote by F M M(n,Λ,D) its closure in the measured Gromov-Hausdorff topology.
                
                
                
                Let (X,μ) ∈ M M(n,Λ,D).
                By <cit.>, there exists a density function ρ_X∈ C^2_∗(X^reg)⊂ C^1,α(X^reg), for any 0<α<1, such that 
                
                    dμ = ρ_X dV_X/Vol(X),
                
                where dV_X is the Riemannian volume element on the regular part X^reg of X. 
                
                Note that ρ_X can be zero at non-orbifold point of X, see <cit.>.
                For (Y,μ̃) ∈ F M M(n,Λ,D), since Y is always a Riemannian manifold due to Theorem <ref>, there exists a strictly positive density function ρ_Y∈ C^2_∗(Y)⊂ C^1,α(Y) on the Riemannian manifold Y such that 
                
                    dμ = ρ_Y dV_Y/Vol(Y).
                
                Let (Y,μ) ∈ F M M(n,Λ,D) and ρ_Y be the density function on Y. Then
                
                    logρ_Y_C^0(Y)≤ C(n,Λ,D,(Y),i_0),
                
                    ρ_Y_C^0,1(Y)≤ C(n,Λ,D,(Y),i_0),
                
                where i_0 is the injectivity radius bound for Y given by (<ref>).
                
                
                Let π:Y→ X=Y/O(n) be the natural projection. Recall <cit.> that
                
                    π_∗(μ)=μ,
                
                which implies that 
                
                    π_∗: L^2_O (Y,μ)→ L^2(X,μ)
                
                is an isometry, where L^2_O (Y,μ) is the O(n)-invariant subspace of L^2 (Y,μ).
                Thus for any u∈ L^2(X,μ) and any open subset U⊂ X, we have
                
                    u_L^2(U,μ)=u∘π_L^2(π^-1(U),μ).
                
                
                The weighted Laplacian Δ_X in a local coordinate (x^1,⋯,x^d) of X^reg has the form
                
                    Δ_X u    =   1/ρ_X |g|^1/2∂/∂ x^j( ρ_X |g|^1/2g^jk∂/∂ x^ku ) 
                       =    g^jk∂^2/∂ x^j ∂ x^k u + lower order terms.
                
                The weighted Laplacian Δ_Y on Y has the same form with ρ_X replaced by ρ_Y. By <cit.>, the L^2(Y,μ)-normalized eigenfunctions of Δ_Y consist of O(n)-invariant part, denoted by ϕ_j^O, and the orthogonal complement, denoted by ϕ_j^⊥. Moreover, ϕ_j=π_∗(ϕ_j^O)=ϕ_j^O ∘π^-1 is an L^2(X,μ)-normalized eigenfunction of Δ_X due to (<ref>), and their corresponding eigenvalues coincide: λ_j=λ_j^O. As a consequence, the eigenfunctions {ϕ_j}_j=1^∞ of Δ_X form an orthonormal basis of L^2(X,μ).
                
                
                
                
                § QUANTITATIVE UNIQUE CONTINUATION
                
                
                We consider the quantitative unique continuation on the Riemannian manifold Y stated in Theorem <ref>, with observation in an open subset V ⊂ Y, where the Laplacian is weighted with the density function ρ_Y as discussed in the previous section. Recall that ρ_Y∈  C^1,α(Y) and ρ_Y>0. Furthermore the density function ρ_Y has uniformly bounded C^0,1-norm by Lemma <ref>. 
                
                
                Consider the wave equation with the weighted Laplacian Δ_Y on Y,
                
                    (∂_t^2-Δ_Y) u= f,
                
                where Δ_Y is given by (<ref>) with ρ_X replaced by ρ_Y.
                In view of the bounded geometric parameters (<ref>)-(<ref>) on Y, we apply the stable unique continuation theorem in <cit.>. 
                This is applicable because Δ_Y in local coordinates has the same principle term as the Laplace-Beltrami operator, and the density ρ_Y only appears in lower order terms with its C^1-norm uniformly bounded.
                However, before we can apply those mentioned results, we need to relax the curvature assumptions in those results where higher order curvature bounds were assumed. In the simple case of manifolds without boundary, one can apply results in <cit.> to smoothen the Riemannian metric.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                This is an improvement of <cit.> for manifolds without boundary in terms of regularity assumptions by removing the assumptions on the higher order curvature bounds in those results. The idea is to construct the desired non-characteristic domains using distance functions associated to smoothened Riemannian metrics. Recall from <cit.> that for sufficiently small ϵ depending only on n,Λ, there exist smooth Riemannian metrics g_ϵ on Y such that
                
                
                    g-g_ϵ_C^0≤ϵ, 
                       |R(Y,g_ϵ)| ≤ 2Λ^2, 
                       ∇^k R(Y,g_ϵ)≤ C(n,Λ) (k+1)!  ϵ^-k/2 ,  ∀  k∈.
                
                The metrics g_ϵ are constructed as the solution of the Ricci flow with the initial metric g.
                The sectional curvature bound (<ref>) is due to the maximum principle.
                
                
                We denote by d_ϵ the distance function of (Y,g_ϵ). 
                The conditions Vol_n(Y,g)≥ c_0 and (<ref>) imply that Vol_n(Y,g_ϵ) ≥ c_0/2 for sufficiently small ϵ.
                
                Then Cheeger's injectivity radius estimate and (<ref>) yield
                
                    inj(Y,g_ϵ) ≥ i_1=i_1(n,Λ,D,c_0).
                
                This shows that the distance function d_ϵ(·,y) of (Y,g_ϵ) is smooth on the ball B^ϵ(y,i_1)∖{y}, where B^ϵ denotes a ball of (Y,g_ϵ). Note that the two distance functions d_Y,d_ϵ are close to each other with respect to the Lipschitz distance: d_L((Y,g), (Y,g_ϵ) )<ϵ, which implies that B(x,i_1/2)⊂ B^ϵ(x,i_1) for sufficiently small ϵ.
                
                
                To obtain the bounds for higher order derivatives of d_ϵ(·,x), we need to work in harmonic coordinates instead of geodesic normal coordinates used in <cit.>, since we do not assume higher order curvature bounds for (Y,g).
                
                For closed manifolds, the C^1,α-harmonic radius can be explicitly estimated in terms of the sectional curvature bound, see <cit.>.
                
                Namely, there exists a constant r_h=r_h(n,Λ,i_1), such that one can find harmonic coordinates in a ball of radius r_h/2 on (Y,g) with the property that
                
                    1/2|v|^2≤∑_i,j=1^n g^ij v_i v_j≤ 2|v|^2  (v∈ℝ^n),    g_ij_C^1,α≤ C(n,Λ,r_h,α),
                
                for any α∈ (0,1).  
                Let us fix a point y∈ Y and a harmonic coordinate (x^1,⋯,x^n) on (Y,g) in a ball of radius r_h/2 centered at y.
                Recall that g_ϵ is the solution of the Ricci flow with initial metric g.
                Then in this harmonic coordinate, the Ricci flow equation, together with curvature bounds (<ref>), and an interpolation argument yield
                
                    (g_ϵ)_ij_C^1,α'≤ C(n,Λ,r_h,α,α'),  ∀ α'∈ (0,α).
                
                
                In the harmonic coordinate (x^1,⋯,x^n) for (Y,g) on B(y,r_h/2), 
                the covariant derivative ∇^ϵ has the following form (e.g. Chapter 2 in <cit.>, p32),
                
                    ((∇^ϵ)^2 d_ϵ(·, y))(∂/∂ x^k,∂/∂ x^l)=∂^2/∂ x^k ∂ x^l d_ϵ(·,y)-∑_i=1^n (Γ^ϵ)_kl^i ∂/∂ x^i d_ϵ(·,y),    k,l=1,⋯,n.
                
                Note that d_ϵ(·,y) is smooth on B(y,r_h/2)∖{y} in the harmonic coordinate (x^1,⋯,x^n) for (Y,g) since g_ϵ is smooth in this coordinate.
                
                
                We know that |∇^ϵ d_ϵ(·,y)|=1, and 
                
                    (∇^ϵ)^2 d_ϵ(·,y)≤ 4 r^-1,  (∇^ϵ)^3 d_ϵ(·,y)≤ C(n,R(Y,g_ϵ)_C^1)r^-2
                
                on B^ϵ(y,r)-B^ϵ(y,r/2) for sufficiently small r depending only on n,Λ,i_1 due to (<ref>). 
                The bound on the second covariant derivative is due to the Hessian comparison theorem.
                
                Then the bound on the third covariant derivative follows by differentiating the Riccati equation in polar coordinates, see e.g. <cit.>.
                
                Thus, in the C^1,α-harmonic coordinate (x^1,⋯,x^n) of (Y,g) chosen above, (<ref>), (<ref>) and (<ref>) give the bound for the C^2,α'-norm of d_ϵ(·,y) in this harmonic coordinate for any α'∈ (0,α):
                
                    d_ϵ(·,y)_C^2,α'≤ C(n,Λ,r_h,α,α',ϵ)r^-2   in  B^ϵ(y,r)-B^ϵ(y,r/2),
                
                for sufficiently small r<r_h/4 depending on Λ,i_1.
                
                Note that the sectional curvature bound (<ref>) independent of ϵ is crucial, as this yields the uniform injectivity radius bound (<ref>) and a uniform choice of r above.
                
                
                Now we can use the method in <cit.> adapted for the simpler case of the closed manifold (Y,g) in C^1,α-harmonic coordinate charts, and then apply <cit.>. 
                
                
                For closed manifolds, the proof is much less technical as the distance function d_ϵ is already smooth within half the injectivity radius.
                Set
                
                    0<h<min{1/10,T/8, i_1/10,r_h/10,π/12Λ}.
                
                Let {z_j} be an h-net in {x∈ V: d_Y(x,∂ V)≥ 2h}, and ρ_0:=min{i_1/2,r_h/2,π/6Λ}. If T≤ρ_0,
                we define the following C^2,α' function
                
                    ψ_j(x,t)=(T-d_ϵ (x,z_j))^2-t^2,
                
                and the domain
                
                    Ω_j^0:={(x,t)∈ Y× [-T,T]: ψ_j(x,t)>8T^2ϵ}-{x: d_ϵ(x,z_j)≤h/2}× [-T,T].
                
                For each z_j, we fix a C^1,α-harmonic coordinate chart (x^1,⋯,x^n) of (Y,g) in the ball B(z_j,ρ_0). By the discussions above d_ϵ(·,z_j) is smooth in this coordinate, and (<ref>) gives a bound on its C^2,α'-norm (with r=h) independent of z_j.
                Hence we see that the C^2,α'-norm of ψ_j is uniformly bounded in the harmonic coordinate chosen above in Ω_j^0 (depending on ϵ,h).
                Next, we show that ψ_j is a non-characteristic function in Ω_j^0.
                We take the gradient of ψ_j with respect to the space variable in the harmonic coordinate (x^1,⋯,x^n) of (Y,g) chosen above,
                
                    ∇_x ψ_j = 2 ( T-d_ϵ(x,z_j) ) ∇ d_ϵ(x,z_j).
                
                Using (<ref>) in the harmonic coordinate satisfying (<ref>), for sufficiently small ϵ, we have
                
                    |∇ d_ϵ(x,z_j)|    =    g^kl (∂_x^k d_ϵ) (∂_x^l d_ϵ) 
                       ≥    (g_ϵ)^kl (∂_x^k d_ϵ) (∂_x^l d_ϵ)-ϵ (∂_x^k d_ϵ) (∂_x^l d_ϵ) 
                       ≥    (1-3ϵ) |∇^ϵ d_ϵ(x,z_j)|^2=1-3ϵ,
                
                where ∇^ϵ is the gradient with respect to the metric g_ϵ.
                Thus in Ω_j^0, the principle symbol
                
                    p ((x,t),∇ψ_j )    =   |∇_x ψ_j|^2-|∂_t ψ_j|^2 
                       ≥    4(T-d_ϵ)^2(1-3ϵ)^2-4t^2 
                       ≥    4ψ_j(x,t) - 24T^2 ϵ >8T^2ϵ,
                
                where p(·,·) stands for the principle symbol of the wave operator with respect to the original metric g.
                Note that the only part changed from <cit.> in the construction above is the distance function d_ϵ, while the harmonic coordinate and the wave operator considered are still with respect to the original metric g.
                
                
                For the case of T>ρ_0,
                define a cut-off function ξ:ℝ→ℝ by
                
                    ξ(x)=(h-x)^3/h^3,  for  x∈ [0,h],
                
                and ξ(x)=0 for x>h. The function ξ(x) on negative numbers can be defined in any way so that ξ(x)≥ 1 for x<0, and ξ(x) is smooth on (-∞,h). 
                This function ξ(x) is of C^2,1 on ℝ and monotone decreasing on [0,+∞). 
                Define the function
                
                    ψ_1,j(x,t)=((1-ξ(ρ_0-d_ϵ(x,z_j)))T-d_ϵ(x,z_j))^2-t^2.
                
                This is a function of class C^2,α' with uniformly bounded norm (<ref>) in the domain
                
                    Ω_1,j^0:={(x,t)∈ Y× [-T,T]: ψ_1,j(x,t)>8T^2ϵ}-{x: d_ϵ(x,z_j)≤h/2}× [-T,T].
                
                Note that in general, the domain characterized by ψ_1,j(x,t)>8T^2 ϵ has two connected components. Here we define Ω_1,j^0 to be the connected component characterized by (1-ξ(ρ_0-d_ϵ(x,z_j)) )T-d_ϵ(x,z_j)>0. Thus, it satisfies that d_ϵ(x,z_j)<ρ_0 in Ω_1,j^0 by the definition of the function ξ, and hence Ω_1,j^0 is contained in one coordinate chart of Y times [-T,T].
                
                A similar argument as above shows that ψ_1,j is non-characteristic on Ω_1,j^0. 
                Then one can follow the method in <cit.> to construct all subsequent ψ_i,j and Ω_i,j^0, Ω_i,j. 
                Note that Theorem 1.2 in <cit.> is applicable to the present case, since the weighted Laplacian Δ_Y has coefficients g^ij and ρ_Y in C^1 with uniformly bounded norms in harmonic coordinates.
                In the end, the unique continuation propagates to following domain
                
                    Ω_ϵ(h) :={(x,t)∈ (Y ∖ V)× [-T,T]: T-|t|-d_ϵ(x,V) > √(ϵ)+h }.
                
                Since the distance functions d_Y,d_ϵ are ϵ-close in the Lipschitz distance,
                it follows that
                
                    Ω_ϵ(h) ⊃{(x,t)∈ (Y ∖ V)× [-T,T]: T-|t|-d_Y(x, V) > √(ϵ)+h +Dϵ}.
                
                Thus, fixing any choice of α,α', the theorem follows by choosing ϵ=h.
                The constants C_4,h_0 explicitly depend on n,T,Λ,D,i_1,r_h,ρ_Y_C^1, the upper bound for Vol_n (Y), and the curvature bounds in (<ref>); the constant C_5 explicitly depends on n, see Appendix A in <cit.>. 
                The volume Vol_n(Y) is bounded above by n,Λ,D.
                
                The curvature bounds in (<ref>) depend on ϵ only polynomially, which can be absorbed into the exponential term exp(h^-C_5).
                
                
                
                
                
                
                By the Sobolev embedding theorem and <cit.>, Theorem <ref> yields the stable unique continuation on the whole domain of influence.
                
                
                Let (Y,g) be a connected, closed, smooth manifold with C^2-smooth Riemannian metric g, satisfying 
                
                    dim(Y)=n,    |R(Y)|≤Λ^2,    diam(Y)≤ D,    Vol(Y)≥ c_0.
                
                Let ρ_Y∈ C^1(Y) be a density function on Y.
                Suppose u∈ H^1(Y×[-T,T]) is a solution of the wave equation (<ref>) with f=0 and ∂_t u (·,0)=0. 
                Let V⊂ Y be a connected open subset with smooth boundary.
                If the norms satisfy
                
                    u(·,0)_H^1(Y)≤ E_1,  u_H^1(V× [-T,T])≤ε_0,
                
                then for 0<h<h_0, we have
                
                    u(·,0)_L^2(Y(V,T))≤ C_4 h^-2/9exp(h^-C_5)E_1/(log(1+hE_1/ε_0)) ^1/6 +C_6 E_1 h^1/3(n+1),
                
                where the domain of influence Y(V,τ):= {y∈ Y: d_Y(y,V) <τ}. The constant C_6 explicitly depends on n,Λ,D,c_0, Vol_n-1(∂ V), i_b(V), where i_b(V) denotes the boundary injectivity radius of V.
                Furthermore, if we additionally assume T≥ r_1 for some r_1>0, then C_6 can be chosen explicitly depending only on n,Λ,D,c_0,r_1 (independent of V).
                
                
                This is a straightforward corollary of Theorem <ref> and the Sobolev embedding, due to the uniform boundedness of the Hausdorff measure of ∂ Y(V,τ) for any τ>0.
                The latter was essentially proved in <cit.>.
                The dependence on Vol_n-1(∂ V),i_b(V) is brought in when T is small. 
                
                
                If T≥ r_1, then the Hausdorff measure of ∂ Y(V,τ) does not depend on V.
                
                
                Note that Proposition 3.10 in <cit.> assumes C^1 curvature bound in the case of manifolds with boundary. For the present case of manifolds without boundary, the sectional curvature bound is enough.
                
                
                
                Later in the next section, we will apply this quantitative unique continuation result to Riemannian manifolds 
                Y ∈ F M M(n,Λ,,D). 
                From discussions in Section <ref>, we know that Y∈ F M M(n,Λ,,D) is a smooth manifold with C^4-smooth Riemannian metric tensor and the sectional curvatures are uniformly bounded (<ref>).
                Note that the constants C_4,C_6 depend on the lower bound for the volume of Y, which is given by Lemma <ref> and (<ref>).
                
                 
                
                On a Riemannian manifold (Y,μ) ∈ F M M(n,Λ,D), we have two measures: μ and the Riemannian volume element dV_Y. 
                
                Since ρ_Y_C^0 and ρ_Y^-1_C^0 are both uniformly bounded by Lemma <ref>, the norms with respect to these two measures are equivalent, considering the fact that Vol(Y) is bounded from below by Lemma <ref> and clearly from above.
                
                From now on, the norms on Y are taken with respect to μ.
                
                
                § DETERMINATION OF FOURIER COEFFICIENTS
                
                
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose Vol_d(X) ≥ v_0 where d= dim(X).
                
                Assume that the interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J for the weighted Laplacian Δ_X are given on a ball B(p,r_0) ⊂ X^reg.
                In this subsection, let 0<η< min{i_0/2,r_0/4} be a fixed small number such that η<inf_x∈ B(p,r_0/2) inj(x), where i_0 is the injectivity radius bound for Y given by (<ref>), and inj(x) for x∈ X^reg is the injectivity radius at x in the Riemannian sense. 
                
                
                Here (Y,μ) ∈ F M M(n,Λ,,D) is the Riemannian manifold stated in Theorem <ref> and π:Y→ X=Y/O(n) is the projection.
                We choose open subsets of B(p,r_0) in the following way. 
                Let {p_i}_i=1^N be a maximal η/2-separated set in B(p,r_0/2).
                
                Let {U_i}_i=1^N be disjoint open subsets containing p_i in B(p,r_0) satisfying
                
                    B(p,r_0/4) ⊂⋃_i=1^N U_i⊂ B(p,r_0/2),   diam(U_i)≤η.
                
                Assume that every U_i contains a ball of radius η/4. 
                Without loss of generality, we assume that every ∂ U_i is smoothly embedded in X and admits a boundary normal neighborhood of width η/10. 
                This is because one always has the choice to propagate the unique continuation from the disjoint balls B(p_i,η/4). An error of order η does not affect the final result. 
                
                
                
                
                Note that with the choice of U_i above, the total number N is bounded above. 
                This is because {B(p_i,η/4)}, and hence {π^-1(B(p_i,η/4))}, do not intersect, and each {π^-1(B(p_i,η/4))} contains a ball of radius η/4 in Y, where π:Y→ X is the natural projection.
                Thus the total number of balls is bounded by geometric parameters including the volume of Y, which in turn is bounded above by n,dim(Y),|R(Y)|,D, i.e.,
                
                    N≤ C(n,dim(Y),Λ,,D) η^-dim(Y).
                
                Recall that the bound for the sectional curvature R(Y) of Y is due to Lemma <ref>.
                Here we have used the condition that η<i_0/2, so that the volume of balls of radius η/4 in Y is bounded below by Cη^dim(Y).
                
                
                Let α=(α_1,⋯,α_N) with α_k∈ [η, D] be a multi-index. We define the domain of influence associated with α by
                
                    X_α:=⋃_k=1^N X(U_k, α_k)=⋃_k=1^N { x ∈ X: d_X(x,U_k) < α_k }.
                
                Denote by 𝒱_J ⊂ C^0,1(X) the function space spanned by the first J orthonormalized eigenfunctions {ϕ_j}_j=0^J.
                Given E_1,ε_1>0, we define
                
                    𝒰(E_1,ε_1):=⋂_k=1^ N {v∈𝒱_J: v∘π_H^2(Y)≤ E_1, W(v ∘π)_L^2(π^-1(U_k)× [-α_k,α_k])≤ε_1 },
                
                where W(v) is the solution of the wave equation
                
                    (∂_t^2 -Δ_Y) W=0    on Y,
                
                with the initial conditions W|_t=0=v and ∂_t W|_t=0=0.
                
                The solution operator of the wave equation commutes with π. Namely, define W(v) to be the solution of the wave equation (∂_t^2 -Δ_X) W=0 on X with the initial conditions W|_t=0=v and ∂_t W|_t=0=0. Then 
                
                    W(v)∘π=W(v∘π).
                
                This is because W(v)∘π satisfies the initial conditions and the wave equation on Y, see <cit.>, and the claim follows from the uniqueness of solution.
                
                The conditions in (<ref>) can be computed using the given interior spectral data on B(p,r_0). Namely, if v=∑_j=0^J v_j ϕ_j, then v∘π=∑_j=0^J v_j ϕ_j^O. 
                Note that ϕ_j^O=ϕ_j ∘π are orthonormal with respect to the measure μ on Y.
                Since λ_j=λ_j^O, then v ∘π_H^2(Y)=∑_j=0^J(1+λ_j^2)v_j^2.
                The wave W(v∘π) on Y is given by
                
                    W(v∘π)(y,t)= ∑_j=0^J v_j cos(√(λ_j) t) ϕ_j^O (y),    y∈ Y,   t≥ 0.
                
                Hence, due to (<ref>),
                
                    W(v ∘π)_L^2(π^-1(U_k)× [-α_k,α_k])=∑_j=0^J v_j cos(√(λ_j) t) ϕ_j _L^2(U_k× [-α_k,α_k]),
                
                where the right-hand side can be computed using the given interior spectra data.
                
                Let u∈ H_O^2(Y) satisfying
                
                    u_L^2(Y)=1,  u_H^2(Y)≤  E_0.
                
                Then for sufficiently small γ∈ (0,N^-2), we can construct a function u_0∈ H^2(X^reg) such that
                
                    u_0 |_X_α=0,     u_0 |_X^c_α+γ=u∘π^-1,    u_0_L^2(X)≤  1,
                
                    u_0∘π_H^2(Y)≤ C_0 E_0 γ^- dim(Y)-2,
                
                where α+γ=(α_1+γ,⋯,α_N+γ). 
                
                
                The proof is similar to Lemma 4.1 in <cit.>. To construct the required partition of unity, one can use the projections of partition of unity on Y. Let {x_j} be a maximal γ/2-separated set in X for γ<i_0/2, where i_0 is the constant given in (<ref>). Then {π^-1(B(x_j,γ/2))} is an open cover of Y due to d_X(x,y)=d_Y(π^-1(x),π^-1(y)).
                We take any C^2 partition of unity {χ_l} subordinate to this open cover of Y, and consider the projection to O(n)-invariant functions,
                
                    χ_l^O:=ℙ_O χ_l.
                
                The projection ℙ_O is defined as
                
                    (ℙ_O v)(y)=1/Vol(O(n))∫_O(n) v(o(y))dℋ_O(n)(o),    y∈ Y,
                
                where dℋ_O(n) denotes the Riemannian measure (or the Haar measure) of O(n).
                By <cit.>, χ_l^O∈ C^2_O(Y) (i.e., the set of O(n)-invariant C^2 functions on Y), and
                
                    ∑_l χ_l^O (y) =1/Vol(O(n))∫_O(n)∑_l χ_l(o(y)) dℋ_O(n)(o)=1.
                
                Note that the summation above is finite for every y∈ Y.
                Thus the partition of unity on X subordinate to {B(x_j,γ/2)} can be taken as 
                
                    χ_l:=χ_l^O∘π^-1∈ C^0,1(X).
                
                Then u_0 can be defined as
                
                    u_0(x):=∑_supp(χ_l)∩ X_α=∅χ_l(x) · (u∘π^-1)(x)  ,   x∈ X.
                
                The first three conditions claimed are clearly satisfied by construction, considering that u∘π^-1_L^2(X)=1 by (<ref>). 
                
                
                
                The last condition follows from the given H^2-bound and
                χ_l^O_C^2≤ Cγ^-2, if we can estimate the number of nonzero χ_l in (<ref>).
                
                It remains to show that the number of nonzero χ_l in (<ref>) for any given x is bounded.
                In fact, knowing a rough bound on the total number of χ_l is already enough for our purpose.
                Since {π^-1(B(x_j,γ/4))} do not intersect and each π^-1(B(x_j,γ/4)) contains a ball of radius γ/4 in Y, the total number of balls is bounded above by C Vol(Y) γ^-dim(Y), which is bounded by n,dim(Y),|R(Y)|,D,γ. Here we have used the condition that γ<i_0/2 in the same way as (<ref>). Compared with <cit.>, the bound we obtained here for the number of nonzero χ_l in the sum (<ref>) (for any given point x) depends on γ, which results in a much rougher H^2-bound. 
                
                
                
                
                Let u∈ H_O^2(Y), and u_0 be constructed as in Lemma <ref>. Suppose
                
                    u_0=∑_j=0^∞ d_j ϕ_j,
                
                where ϕ_j are the orthonormalized eigenfunctions of Δ_X.
                We consider
                
                    u_0:= u_0∘π=∑_j=0^∞ d_j ϕ_j^O,
                
                where ϕ_j^O are the eigenfunctions of the O(n)-invariant component Δ_O of Δ_Y. 
                Then we define 
                
                    u_J:=∑_j=0^J d_j ϕ_j,  u_J:=u_J∘π=∑_j=0^J d_j ϕ_j^O.
                
                Let E_1=C_0 E_0 γ^- dim(Y)-2. For any ε_1>0, there exists sufficiently large J such that u_J∈𝒰(E_1,ε_1), where the set 𝒰(E_1,ε_1) is defined in (<ref>).
                
                The first norm condition in (<ref>) is satisfied due to (<ref>).
                
                The second norm condition is a consequence of Lemma <ref>, finite speed propagation of the wave equation on Y, and the energy estimates, see <cit.>. 
                
                Due to Theorem <ref> and d_X(x,y)=d_Y(π^-1(x),π^-1(y)), the domains of influence on X and Y match. More precisely, define
                
                    Y(U,τ):= {y∈ Y: d_Y(y,U) <τ}.
                
                Then
                
                
                    Y(π^-1(U),τ)=π^-1( X(U,τ) ).
                
                
                We prove (<ref>) as follows. Let y∈ Y such that d_Y(y, π^-1(U))<τ. Let x=π(y). Then d_Y(π^-1(x), π^-1(U))<τ by definition, which yields that d_X(x,U)<τ. This shows that Y(π^-1(U),τ) ⊂π^-1( X(U,τ) ). 
                
                
                On the other hand, let z=π^-1(x) such that d_X(x,U)<τ, which gives d_Y(π^-1(x), π^-1(U))<τ. Thus there exists z_0∈π^-1(x) such that d_Y(z_0, π^-1(U))<τ. 
                Since X=Y/O(n), there exists o∈ O(n) such that o(z)=z_0.
                Since O(n) acts isometrically on Y, we have
                d_Y(z,π^-1(U))=d_Y(z,o^-1(π^-1(U)))<τ. 
                This proves the other direction π^-1( X(U,τ) ) ⊂ Y(π^-1(U),τ).
                
                
                
                
                
                
                
                
                Let u∈ H_O^2(Y) satisfying
                
                    u_L^2(Y)=1,  u_H^2(Y)≤  E_0.
                
                Set u=u∘π^-1∈ L^2(X).
                Let α=(α_1,⋯,α_N), α_k∈ [η, D] be given, and X_α be defined in (<ref>). 
                
                Then for any ε>0, there exists sufficiently large J,
                
                such that by only knowing the first J interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J of X and the first J Fourier coefficients {a_j}_j=0^J of u, we can find {b_j}_j=0^J and u^a=∑_j=0^Jb_j ϕ_j, such that
                
                    u^a-χ__X_αu_L^2(X,μ) <ε,
                
                where χ denotes the characteristic function.
                
                
                
                We consider the following minimization problem in 𝒰 defined in (<ref>).
                
                Let u_min∈𝒰 be the solution of the minimization problem
                
                    u_min-u_L^2(X)=min_w∈𝒰(E_1,ε_1)w-u_L^2(X),
                
                with parameters E_1,ε_1 to be determined.
                Since w∈𝒰, the norm conditions in (<ref>) yield the following by interpolation and energy estimate:
                
                    W(w ∘π)^2_H^1(π^-1(U_k)× [-α_k,α_k])   ≤   W(w ∘π)_L^2(π^-1(U_k)× [-α_k,α_k])W(w ∘π)_H^2(Y× [-α_k,α_k])
                       ≤   ε_1 · C(D) w∘π_H^2(Y)≤ C(D) E_1 ε_1.
                
                Denote U_k=π^-1(U_k).
                Then we apply our quantitative unique continuation result (Proposition <ref>), in view of the bounded geometric parameters (<ref>)-(<ref>) for (Y,μ) ∈ F M M(n,Λ,,D) and Lemma <ref>.
                
                
                By Lemma <ref>, Remark <ref>, (<ref>) and (<ref>), we have
                
                    w_L^2(X(U_k,α_k)) = w∘π_L^2(Y(U_k,α_k))< ε_2=ε_2(h,E_1,η,γ,ε_1),
                 
                for all k=1,⋯,N.
                
                
                
                
                Hence,
                
                    w_L^2(X_α) < Nε_2.
                
                Then for any w∈𝒰 and in particular for w=u_min,
                
                    w-u_L^2(X)^2    =   w-u_L^2(X_α)^2+  w-u_L^2(X_α^c)^2 
                       >   u^2_L^2(X_α)-2Nε_2 +  w-u^2_L^2(X_α^c).
                
                
                On the other hand the following estimate holds for u_J:
                
                    u_J-u_L^2(X)^2    ≤    (u_J-u_0_L^2(X) + u_0-u_L^2(X))^2 
                       ≤   u_J-u_0^2_L^2(X)+4u_J-u_0_L^2(X)  + u_0-u^2_L^2(X)
                       ≤    C λ_J^-1/2γ^-dim(Y)-2+ u^2_L^2(X_α)+ u_0-u^2_L^2(X_α+γ-X_α) ,
                
                where
                
                
                    u_0-u_L^2(X_α+γ-X_α)≤u_L^2(X_α+γ-X_α) < C γ^1/2(n+1).
                
                The first inequality above is due to the construction of u_0, cf. (<ref>).
                The last inequality above is due to (<ref>), Sobolev embedding (on Y) and
                
                    μ(X_α+γ-X_α)    =   μ(π^-1(X_α+γ)-π^-1(X_α)) 
                       =   ρ_Y/Vol(Y)Vol(π^-1(X_α+γ)-π^-1(X_α)) 
                       =   ρ_Y/Vol(Y)Vol(Y_α+γ-Y_α) < C N γ <C γ^1/2,
                
                where we have used <cit.>, (<ref>), (<ref>), the upper bound for ρ_Y, the lower bound for Vol(Y) due to Lemma <ref> and γ<N^-2. 
                In the above we denote the difference of sets by A-B:=A∩ B^c for two sets A,B.
                The domain Y_α is defined by
                
                    Y_α:=⋃_k=1^N Y(π^-1(U_k), α_k)=⋃_k=1^N { y ∈ Y: d_Y(y,π^-1(U_k)) < α_k }.
                
                The fact that 
                
                    Vol( Y(π^-1(U_k),α_k+γ)-Y(π^-1(U_k),α_k) ) <Cγ
                
                follows from <cit.>. Note that here the constant C depends on η due to the condition α_k≥η, the same as the special case of Proposition <ref>.
                
                
                
                Hence,
                
                    u_J-u_L^2(M)^2 < C λ_J^-1/2γ^-dim(Y)-2+ u^2_L^2(X_α)+ C γ^1/n+1.
                
                For sufficiently large J, we have u_J∈𝒰(E_1,ε_1) by Lemma <ref> with E_1=C_0 E_0 γ^-dim(Y)-2. This indicates that the minimizer u_min also satisfies
                
                    u_min-u_L^2(X)^2 < C λ_J^-1/2γ^-dim(Y)-2+ u^2_L^2(X_α)+ C γ^1/n+1.
                
                
                Combining the two inequalities (<ref>) and (<ref>), we have
                
                    u_min-u^2_L^2(X_α^c) < 2Nε_2 +C(Λ) λ_J^-1/2γ^-dim(Y)-2+ Cγ^1/n+1.
                
                The fact that u_min_L^2(X_α)< Nε_2 implies that
                
                    χ_X_αu-(u-u_min)^2_L^2(X)   =   u_min-χ_X_α^cu^2_L^2(X)
                       =   u_min-χ_X_α^cu^2_L^2(X_α^c)+u_min^2_L^2(X_α)
                       <    2Nε_2 +C(Λ) λ_J^-1/2γ^-dim(Y)-2+Cγ^1/n+1 +N^2ε_2^2 .
                
                
                Now observe that the Fourier coefficients of u_min is solvable given interior spectral data because it is equivalent to a polynomial minimization problem in a bounded domain in ℝ^J, see Remark <ref>. Suppose we have found a minimizer u_min=∑_j=0^J c_j ϕ_j. Since the first J Fourier coefficients of u are given as a_j, we can replace the function u-u_min in the last inequality by ∑_j=0^J a_j ϕ_j-u_min and the error in L^2-norm is controlled by λ_J^-1/2. Hence by the Cauchy-Schwarz inequality, we obtain
                
                    χ_X_αu-∑_j=0^J(a_j-c_j)ϕ_j ^2_L^2(X) < 4Nε_2+2N^2ε_2^2 +C(Λ) λ_J^-1/2γ^-dim(Y)-2+ Cγ^1/n+1,
                
                which makes u^a :=∑_j=0^J b_j ϕ_j with b_j=a_j-c_j our desired function. 
                
                For choices of parameters, making each term of the right-hand side of (<ref>) small gives choices of ε_2, γ, λ_J in this order, where η and hence N are fixed. The choice of γ gives E_1, and together with ε_2 gives choices of h, then ε_1 by Proposition <ref>, which at last determines J by Lemma <ref>.
                
                
                We apply Proposition <ref> to u being constant function 1 (i.e., the first eigenfunction ϕ_0=1 since μ(X)=1), 
                
                and thus we can approximate any μ (X_α) with explicit estimates.
                
                Let α=(α_1,⋯,α_N), α_k∈ [η, D] be given, and X_α be defined in (<ref>). 
                Then for any ε>0, there exists sufficiently large J, such that by only knowing the first J interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J of X, we can compute a number μ^a (X_α) satisfying 
                
                    |μ^a (X_α)-μ (X_α) |<ε.
                
                
                § RECONSTRUCTION OF INTERIOR DISTANCE FUNCTIONS
                
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose Vol_d(X) ≥ v_0 where d= dim(X).
                
                Given an open subset U⊂ X, the interior distance functions r_x: U→ℝ corresponding to x∈ X are defined by r_x (z)=d(x,z),  z∈ U. Denote ℛ_U(X)={r_x: x∈ X}, see e.g. <cit.> for basic properties. 
                Assume that the finite interior spectral data for the weighted Laplacian Δ_X are given on a ball B(p,r_0) ⊂ X^reg.
                In this subsection, we consider interior distance functions on B(p,r_0/2), denoted by B in short.
                Let 0<η< min{i_0/2,r_0/4} be such that η<inf_x∈ B(p,r_0/2) inj(x), and {p_i}_i=1^N be an η/2-net in B=B(p,r_0/2), where the total number N is bounded by (<ref>).
                Let {U_i}_i=1^N be disjoint open neighborhoods of p_i satisfying (<ref>).
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                Let β=(β_1,⋯,β_N), β_i∈ℕ be a multi-index. Define
                
                    X_β^∗:= ⋂_i: β_i>0{ x∈ X: d(x,U_i)∈ [β_i η-2η, β_iη) }.
                
                
                
                
                
                Under the assumptions of Lemma <ref>, we have
                
                    |μ^a (X_β^∗)-μ (X_β^∗) |<2^Nε.
                
                
                Let η∈ (0,i_0/2). We choose all multi-indices β with β_i∈ℕ, β_i ∈ (0, D/η] for all i such that 
                
                    μ^a (X_β^∗) > c_∗η^(Y),
                
                where c_∗ is a uniform constant  depending on n,Λ,D, dim(X),v_0. 
                For all such multi-indices β, define piecewise constant functions
                
                    r_β(z) := β_i η,   if z∈ U_i.
                
                Denote by ℛ^∗ the set of all such functions r_β.
                
                
                For any η>0, there exists a constant J_0 such that for any J>J_0, one can construct a finite set ℛ^∗ from the first J interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J of X such that
                
                    d_H(ℛ^∗,ℛ_B(X) ) <4η,
                
                where B=B(p,r_0/2) and d_H is the Hausdorff distance in L^∞(B).
                
                
                
                
                
                
                
                
                For any x∈ X, it is clear that there exists β_i∈ℕ,  β_i>0 such that 
                
                    B(x,η/2) ⊂ X_β^∗.
                
                It is clear that such r_β is an approximation of r_x.
                Indeed, by definition (<ref>), the fact that x∈ X_β^∗ gives |d(x,U_i)-β_iη|<2η for all i.
                Since diam(U_i)≤η, this shows for all i,
                
                    |d(x,z)-β_iη|<3η,  ∀  z∈ U_i.
                
                
                Furthermore, since π^-1(B(x,η/2)) contains a ball of radius η/2 in Y, 
                
                    μ (B(x,η/2)) = μ (π^-1(B(x,η/2)))    ≥   μ (B_Y(η/2)) 
                       =   ρ_Y/ Vol(Y) Vol(B_Y(η/2)) 
                       ≥    c η^(Y),
                
                where the last inequality used Lemma <ref>.
                We choose 
                ε_∗=c η^(Y)/2 and ε=2^-Nε_∗ in Lemma <ref>. Hence,
                
                    μ^a (X_β^∗) > μ(X_β^∗)-ε_∗≥μ(B(x,η/2)) -ε_∗≥ε_∗.
                
                
                
                On the other hand, let β such that μ^a( X_β^∗) >ε_∗, where ε_∗ and ε are chosen as above. By Lemma <ref>, we see that 
                μ( X_β^∗)>0, which means that X_β^∗ is nonempty. Choosing any x'∈ X_β^∗, the same argument as before shows that r_x' is an approximation of r_β.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                We say that an element r_β∈ℛ^∗ and a point x∈ X correspond to each other if r_β-r_x_L^∞(U)<4η.
                
                Observe that ℛ^∗ determines the distances on the η/2-net {p_i}_i=1^N within error 8η. Namely, ℛ^∗ determines
                
                    D^a(p_i,p_k):= inf_r_β∈ℛ^∗(r_β(p_i) + r_β(p_k) ),
                
                and it satisfies (see <cit.>)
                
                    |D^a(p_i,p_k)-d(p_i,p_k)|<8η,   for all   i,k=1,⋯,N.
                
                While a space X∈ M M(n,Λ,D) has curvature bounded from below in the sense of Alexandrov, it does not have curvature bounded from above in general due to Fukaya's counterexample <cit.>.
                
                In the case of orbifolds, i.e., (X)=n-1, with a lower volume bound, it is possible to have an upper bound for the sectional curvature of the regular part, if the metric on the regular part is of class C^4.
                For this reason, we impose a bound on the covariant derivatives of the curvature tensor (<ref>), and consider the  class M M(n,Λ,,D).
                
                
                
                Let X ∈ M M(n,Λ,,D). Suppose (X)=n-1 and Vol_n-1(X)≥ v_0. Then the sectional curvature of X^reg is bounded above by C(n,Λ,,D,v_0).
                
                This can be argued similar to Lemma <ref>.
                Suppose the sectional curvature is unbounded at p_k∈ X_k, then <cit.> implies that the isotropy group G_p_0 at lim_k→∞ p_k=p_0 ∈ X_0=lim_k→∞ X_k has positive dimension.
                However, since a uniform lower volume bound is assumed, the sequence of orbifolds X_k does not collapse (e.g. <cit.>) and hence X_0 is an orbifold. So the isotropy group of orbifold X_0 should have zero dimension, which is a contradiction.
                
                In the orbifold case, we apply the method in <cit.> to reconstruct the Riemannian metric of X^reg locally and further a global approximation. 
                Away from the singular set S, Cheeger's injectivity radius estimate (<cit.>) is still valid, depending on the lower volume bound and sectional curvature bound of
                X^reg (Lemma <ref>):
                
                    i(x)    ≥   min{π/√(|R(X^reg)|)),   c(n,D,|R(X^reg)|,v_0),   d(x,S) }
                       =   min{ c(n,Λ,,D,v_0),   d(x,S) } ,   x∈ X^reg.
                
                
                With sectional curvature bound (Lemma <ref>) and injectivity radius bound (<ref>), an approximation of interior distance functions as given in Lemma <ref> gives the approximate Riemannian structure locally at points bounded away from the singular set. 
                
                Namely, if we a priori know a point x_0∈ X^reg bounded away from the singular set S, we can approximately reconstruct the Riemannian metric at x_0.
                
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose (X)=n-1 and Vol_n-1(X) ≥ v_0. Let σ∈ (0,1), and r_0>0 such that B(p,r_0) ⊂ X^reg.
                Then for any sufficiently small η,
                
                there exists a constant J_0=J_0(n,Λ,,D,v_0,η) such that the following holds.
                
                Let x_0∈ X^reg satisfying d(x_0,S)>σ, and
                let r_0∈ℛ^∗ be a corresponding element to x_0 in the sense of Lemma <ref>. 
                
                
                Then there exists a basis {v_k}_k=1^n-1 in the tangent space T_x_0 X such that
                one can calculate numbers g_kl directly from the finite interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J for J > J_0 such that
                
                    |g_kl - g_kl(x_0) | < C(σ) η^1/8,
                
                where g_kl(x_0) is the metric components in the Riemannian normal coordinates associated to the basis {v_k}_k=1^n-1. The constant C(σ) depends on n,Λ,,D,v_0,r_0,σ, and C(σ)→∞ as σ→ 0.
                
                
                
                
                The idea is to apply Theorem 2.4 in <cit.>, but the issue here is that we need to stay away from the singular set. 
                Lemma <ref> shows that ℛ^∗ is a 4η-approximation of the interior distance functions on B(p,r_0/2).
                To construct the Riemannian metric at x_0, we need to use the data in the η^1/4-neighborhood (w.r.t. the ℓ^∞-norm) of r_0 in ℛ^∗, as formulated in <cit.>.
                For any r_ℓ∈ℛ^∗ satisfying 
                
                    r_ℓ-r_0_ℓ^∞<η^1/4,
                
                applying <cit.> gives
                
                    d(x_ℓ,x_0)< Cη^1/8,
                
                where x_ℓ is any point corresponding to r_ℓ.
                Note that <cit.> is valid for Alexandrov spaces with curvature bounded from below.
                This indicates
                
                    d(x_ℓ,S)>σ-Cη^1/8.
                
                In particular, x_ℓ∈ X^reg if η<Cσ^8. Due to the convexity of the regular part X^reg, any minimizing geodesic connecting x_0,x_ℓ or connecting points in B(p,r_0) does not intersect with the singular set.
                Hence, together with Lemma <ref> and (<ref>), Theorem 2.4 in <cit.> is applicable.
                The (explicit) dependence of C(σ) on σ traces back to <cit.>, which blows up as σ→ 0.
                
                
                § FINDING SINGULAR SET AND RECONSTRUCTION OF METRIC
                
                
                To reconstruct the metric structure from the given finite interior spectral data, one needs to determine if a corresponding point of an element in ℛ^∗ constructed in Definition <ref> is bounded away from the singular set S and then applies Proposition <ref>.
                We propose the following procedure.
                
                 Step 1. The first step is to find out the approximate location of the singular set S.
                We follow the notations is <cit.>.
                Let x∈ X^reg, denote by i(x,ξ) the distance from x to the cut point along the geodesic γ_x,ξ from x with the unit initial vector ξ. 
                
                Let ρ<i(x,ξ) and y=γ_x,ξ(ρ)∈ X^reg, and define
                
                    N(x,ξ;ρ,s,ε):=B(y,s+ε) ∖ B(x,ρ+s).
                
                Since d(x,y)=ρ<i(x,ξ) in our setting, the minimizing geodesic connecting x,y is unique, so we can also use the notation 
                N(x,y;ρ,s,ε)
                 to denote the set (<ref>).
                 
                Recall <cit.> that one can find the cut locus by examining if the measure of N(x,ξ;ρ,s,ε) is zero, when the complete spectral data are known. To study stability, we need a quantitative version of the lemma.
                Before reaching the cut locus, i.e., if ρ+s+ε≤ i(x,ξ), it follows from the triangle inequality that the set N(x,ξ;ρ,s,ε) contains a ball of radius ε/2 centered at γ_x,ξ(ρ+s+ε/2). Hence the same volume estimate as (<ref>) applies:
                
                    μ(N(x,ξ;ρ,s,ε))≥ c_∗ε^(Y),   if ρ+s+ε≤ i(x,ξ).
                
                
                If ρ+s>i(x,ξ), <cit.> says that the set N(x,ξ,ρ,s,ε) is empty for sufficiently small ε. The following lemma gives a uniform choice of ε, depending on the space X.
                
                
                Let ρ≥ρ_0>0, and let τ>0 be fixed.
                Then there exists a constant ε=ε(X,ρ_0,τ)>0 such that
                if a point x∈ X^reg and a unit vector ξ∈ T_x X satisfy ρ<i(x,ξ) and
                ρ+s≥ i(x,ξ)+τ, then N(x,ξ;ρ,s,ε) =∅ for all ε≤ε.
                
                The proof is done by contradiction. Suppose that there exist ε_j→ 0, and points and numbers x_j,y_j,ρ_j,s_j satisfying the conditions above, such that N(x_j,y_j;ρ_j,s_j,ε_j)≠∅. This means that there exist points q_j such that q_j∈ B(y_j,s_j+ε_j) and q_j∉ B(x_j,ρ_j+s_j), i.e.,
                
                    d(y_j,q_j)<s_j+ε_j,    d(x_j,q_j)≥ρ_j+s_j.
                
                Picking converging subsequences, denote the limit of subsequences of x_j,y_j,ρ_j,s_j by x,y,ρ,s as j→∞. Then,
                
                    d(y,q)≤ s,    d(x,q)≥ρ+s.
                
                Hence d(x,q)≤ d(x,y)+d(y,q)≤ρ+s, which shows
                d(x,q) = ρ+s and d(y,q)=s.
                In particular, d(x,q)=d(x,y)+d(y,q).
                If q∈ X^reg, then the convexity of the regular part means that minimizing geodesics between these points are in the regular part. Then the equality above implies that q is on the extension of the geodesic γ passing x,y, and the geodesic γ(t) starting from x is minimizing at least until t=ρ+s. 
                
                On the other hand, if q∈ S is a singular point, then <cit.> implies that any minimizing geodesic from y to q only intersects with S at the end point q, see <cit.>. 
                
                Then the same argument shows again that 
                q is on the extension of the geodesic γ passing x,y, and the geodesic γ(t) starting from x is defined and minimizing at least until t=ρ+s.
                
                
                However, due to the fact that the limit of minimizing geodesics is minimizing in complete length spaces and the limit of cut points is cut point, the cut distance is continuous in X. Note that one needs to take the closedness of the singular set into account. Hence ρ_j+s_j≥ i(x_j,ξ_j)+τ gives ρ+s≥ i(x,ξ)+τ for some subsequences, where τ>0 is fixed. 
                Note that the condition ρ≤ i(x,ξ) combined with the above indicates that s≥τ>0.
                Since geodesics do not minimize past the first cut point, γ(t) cannot minimize at t=ρ+s, which is a contradiction.
                
                
                
                Note that the choice of ε above depends on the space X and σ.
                To determine singular points, one can replace x by x'=γ_x,ξ(ρ/2), and repeat the procedure above. If the cut locus from x' does not move further, then it is a singular point; otherwise, it is not a singular point. 
                
                
                For y< i(x), denote by cut_x,y the cut point with respect to x along the extension of the minimizing geodesic from x to y.
                
                
                
                Let ρ≥ρ_0>0, d(x,S)>2ρ, x'=γ_x,ξ(ρ/2) and y=γ_x,ξ(ρ).
                
                Let σ>0 be fixed.
                Then there exist constants τ,ε, depending on X,ρ_0,σ, such that the following holds.
                
                Let s be the smallest number such that N(x,y;ρ,s,ε)=∅ for some 0<ε≤ε.
                If N(x',y;ρ/2,s+τ/2,ε)=∅ for some 0<ε≤ε, then d(cut_x,y,S)<σ.
                
                
                
                
                
                
                
                First, we claim the following observation:
                there exists τ=τ(X,ρ_0,σ)>0 such that if d(cut_x,y,S)≥σ, then d(cut_x',y, cut_x,y)≥τ.
                This can be seen from a straightforward argument by contradiction, using the convexity of the regular part and the fact that the limit of cut points is cut point.
                
                Let ε be the constant determined in Lemma <ref> with the fixed parameter τ/10.
                One can assume that ε<τ/10.
                By (<ref>) and Lemma <ref>, the condition that s is the smallest number satisfying N(x,y;ρ,s,ε)=∅ for some ε≤ε implies 
                
                    i(x,ξ)-ε < ρ+s ≤ i(x,ξ)+τ/10,
                
                
                
                which shows
                
                    ρ+s+τ/2+ε <  i(x,ξ)+τ.
                
                This implies that d(cut_x',y, cut_x,y)<τ and thus the lemma. This is because if the opposite is true: d(cut_x',y, cut_x,y)≥τ, then (<ref>) and
                
                    ρ/2+s+τ/2+ε< i(x,ξ)-ρ/2+τ = d(x', cut_x,y)+τ≤ d(x', cut_x',y),
                
                shows that the measure of N(x',y;ρ/2,s+τ/2,ε)≠∅ for every 0<ε≤ε is bounded below, which is a contradiction.
                
                
                
                One can easily check that Lemma <ref> is still valid if we can only test for the approximate measure μ^a, thanks to (<ref>) and Lemma <ref> (with N=2).
                This requires us to choose sufficiently small ε depending on X,σ, and hence sufficiently small η for finer slicing procedures.
                
                 Finding singular points. Let σ>0 be fixed.
                Let τ,ε>0 be the constants determined in Lemma <ref>, and choose an arbitrary ε≤ε.
                Then we choose sufficiently large J such that the interior spectral data {λ_j,ϕ_j|_B(p,r_0)}_j=0^J determine μ^a such that
                |μ^a (X_β^∗)-μ (X_β^∗) |<c_∗ε^dim(Y)/4 in Lemma <ref>, so that μ^a has enough accuracy to distinguish N(x,y;ρ,s,ε) being empty or not, thanks to (<ref>).
                We range x,y over all pair of points from the η/2-net {p_i}_i=1^N in B(p,r_0/2) ⊂ X^reg. 
                For each pair x,y, choose a point x' in the net such that d(x',x)≤ d(x,y)/2 +η and d(x',y)≤ d(x,y)/2 +η.
                
                This finds the midpoint between x and y in the net with error η.
                
                Recall that the distances between points in the net can be determined within error 8η, see (<ref>).
                We search for x,y,ρ,s with d(x,y)=ρ<i(x) according to the following criteria:
                
                
                  (1)s is the smallest number such that μ^a(N(x,y;ρ,s,ε) )<c_∗ε^dim(Y)/2;
                
                
                  (2) choose from s chosen in (1) such that μ^a(N(x',y;ρ/2,s+τ/2,ε))<c_∗ε^dim(Y)/2.
                 We denote the data x,y,ρ,s obtained above by
                
                    d_S^σ(x,y):=(ρ+s,s ).
                
                
                Note that the criterion (2) above uses x' as a point in the net η-close to the midpoint of x and y, while in Lemma <ref>, x' is the actual midpoint. This can be handled by choosing sufficiently small η, considering that the measure of the set N(x,y;ρ,s,ε) is continuous with respect to x, and one can choose sufficiently small η so that the difference of measure by perturbing x by η is less than c_∗ε^dim(Y)/8.
                
                
                 Step 2. The next step is to approximately identify the singular set in the data ℛ^∗.
                In the case of finite spectral data, the points x,x',y in Lemma <ref> need to be chosen from the finite number of points {p_i}_i=1^N, the η/2-net in B(p,r_0/2)⊂ X^reg. 
                Without loss of generality, assume that r_0 is less than the injectivity radius bound given in (<ref>).
                Note that the procedure described in Lemma <ref> stops at cut points.
                In order to find the singular set S, the points have to satisfy that minimizing geodesics passing pair of points can reach the singular set. 
                
                
                
                
                
                For any δ>0, there exists η=η(X,δ)>0 such that for any η-net Γ in B(p,r_0/2), the δ-neighborhood of ∪_x,y∈Γcut_x,y contains the singular set S.
                
                The proof is a straightforward argument by contradiction. If not, assume that there exist δ>0, and η_j→ 0, η_j-nets Γ_j in B(p,r_0/2) such that there exist z_j∈ S such that for any q_j∈∪_x,y∈Γ_jcut_x,y, one has d(z_j, q_j)>δ. 
                Pick any point x_j∈Γ_j ∩ B(p,r_0/4), and consider the minimizing geodesics [x_j z_j] from x_j to z_j.
                Then one can find converging subsequences x_j→ x, z_j→ z∈ S, and converging minimizing geodesics [x_j z_j] → [xz]. 
                Pick any y∈ B(p,r_0/2)-B(p,3r_0/8) so that y∈ [xz], and thus cut_x,y=z.
                Since Γ_j are η_j-nets with η_j→ 0, there exist y_j∈Γ_j such that y_j→ y. Since the limit of cut points is cut point, then for sufficient large j, d(z_j,cut_x_j,y_j) is small than any given δ, which is a contradiction.
                
                Let τ>0 be fixed. Then there exist constants η,ε, depending on X,τ, such that the following holds.
                
                Let Γ be any η-net in B(p,r_0/2)⊂ X^reg.
                For any z∈ S, there exist x,y∈Γ such that d(cut_x,y,z)<τ.
                Moreover, let s be the smallest number such that N(x,y;ρ,s,ε)=∅ for some 0<ε≤ε.
                Then N(x',y;ρ/2,s+τ/2,ε)=∅ for all ε≤ε.
                
                Let τ>0 be fixed. There exists δ=δ(X,τ)>0 such that for any x,y∈ X satisfying d(cut_x,y,S)<δ, then d(cut_x',y, cut_x,y)<τ/10.
                One can assume that δ<τ.
                Then by Lemma <ref>, there exists η=η(X,δ)>0 such that the δ-neighborhood of ∪_x,y∈Γcut_x,y contains S, where Γ is any η-net in B(p,r_0/2).
                This means that for any z∈ S, there exists x,y∈Γ such that d(cut_x,y,z)<δ<τ.
                Hence we have d(cut_x',y, cut_x,y)<τ/10.
                Recall from (<ref>) that ρ+s>i(x,ξ)-ε. Hence, for all ε≤ε<τ/10,
                
                    ρ/2+s+τ/2>i(x,ξ)-ε-ρ/2+τ/2 =d(x', cut_x,y)-ε+τ/2>d(x', cut_x',y)+τ/10,
                
                where ε is the constant determined in Lemma <ref> with the fixed parameter τ/10.
                Thus, Lemma <ref> shows that N(x',y;ρ/2,s+τ/2,ε)=∅ for all ε≤ε.
                
                
                
                Let σ∈ (0,1) be fixed. Setting the fixed parameter τ in Lemma <ref> to be the constant determined in Lemma <ref>, one determines a constant η.
                One can assume τ,10η≤σ without loss of generality.
                
                Let Γ={p_i}_i=1^N be any η/2-net in B(p,r_0/2)⊂ X^reg.
                
                
                
                
                
                With the data d_S^σ obtained in Step 1 (<ref>), we search for all multi-indices β∈ℛ^∗ evaluated on all pair of points x,y∈Γ in the data d_S^σ such that 
                
                    (r_β(x),r_β(y) )-d_S^σ(x,y) _ℓ^∞ <3σ.
                
                Denote by S^∗ all such multi-indices in ℛ^∗ according to (<ref>).
                Here we abuse the notation in Definition <ref> also denoting the set of multi-indices by ℛ^∗.
                We shall show that S^∗ has the following property.
                
                ∙Claim that Step 1 and (<ref>) find S^∗⊂ℛ^∗ and S⊂ X such that 
                
                    d_H(S, S)<C_sσ^1/2,  S:= ⋃_β∈ S^∗X_β^∗,
                
                for some constant C_s>1, where X_β^∗ is defined in (<ref>).
                We prove this claim as follows.
                For any z∈ S, there exists x,y∈Γ such that they generate admissible data according to criteria (<ref>) due to Lemma <ref>. Moreover, we have d(cut_x,y,z)<τ≤σ.
                Recall (<ref>) that |ρ+s-i(x,ξ)|<τ where (ρ+s,s)=d_S^σ(x,y), which implies 
                
                    |ρ+s-d(x,cut_x,y)|<τ,    |s-d(y,cut_x,y)|<τ.
                
                Hence by the triangle inequality,
                
                    |ρ+s-d(x,z)|<2σ,    |s-d(y,z)|<2σ.
                
                By definition (<ref>), one can find a multi-index β∈ℛ^∗ such that z∈ X_β^∗, i.e., r_z-r_β_L^∞(Γ)<4η, and thus this multi-index β satisfies (<ref>). This proves one direction of the claim (<ref>).
                
                For the other direction of the claim, take any point z'∈ X_β^∗ for some β∈ S^∗, i.e., satisfying r_z'-r_β_L^∞(Γ)<4η.
                By definition of S^∗, there exist p_1,p_2∈Γ in the data d_S^σ such that (<ref>) is satisfied.
                Then by triangle inequality,
                
                    (r_z'(p_1),r_z'(p_2) )-d_S^σ(p_1,p_2) _ℓ^∞ <3σ+4η.
                
                Since d(p_1,cut_p_1,p_2)=d(p_1,p_2)+d(p_2,cut_p_1,p_2), combining with (<ref>), we have
                
                    |d(z',p_1)-d(z',p_2)-d(p_1,p_2) |<9σ,
                
                where we have used η≤σ/10.
                
                Thus, applying <cit.> gives
                
                    d(z',cut_p_1,p_2) ≤ |d(cut_p_1,p_2,p_2)-d(z',p_2)|+C σ^1/2 < Cσ^1/2+4σ+4η.
                
                On the other hand, by Step 1 criteria (<ref>) and Lemma <ref>, there exists z∈ S such that d(z, cut_p_1,p_2)<σ. Thus, we have
                
                    d(z',z)≤ d(z',cut_p_1,p_2)+d(z,cut_p_1,p_2) < Cσ^1/2+5σ+4η,
                
                which concludes the proof of the claim (<ref>).
                
                
                 Step 3. We search for all β∈ℛ^∗ such that
                
                    η·β-β_S_ℓ^∞> 3C_s σ^1/2,   for all β_S∈ S^∗,
                
                where S^∗ is chosen in Step 2 according to (<ref>).
                
                Denote by Γ={p_i}_i=1^N the η/2-net in B(p,r_0/2).
                For any point x_0∈ X_β^∗ corresponding to such a multi-index β, and for any z'∈ X_β_S^∗ for an arbitrary β_S∈ S^∗, we know
                
                    r_x_0-r_z'_ℓ^∞(Γ)   ≥   r_β-r_β_S_ℓ^∞-r_x_0-r_β_ℓ^∞-r_z'-r_β_S_ℓ^∞
                       >   η·β-β_S_ℓ^∞-8η > 3C_s σ^1/2-8η.
                
                This implies that there exists a point in Γ, say p_1∈Γ, such that 
                
                    |r_x_0(p_1)-r_z'(p_1)|>3C_s σ^1/2-8η,
                
                which yields d(x_0,z')>3C_s σ^1/2-8η by the triangle inequality.
                Since z' is chosen arbitrarily in S, it follows from (<ref>) that 
                
                    d(x_0,S)≥ d(x_0,S)-C_sσ^1/2 > 3C_s σ^1/2-8η-C_sσ^1/2 >σ^1/2,
                
                which gives a lower bound away from the singular set S. Thus, we can apply Proposition <ref> to recover the Riemannian metric locally near x_0. Note that the constant in Proposition <ref> depends on σ, and goes to infinity as σ→ 0.
                
                
                 Step 4. The last step is to show that the set of points not corresponding to multi-indices picked up in Step 3 is small in σ.
                Namely, for all other multi-indices β not picked up in Step 3, there exists β_S∈ S^∗ such that 
                
                    η·β-β_S_ℓ^∞≤ 3C_sσ^1/2.
                
                This means that for any x∈ X_β^∗ and any z∈ X_β_S^∗∖ S ⊂S, we have for all i=1,⋯,N,
                
                    |d(x,p_i)-d(z,p_i)|<3C_sσ^1/2+8η,
                
                where {p_i}_i=1^N is the η/2-net in B(p,r_0/2).
                Take any point in {p_i}_i=1^N, say q_1, and consider the minimizing geodesic [q_1 z]. Since {p_i} is an η/2-net, one can find another point in {p_i}, say q_2, such that 
                
                    |d(z,q_1)-d(z,q_2)-d(q_1,q_2)|<η.
                 Then the triangle inequality shows
                
                    |d(x,q_1)-d(x,q_2)-d(q_1,q_2)|<6C_sσ^1/2+17η.
                
                Thus, applying <cit.> gives d(x,z)<Cσ^1/4.
                Then by (<ref>) we have
                
                    d(x,S)≤ d(x,z)+d(z,S) < Cσ^1/4+C_sσ^1/2.
                
                From the constructions in Step 1-4 above, we arrive at the following result.
                
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose (X)=n-1 and Vol_n-1(X) ≥ v_0.
                Let σ∈ (0,1), and U⊂ X^reg be an open subset containing a ball B(p,r_0).
                Then for sufficiently small η, there exists a constant J_0 such that the following holds.
                
                
                
                
                (1) The finite interior spectral data {λ_j,ϕ_j|_U}_j=0^J for J > J_0 determines a C_0η^1/2-net {x_i}_i=1^I in X, where the constant C_0 depends on n,Λ,D.
                
                
                (2) The finite interior spectral data also determines a subset of the net, say {x_i}_i=1^I_0 for some I_0<I, with the following property:
                there exists a basis {v_k}_k=1^n-1 in the tangent space T_x_i X such that
                
                one can calculate numbers g_kl directly from the finite interior spectral data {λ_j,ϕ_j|_U}_j=0^J for J > J_0 such that
                
                    |g_kl - g_kl(x_i) | < C'_1(σ) η^1/8,   for   i=1,⋯,I_0,
                
                where g_kl(x_i) is the metric components at x_i in the Riemannian normal coordinates associated to the basis {v_k}_k=1^n-1.
                The index I_0 depends on σ,η.
                
                The constant C'_1(σ) depends on n,Λ,,D,v_0,r_0,σ, and C'_1(σ)→∞ as σ→ 0.
                
                
                (3) The remaining part {x_i}_i=I_0+1^I of the net satisfies d(x_i,S)<C_3σ^1/4 for i=I_0+1,⋯,I, where the constant C_3 depends on n,Λ,D.
                
                
                
                
                
                By Proposition <ref>, the finite interior spectral data for sufficiently large J determine a 4η-approximation ℛ^∗ of the interior distance functions. 
                We consider the set of points {x_i}_i=1^I in X corresponding to the elements in ℛ^∗.
                Then applying <cit.> yields that {x_i}_i=1^I is a C_0η^1/2-net in X. Note that the Proposition 5.2(2) in <cit.> applies to Alexandrov spaces with curvature bounded from below. 
                This proves the first statement.
                The second statement follows from our constructions in Step 1-4 above and Proposition <ref>, where the subset {x_i}_i=1^I_0 of the net corresponds to elements in ℛ^∗ picked up in Step 3 (<ref>).
                The third statement is proved in Step 4.
                
                Let (X,μ) ∈ M M(n,Λ,,D) and p∈ X^reg. Suppose (X)=n-1 and Vol_n-1(X) ≥ v_0.
                Let σ∈ (0,1), and U⊂ X^reg be an open subset containing a ball B(p,r_0).
                Then for sufficiently small η, there exists a constant J_0 such that the following holds.
                
                The finite interior spectral data {λ_j,ϕ_j|_U}_j=0^J for J > J_0 determines a subset {x_i}_i=1^I_0 of a C_0η^1/2-net {x_i}_i=1^I in X and numbers d_i,i' for i,i'∈{1,⋯, I_0}, such that
                
                    |d_i,i'-d(x_i,x_i')|< C_1(X,σ) η^1/8,   for   i,i'=1,⋯,I_0.
                
                Moreover, the remaining part {x_i}_i=I_0+1^I of the net is in the C_3σ^1/4-neighborhood of the singular set of X.
                
                Theorem <ref> gives an approximation of the Riemannian metric in normal coordinates near {x_i}_i=1^I_0, where {x_i}_i=1^I_0 is the subset of the net in X satisfying d(x_i,S)>σ^1/2 determined in Step 3. Furthermore, since x_i(i=1,⋯,I_0) is bounded away from the singular set, given an element r∈ℛ^∗ close to the element corresponding to x_i (w.r.t. the ℓ^∞-norm), following the method in <cit.> also gives an approximation of the coordinate of r. Hence one obtains an approximation of d(x_i_1,x_i_2) when x_i_1,x_i_2 are close to each other, see <cit.>.
                
                For global approximations, the issue is that a shortest path between two points x,y∈ X^reg may come close to the singular set S.
                Suppose that x,y∈ X^reg satisfy d(x,S)>σ^1/2 and d(y,S)>σ^1/2.
                Observe that there exists a constant σ=σ(X,σ) such that any shortest path between x,y is bounded away from S by σ. This follows from a straightforward compactness argument in the space X, using the convexity of the regular part.
                
                
                Now applying Theorem <ref> while replacing the fixed parameter σ there by (σ/2C_3)^4, one sees from Theorem <ref>(3) that near any shortest path between x,y, there exist points in the net (depending on σ) at which one can construct the Riemannian metric locally as in Theorem <ref>(2).
                The constant C_1 now depends on σ which depends on X,σ.
                Then following the global constructions in <cit.> gives the estimate for the distances between points in the net bounded away from S by σ^1/2.
                
                The points x_i in the net satisfying d(x_i,S)>σ^1/2 can be determined by the given interior spectral data,
                
                using Step 3 (<ref>). 
                
                However, note that this procedure does not guarantee to find all points x_i satisfying d(x_i,S)>σ^1/2.
                The remaining part of the net is in the C_3σ^1/4-neighborhood of the singular set, as argued in Step 4.
                
                
                
                
                At last, we prove Theorem <ref>.
                
                
                This is a direct consequence of Theorem <ref>.
                The space X∖ S_σ,δ is the C_0η^1/2-neighborhood of {x_i}_i=1^I_0, and the finite space is X={x_i}_i=1^I_0.
                Using the numbers d_i,i' determined in Theorem <ref> from the given interior spectral data, one can construct a metric d_X on X such that 
                
                    |d_X(x_i,x_i')-d_i,i'|<2C_1 η^1/8,   for   i,i'=1,⋯,I_0.
                
                The existence of such a metric d_X is guaranteed by Theorem <ref>.
                As the space X is finite, this is equivalent to finding a solution of a system of linear inequalities with I_0(I_0-1)/2 variables.
                
                
                Then (X,d_X) is a 4C_1η^1/8-Gromov-Hausdorff approximation of X∖ S_σ,δ, equipped with the restriction of the metric d of X on X∖ S_σ,δ.
                Moreover, S_σ,δ is contained in the C_0η^1/2-neighborhood of the remaining part {x_i}_i=I_0+1^I of the net in X, which is at most C_3σ^1/4+C_0η^1/2 distance away from the singular set.
                For the relation between parameters η and δ, one needs an explicit computation from the very beginning of the reconstruction.
                The three logarithms trace back to the logarithmic and exponential dependence of constant in Theorem <ref>, and to Lemma <ref> where all N points {p_k}_k=1^N are used in the reconstruction, see Lemma <ref> and (<ref>).
                One can refer to similar computations in the Appendix of <cit.>.
                
                When we are given the complete interior spectral data {λ_j,ϕ_j|_U}_j=0^∞,
                where U⊂ X is an open and non-empty set, the uniqueness result follows from the fact that the interior spectral data
                determines the heat kernel H(x,y,t) for t>0, x,y∈ U, and applying the result of <cit.> that the  
                inverse problem for the heat kernel is uniquely solvable.
                
                
                § APPLICATIONS OF GEL'FAND'S PROBLEM FOR COLLAPSING MANIFOLDS
                
                
                 §.§ Manifold learning
                
                
                Gel'fand's  problem <ref> is encountered in manifold learning, also called dimensionality reduction in data science, see e.g. <cit.>.
                In manifold learning one is given a point cloud, that is, a finite subset {y_j}_j=1^N⊂^n of 
                 N points that are close to a d-dimensional submanifold M_0 in an n-dimensional Euclidean  space, where n is
                 much larger than d<cit.>. Then, the goal in learning is to find a manifold M⊂^m 
                 that is diffeomorphic to M_0 so that d<m<n. In a more general problem, called the geometric Whitney
                 problem <cit.>, one is given a possibly discrete metric space (Y,d_Y) (that is, not necessarily a subset of an Euclidian space), and the task is to construct a smooth Riemannian manifold M
                that is close to Y in the Gromov-Hausdorff sense.
                 The manifold learning problem is encountered in situations where a given data set 
                 is assumed to have an intrinsic structure close to a smooth manifold.
                 An example of the above problems is encountered in
                the cryogenic electron-microscopy (Cryo-EM), see <cit.>, that is, a method to find
                the structure of large 
                molecules or viruses.
                This method is applied to a large number of identical  
                samples of the  molecules 
                frozen in vitreous water.
                By applying electron-microscopy to  a thin slice of the ice containing the molecules, 
                 one obtains 
                a large number of images of the molecule, all taken from unknown directions.
                  For example, the images of the COVID-19 viruses (including their  spike proteins)
                 were obtained using this technique. The 2017 Nobel Prize in Chemistry was awarded for developing methods for this problem.
                 Mathematically, 
                each measured image of individual viruses consists of m× m pixels that correspond to points in ^n where n=m^2.
                The number of images of the virus particles, N, is quite large and all images are taken of identical objects (e.g. a virus) but
                from unknown directions. Mathematically speaking, the unknown directions correspond to
                an element in the 3-dimensional rotation group  SO(3). Thus the set of images
                can be regarded as a set of N data points y_1,…,y_N in a high dimensional Euclidean
                space ^n that lie approximately on an orbit
                of the 3-dimensional rotation group  SO(3). The individual images are 
                very noisy
                and thus the imaging problem corresponds to finding an approximation of a  submanifold M_0 in ^n that
                is diffeomorphic to SO(3), when we are given y_j=x_j+ξ_j, j=1,2,…,N,  where x_j are randomly sampled
                points on M_0 and ξ_j∼ N(0,σ^2I) are Gaussian measurement errors.
                In addition to Cryo-EM, manifold learning is widely applied in statistics, visualization of data, astronomy  and chemistry, see
                <cit.>.
                
                
                
                An extensively studied method in manifold learning is the diffusion maps (or the spectral embedding) introduced
                 by Coifman and Lafon, see <cit.>.
                In this method one considers the data points X={x_j}_j=1^N⊂^n, which are sampled from a manifold M_0 as the nodes of a graph, and uses Euclidean
                distances between these points to
                compute a kernel function k_t:X× X→ that approximates the heat kernel of the graph Laplacian
                of X at a small time t. The values k_t(x_j,x_j') can be considered as an
                approximation the values of the heat kernel H(x_j,x_j',t) of the manifold M_0 at the points 
                x_j,x_j'∈ M_0 at the time t. The values of the kernel h_t(x_j,x_j')=H(x_j,x_j',t) at the times t∈{t_1,t_2,…,t_l}⊂_+
                correspond to the pointwise heat data, PHD, considered in <cit.>.
                
                
                Then, in the diffusion map 
                  algorithm one continues by computing the first J eigenfunctions ϕ_j(x)
                of the integral operator having the kernel h_t(x,x'), x,x'∈ X, and uses the eigenfunction evaluation map,
                defined as
                
                    Φ^(K,J):X→^J,
                          Φ^(K,J)(x)=(ϕ_K(x),ϕ_K+1(x),…,ϕ_K+J-1(x) ),
                
                to construct the set Φ^(K,J)(X)⊂^J that  approximates the Riemannian manifold M_0. 
                Closely related to the present paper and the part I of the paper <cit.> is the problem of reconstructing an approximation
                of the Riemannian manifold M_0 in the case when we are given a spatially restricted data. In this case, one has only the local finite interior spectral data,
                that is,  λ_j and ϕ_j(x_k),
                where j≤ J, k≤ K and J,K∈ℤ_+, or the
                local pointwise heat data
                h_t_ℓ(x_k,x_k')=H(x_k,x_k',t_ℓ), where  ℓ≤ L and the sample points {x_k: k=1,2,…,K} are an -dense set of points in a possibly small ball B(x_0,r) of M_0, not an -dense set of points in the whole manifold M_0. The points in the ball B(x_0,r) can be considered as spatially restricted marker points.
                
                The collapsing of manifolds is encountered in manifold learning when
                the data set has several degrees of freedom in different scales <cit.>. As an motivational example,
                one can consider a d-dimensional submanifold M_0⊂^n that is an embedded image of N_0× S^1
                of a map f:N_0× S^1→^n, where N_0 is a (d-1)-dimensional manifold and 
                
                f(z,θ)=f_0(z)+(cosθ ν_1(z)+sinθ ν_2(z) ),   z∈ N_0, θ∈ [0,2π],
                
                where f_0:N_0→^n is a smooth embedding, and ν_1(z) and ν_2(z) are two orthogonal unit
                vectors normal to f_0(N_0) at the point f_0(z). The image f(N_0× S^1) models
                a data set that is determined by a (d-1)-dimensional variable z∈ N_0 and an additional 
                degree of freedom θ∈ S^1 that has a small-scale effect. In the example above on  Cryo-EM,
                such  an almost collapsed case is encountered when a large protein (e.g. a virus)
                is connected to a small part by a chemical bond that allows a rotation. 
                An example of such a bond is the peptide amide linkage between two parts of the molecule, 
                see <cit.> for an overview on the degree of freedom in
                the 3-dimensional structure of proteins.
                In this case, the problem of finding the manifold for the possible EM-images of the molecule means finding a manifold M diffeomorphic to SO(3)× S^1
                that is almost collapsed to a lower dimensional manifold.
                In data science, the collapsing of manifolds is closely related to multi-scale models
                where the intrinsic dimension of the data set is modelled by a function that depends on the
                scale (i.e. the size of window) in which the data set is approximated by a smooth structure, see <cit.>.
                Next, we consider a numerical example in the case when we are given the heat kernel data, or distance data, on the whole manifold.
                
                
                
                
                  §.§.§ A numerical example on manifold learning with complete
                Gel'fand data
                
                
                As a simple example, let us consider the case when the 
                data are given on the whole manifold, that is, on the ball B=B_M(x_0,R_0)
                where R_0>(M) so that B=M.  In this case, we demonstrate numerically two methods to learn
                a manifold, where one method uses the complete Gel'fand data and the other method is based on an approximation of the heat kernel.
                The methods developed in this paper for the local Gel'fand data are not yet numerically implemented or analyzed.
                
                We consider  a manifold M that is the torus M_R,r=𝕊^1_R×𝕊^1_r, endowed with the flat product metric, where 𝕊^1_R is the circle of radius R. In the numerical test, the larger radius is R=10 and the smaller radius is r=1/2.
                Due to visualization purposes, we also consider an embedded
                2-dimensional torus
                in ^3,  see Figure <ref> (left), given by
                      M^emb_R,r=
                
                      {((R+r cos(s_1))cos(s_2),(R+rcos(s_1)) sin(s_2), rsin(s_1) )∈^3: s_1,s_2∈ [0,2π] }.
                
                As the ratio R/r is large, the intrinsic metric of the torus M^emb_R,r inherited from ^3 is relatively close to the flat metric of M_R,r.
                
                
                
                We consider data points X={x_j∈ M_R,r: j=1,2,…,N} with N= 2048, where 
                x_j are points sampled from the torus M_R,r. The points x_j could be
                sampled randomly on M_R,r, but to make the visualization of the situation
                clearer, we sample x_j  randomly on a closed geodesic γ_x_0,ξ_0,
                where, in the (s_1,s_2)-coordinates, 
                
                
                x_0=(0,0)
                and ξ_0=1/m_s_1+_s_2 with m=35. In other words, the sample points are x_j=γ_x_0,ξ_0(τ_j), where τ_j are    independent samples from
                the uniform distribution on [0,2π]. The geodesic γ_x_0,ξ_0 on M_R,r is visualized 
                in Figure <ref> (left)
                as a helical curve γ^emb on M^emb_R,r,
                
                    γ^emb(s)=(
                    (R+r cos(ms))cos(s),(R+rcos(ms))sin(s), rsin(ms))∈^3.
                
                
                
                 
                
                 We consider the following two cases:
                
                  (C1) (Two-dimensional representation of the manifold using the heat kernel)
                We assume that we are given the Gel'fand data, that is, the values of the heat kernel h_t(x_i,x_j),i,j=1,…,N,t=100.
                
                 We compute a 2-dimensional representation of the manifold using
                 the heat kernel of the  manifold. The eigenfunctions of the Laplacian
                 operator coincide with the eigenfunctions of the integral operator defined by the heat kernel h_t(x,y),
                and we use  these eigenfunctions to represent the manifold M_R,r  in ^2. 
                 
                 
                 
                  (C2) (Two-dimensional representation of the manifold using an exponential kernel)
                
                
                
                
                
                
                We assume that we are given the geodesic distances  d_M_R,r(x_i,x_j),i,j=1,…,N, on the flat torus M_R,r.
                
                Using these data we compute the exponenal kernel  k_t(x_i,x_j) defined in (<ref>),
                that is, 
                 an approximation of the heat kernel of the  manifold. Using the eigenvectors of the matrix
                 associated to  a weighted version of the exponential kernel, we compute
                 a 2-dimensional representation of the manifold M_R,r  in ^2. 
                
                
                
                
                
                Observe that when r is small, the torus M_R,r
                is `almost collapsed' to the limit space that is the circle M_R=𝕊^1_R. 
                
                
                In the case (C1) we use a version of the diffusion map algorithm of
                 Coifman and Lafon, see <cit.>, which we modify below
                 so that it uses 
                the heat kernel of the manifold.
                
                
                
                As data, we use the point values of the heat kernel h_t(x_i,x_j) of the Riemannian manifold (M_R,r,g),
                where g=ds_1^2+ds_2^2 is the flat product metric of the torus and t=100. 
                
                On the flat torus M_R,r, the values of the heat kernel are computed using the formula 
                      h_t(x_i,x_j)=
                
                      ∑_ℓ_1,ℓ_2∈1/4π texp(-1/4t((s_1(x_i)-s_1(x_j)-Rℓ_1 )^2+(s_2(x_i)-s_2(x_j)-rℓ_2)^2 )),
                
                where 
                
                
                
                s_1(x) and s_2(x) denote the s_1 and s_2
                coordinates of the point x∈ M_R,r. 
                (We note that the flat torus was chosen to be our example as its heat kernel can be easily computed 
                using formula (<ref>).)
                First, we define the weighted kernel
                h̃_t(x_i,x_j)=1/h(x_i)h_t(x_i,x_j),   
                h(x_i)=∑_j=1^N h_t(x_i,x_j).
                
                Second, 
                 we compute the singular value decomposition 
                of the matrix (h̃_t(x_i,x_j))_i,j=1^N,
                given by h̃_t(x_i,x_j)=∑_p=1^N ω_pϕ_p(x_i)ψ_p(x_j),
                where v_p=(ϕ_p(x_i))_i=1^N and w_p=(ψ_p(x_j))_j=1^N are orthogonal
                vectors in ^N, and ω_p+1≥ω_p≥ 0.
                We define the map Φ^(K,J), see
                (<ref>), using  the right eigenfunctions ϕ_p(x),
                where p=K,K+1,…, K+J-1,  in the singular value decomposition. 
                Finally, a low-dimensional approximation
                of the manifold M_R,r is computed by 
                the image of the 2-dimensional eigenfunction map
                Φ^(2,2)(X)⊂^2. The result is shown in Figure <ref> (center).
                
                
                
                In the case (C2) we implement a simplified version of the algorithm used above in the case (C1), where
                the heat kernel h_t:X× X→ is replaced by an exponential kernel function k_t:X× X→.
                We assume that we are given the geodesic distances d_M_R,r(x_i,x_j), i,j=1,…,N, and 
                we compute the exponential kernel using the formula
                  
                 k_t(x_i,x_j)=1/4π texp(-d_M_R,r(x_i,x_j)^2/4t)
                 
                  with t=100. 
                  We note that this algorithm is a modified version of the classical diffusion map algorithm where we 
                  use the geodesic distances on M_R,r instead of  the Euclidean distances on M^emb_R,r⊂^3.
                  The kernel  k_t(x_i,x_j) can be considered as an  approximation of the heat kernel of the manifold, as due to Varadhan's formula, (<ref>) is the leading order asymptotics of the heat kernel as t→ 0.
                Then, the diffusion map algorithm  is performed by computing the right eigenfunctions ϕ_j(x) in the singular value decomposition 
                of the matrix (k̃_t(x_i,x_j))_i,j=1^N, where
                k̃_t(x_i,x_j)=1/k(x_i)k_t(x_i,x_j),   
                k(x_i)=∑_j=1^N k_t(x_i,x_j).
                
                 Finally, a low-dimensional approximation
                of the manifold M_R,r is given by 
                the image of the 2-dimensional eigenfunction map
                Φ^(2,2)(X)⊂^2, see Figure <ref> (right).
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                In both cases (C1) and (C2), the image of the 2-dimensional eigenfunction map
                Φ^(2,2)(X)⊂^2 in Figure <ref> is topologically close to a circle, that is, M_rep:=Φ^(2,2)(X) is an approximation
                of the limit space M_R. 
                
                We point out that in the earlier sections of the paper we have studied the local version of
                Gel'fand's  problem <ref>, where 
                the heat kernel H(x_j,x_j',t_ℓ) are  given at points x_j that do not fill
                the whole manifold M but only fill a possibly small metric ball B=B(x_0,R_0)⊂ M with R_0<(M), that is, {x_j: j=1,…,N}⊂ B.
                The data missing from the points x∈ M∖ B(x_0,R_0) are compensated by measuring
                the heat kernel at several times t_ℓ>0, or alternatively, measuring
                a large number eigenfunctions ϕ_j on B(x_0,R_0) and the corresponding eigenvalues,
                see Theorem <ref>.
                
                
                
                
                
                  
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                 §.§ Collapsing manifolds in physics
                
                
                
                In modern quantum field theory,
                in particular string theory, one often models  the Universe
                as a high-dimensional, almost collapsed manifold.
                This type of considerations started from the Kaluza-Klein theory in 1921
                in which the 5-dimensional Einstein equations are considered on ^4× S^1(), that is, the
                Cartesian product of the standard 4-dimensional space-time with the Minkowski metric and 
                a circle S^1() of radius . As → 0, the 5-dimensional Einstein equation yields
                a model containing both the Einstein equation and Maxwell's equations.
                In this subsection we shortly review this and discuss its relation to the considerations
                presented in the main text of this paper.
                
                
                
                
                First, let us begin with a basic example of collapsing, as illustrated in Figure <ref>.
                We start with a product manifold D^2 × [0, ϵ], where D^2 is the unit disk in ^2 with the Euclidean metric.  Consider the action of the finite group _m,  
                 m ≥ 2, on D^2 by a rotation of angle 2 π/m around the origin. Then we define M_ϵ by identifying points 
                 ( x, 0) ∈ D^2 ×{0} with (e^2π i/m· x, ϵ) ∈ D^2 ×{ϵ} in the manifold D^2 × [0, ϵ],
                 where e^2π i/m· x stands for the rotation by the angle 2 π/m.
                 Observe that this gives rise to closed vertical geodesics of length m ϵ except for the points corresponding to the origin.
                
                 When ϵ→ 0, M_ϵ collapses to a 2-dimensional orbifold X=D^2/_m, which has conic singular point only at the origin.
                
                
                
                In the Kaluza-Klein theory, one starts with a 5-dimensional
                manifold N=^4× S^1 with the metric ĝ=ĝ_jk,
                j,k=0,1,2,3,4, which is a Lorentzian metric of the type (-,+,+,+,+). The "hat"
                marks the fact that ĝ is defined on a 5-dimensional manifold.
                 We use on N the coordinates x=(y^0,y^1,y^2,y^3,θ),
                where θ is consider as a variable having values on [0,2π].
                Let us start with a background metric (or the non-perturbed metric)
                
                    ĝ_jk(x)dx^jdx^k=
                    η_νμ(y,θ)dy^ν dy^μ +^2dθ^2,
                
                where >0 is a small constant and ν,μ,j,k are summation indices taking
                values ν,μ,j,k ∈{ 0,1,2,3} . Here, we consider
                the case when [η_νμ(y,θ)]_ν,μ=0^3=(-1,1,1,1).
                Next we consider perturbations of the background metric ĝ in the following form
                
                    ĝ_jk(x)dx^jdx^k=e^-σ/3(e^σ (dθ +κ A_μ dy^μ)^2+g_νμdy^ν dy^μ),
                
                where κ>0 is a constant, σ=σ(y,θ) is a function close to the constant c():=3 log,
                the 1-form A=A_μ(y,θ) dy^μ is small and 
                 g_νμ(y,θ) is close to e^σ(y,θ)/3η_νμ(y,θ).  
                
                Next, assume that ĝ_jk satisfies the 5-dimensional Einstein equations
                
                
                
                
                
                
                
                    Ric_jk(ĝ)-1/2(ĝ^pqRic_pq(ĝ)) ĝ_jk=0   in ^4× S^1,
                
                and write g=g_μν(y,θ),A=A_μν(y,θ) and σ=σ_μν(y,θ) in terms of Fourier series,
                
                g_μν=∑_m=-∞^∞ g_μν^m(y)e^imθ,  
                A_μ=∑_m=-∞^∞ A_μ^m(y)e^imθ,  σ=∑_m=-∞^∞σ^m(y)e^imθ.
                
                Here, the functions A_μ^m(y) and σ^m(y)  correspond to some physical fields. When  
                is small, the manifold N
                is almost collapsed in the S^1 direction and all these functions
                with m≠0 corresponds to physical fields (or particles) of a very high energy which do not
                appear in physical observations with a realistic energy. Thus
                 one considers only the terms m=0 and after suitable approximations
                  (see <cit.> and <cit.>) 
                one observes that  
                 the matrix g_μν^0(y), considered as a Lorentzian metric
                on M=^4, the 1-form A^0(y)=A_μ^0(y)dy^μ,
                and the scalar function ϕ(y)=1/√(3)σ^0(y)
                
                
                
                
                
                
                
                
                
                
                satisfy
                      Ric_μν(g^0)-1/2((g^0)^pqRic_pq(g^0))  g^0_μν=T_μν   in ^4,
                      T_μν= κ^2 e^√(3) ϕ((g^0)^βF_νF_βμ-
                1/4 (g^0)^β(g^0)^γδF_γF_βδg^0_μν)
                
                      +∇_μϕ ∇_νϕ-
                1/2((g^0)^β∇_ϕ ∇_βϕ)g_μν,
                      d(*F)=0   in ^4,
                      □_g^0ϕ=
                1/4κ^2 e^√(3) ϕ(g^0)^βF_νF_βμ   in ^4,
                
                where T is called the stress energy tensor, * is the Hodge operator with respect to  the metric g^0_μν(y),
                ∇_μϕ=ϕ/ x^μ denotes the partial derivative,
                □_g^0 is the Laplacian 
                 with respect to the Lorentzian metric g^0_μν(y)
                 (i.e. the wave operator),
                and
                F_μν(y) is the exterior derivative of the 1-form
                A^0(y)=A_μ^0(y)dy^μ,
                F=dA^0.
                Physically, if we write F=E(y)∧ dy^0+B(y), then E corresponds
                to the electric field and B the magnetic flux  (see <cit.>). 
                
                
                
                
                
                
                
                
                 
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                In the above, equation  (<ref>) for  F=dA^0 are 
                the 4-dimensional formulation of
                Maxwell's equations, and (<ref>) is 
                a scalar wave equation corresponding  a mass-less scalar field that interacts
                with the A^0 field.
                Equations (<ref>)-(<ref>) are the 4-dimensional Einstein equations in a curved space-time
                with stress-energy tensor T, which corresponds to the stress-energy of the electromagnetic
                field F and scalar field ϕ. Thus Kaluza-Klein theory unified  the 4-dimensional
                Einstein equation and Maxwell's equations. However, as the
                scalar wave equation did not correspond to particles observed in physical
                experiments, the theory was forgotten for a long time due to
                the dawn of quantum mechanics. Later in 1960-1980, it was re-invented in the creation of string theories
                when the manifold S^1 was replaced by a higher dimensional manifolds, see <cit.>.
                However,  the Kaluza-Klein theory is still considered as an interesting
                simple model close to string theory suitable for testing
                ideas.
                
                
                
                To consider the relation of the Kaluza-Klein model to the main text in the paper, let
                us consider 5-dimensional Einstein equations 
                 with some matter model on  manifolds N_=× M_, >0,
                 where ({t}× M_,  g_(·,t)), t∈, are compact 4-dimensional Riemannian manifolds.
                 Let ĝ_=-dt^2+g_(·,t) be the background
                metric on N_ and  assume that the metric  g_(·,t) is independent of the variable
                t. Moreover, assume that we can make small perturbations to the matter fields
                in the domain _+× M_
                that cause the metric to become a small perturbation
                 ĝ_(·,t;h)=-dt^2+g_(·,t;h) of  the metric  g_(·,t),
                 where h>0 is a small parameter related to the amplitude of the perturbation.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                By representing tensors  g_(t,x;h) for all h at appropriate coordinates (the so-called wave
                gauge coordinates), one obtains that
                the tensor g̃(t,x)=_h ĝ_(t,x;h)|_h=0
                satisfies the linearized Einstein equations, that is, a wave equation
                      □̂g̃_jk(t,x)+b_jk^lpq(t,x)∇̂_lg̃_pq(t,x)+
                c_jk^pq(t,x)g̃_pq(t,x)=T̃_jk(t,x)  on M_×,
                
                      g̃_pq(t,x)=0  for t<t_- for some t_-∈,
                see <cit.>, where □̂=□^ĝ
                is the wave operator, ∇̂=∇^ĝ is the 
                covariant derivative
                with respect to the metric ĝ, and T̃ is a source
                term corresponding to the perturbation of the stress-energy tensor.
                We mention that in realistic physical models,  T̃ should satisfy
                a conservation law but we do not discuss this issue here.
                
                Let us now consider a scalar equation analogous to (<ref>)  for 
                a  real-valued function  U_(t,x) on × M_,      □̂U_(t,x)+B_^ν(x)∇_νU_(t,x)+
                C_(x)U_ (t,x)=F_(t,x)  on × M_,
                
                      U_(t,x)=0  for t<0.
                
                We will apply to the solution U_(t,x)  the wave-to-heat transformation 
                
                    ( Tf)(t)=1/2π∫_× e^-ξ^2t+it'ξf(t') dt'dξ
                
                 in the time variable, and denote u_(t,x):=( TU_)(t,x) and
                  f_(t,x):=( TF_)(t,x). 
                
                
                
                Then u_ satisfies the heat equation 
                      (/ t+Δ_g_+B_^ν(x)∇_ν+C_(x) )u(t,x)=f_(t,x)  on M_×_+,
                         u_|_t=0=0,
                
                where Δ_g_ is the 3-dimensional (nonnegative definite) Laplace-Beltrami operator
                on (M_,g_). 
                
                
                Then, if we 
                can control the source term F_ and measure the field U_
                for the wave equation (with a measurement error), we can also produce many sources
                f_ for the heat equation and compute the corresponding fields
                u_. In this paper we have assumed that we are given the values
                of the heat kernel, corresponding to measurements with point sources, at the δ-dense points in the
                subset (δ,δ^-1)× M_ of the space-time with some error.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                Due to the above relation of the heat equation to the wave equation
                and the hyperbolic nature of the linearized Einstein equation, the inverse problem for
                the pointwise heat data can be considered as a (very much) simplified
                version of the question: if the observations of the small perturbations
                of physical fields in the subset ×Ω of an almost stationary, almost collapsed
                universe × M_ can be used to find the metric of 
                × M_ in a stable way.
                 As M_ can be considered as a S^1-fiber bundle π:M_→ M_0 on a 3-manifold M_0, it is interesting
                to ask if the  measurements at the δ-dense subset, where δ is much larger than , 
                can be used to determine e.g. the relative
                volume of the
                almost collapsed fibers π^-1(y), y∈ M_0. Physically,
                this means the question if 
                 the macroscopic measurements be used to find information on
                the possible changes of the parameters of the almost collapsed structures of the universe
                in the space-time. 
                We emphasize that the questions discussed here in this subsection are
                not related to the practical testing of string theory, but more to the philosophical
                question: can the properties of the almost collapsed structures in principle be observed
                using macroscopic observations, or not.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                 
                
                
                
                
                
                 
                
                
                
                § REMARKS ON THE SMOOTHNESS IN CALABI-HARTMAN <CIT.> AND MONTGOMERY-ZIPPIN
                <CIT.>
                
                
                The goal of this section is to prove a generalization of Montgomery-Zippin's Theorem <cit.>, Proposition <ref>, which is used in part I of the paper <cit.>.
                Let N be either the manifold M or G × M, where G is a Lie group
                 of isometries acting on M.
                We consider Zygmund classes C^s_*(N), s>0. To define these spaces,
                we cover N by a finite number of coordinate charts, (U_a, Φ_a) with, e.g.,
                Φ_a(U_a)= Q_2 r, where Q_r is a cube with side r,
                and assume that 
                
                    ⋃_a=1^J (Φ_a^-1(Q_r) ) =N.
                
                Then the definition of the norm in C^s_*(N)  is analogous to
                the definition of these spaces in Euclidean space, cf. <cit.>.
                
                
                Let ŝ∈_+, s̃∈ (0, 1] and s=ŝ +s̃.
                 We say that f: N → is in C^s_*(N), if for
                some δ <r,
                f_C^s_*(N) = f_C(N)+ ∑_a=1^J [
                sup_x ∈ Q_rsup_|h| < δ∑_|β| ≤ŝ1/|h|^s̃ |(Δ ^2_h ^βf_a)(x) | ] < ∞,
                
                where f_a:=(φ_a  f)∘Φ_a^-1 and φ_a∈ C^∞_0(U_a) are functions
                for which ∑_a=1^J φ_a(x)=1, and
                
                    (Δ ^2_h f_a)(x)= f_a(x+h)+f_a(x-h)-2f_a(x).
                
                If condition (<ref>) is satisfied, it defines the norm of f in 
                C^s_*(N). We will also consider maps F:N→ M and denote F∈ C^s_*(N;M) if the coordinate representation of F in the local
                coordinates are C^s_*-smooth.
                
                
                Note that even though the norm (<ref>) depends on the local coordinates used,
                the smooth partition of unity φ_a, and on δ, the resulting norms are
                equivalent. Also, when s ∉_+, C^s_*(N) coincides with the Hölder spaces C^ŝ, s̃(N).
                
                By <cit.>, the norm involving only terms with the finite differences along the coordinate axis x^j
                of the partial derivatives  along the same coordinate axis, namely
                
                f_C^s_*(N)^(1)=f_C(N)+ ∑_a=1^J  [
                sup_x ∈ Q_rsup_0<ρ < δ∑_j=1^n 1/ρ^s̃ |(Δ ^2_ρ, j ^ŝf_a/( x^j)^ŝ)(x) | ],
                
                is equivalent to (<ref>), where
                
                    (Δ ^2_ρ, j h)(x)= h(x+ρ e_j)+h(x-ρ e_j)-2h(x).
                
                
                Let  (M^i,h^i), i=1,2, be Riemannian manifolds and let B^i⊂ M^i
                be metric balls. We will use estimates presented in  <cit.> for a
                map F:M^1→ M^2 whose restriction to the ball B^1  defines 
                an isometry F: (B^1, h^1) → (B^2, h^2). We note that the constants in
                these estimates  depend only
                on the norms of h^i,   i=1,2, in the appropriate function classes C^s_* in
                some larger balls containing B^i and the radii of these balls.
                Thus, if a Lie group G has isometric actions on manifold M 
                and F=F_g:M→ M is the action of the group element g ∈ G, then in the case
                where M and G are compact, we obtain uniform estimates by covering the manifold and the Lie group
                with finite number of balls. In the case where M is compact but G  is not, we observe
                that the estimates on a finite collection of balls covering M and one ball B⊂ G for which ∪_g∈ G(g B)=G yield uniform estimates 
                 on the space G× M.
                
                Let (M, h) be a compact Riemannian manifold with a metric h ∈ C^s_*(M),s>1, i.e., in suitable local coordinates the elements h_jk(x) of the metric tensor h are in C^s_*(M).
                Let G be a Lie  group of transformations acting on M as isometries, i.e., for
                any g ∈ G  the action of g, denoted by F_g:M → M, is an isometry.
                Let us denote
                
                    F: G × M → M,    F(g, x)= F_g(x).
                
                We will use local coordinates (x^j)_j=1^n=(x^1,…,x^n) of M and 
                (g^)_=1^p=(g^1,…,g^p)  of G. Below, we will use  Latin indices i,k,l for 
                coordinates on M and Greek indices ,β,γ on coordinates on G.
                Thus e.g. for a function f:G× M→ we often denote
                  _ f(g,x)= f/ g^(g,x) and _j f(g,x)= f/ x^j(g,x).
                
                
                Next, let us assume that
                 h∈ C^s_*(M),s>1. Our aim is to prove that
                 F:G× M→ M is in C^s+1_*(G× M;M), that is, in local coordinates
                 the components F^m(g,x) of F(g,x) are in  C^s+1_*(G× M).
                 Let us start with the case when 1 <s ≤ 2.
                 Let g∈ G be fixed for the moment and denote F=F_g. Then F ∈ C^2(M;M) due to <cit.>.
                Using local coordinates of M in sufficiently small balls B_1 and B_2 satisfying  F(B_1)⊂ B_2, we have by 
                 <cit.>,
                Γ^(1), p_i j(x)  F^m/ x^p -Γ^(2), m_p q(F(x))  F^p/ x^i  F^q/ x^j= ^2 F^m/ x^i  x^j,
                
                where Γ^(1), p_i j(x) and  Γ^(2), p_i j(x) are the Christoffel symbols of the metric h in 
                B_1 and B_2, correspondingly.
                As h ∈ C^s_* and F ∈ C^2, we see easily that all terms in the left side of formula (<ref>), except maybe  the term Γ^(2), m_p q(F(x)), are in
                C^s-1_*. Next we consider
                this term and will show that since 
                Γ^(2), m_p q∈ C^s-1_* and F ∈ C^2, their
                composition satisfies 
                Γ^(2), m_p q∘ F ∈ C^s-1_*(B_1).
                
                To show this, observe that since F ∈ C^2, we have
                |F(x+y/2) -1/2 F(x)-1/2 F(y)| ≤ C|x-y|^2.
                
                Let (s-1)/2< t <s-1. Then Γ^(2), m_p q∈ C^t(B_2), and for x,y∈ B_1 we have
                
                    |Γ^(2), m_p q(F(x+y/2))- 
                    Γ^(2), m_p q( 1/2 F(x)+1/2 F(y)) | ≤ c |x-y|^2t≤ c'|x-y|^s-1.
                
                Moreover, we see that
                      |Γ^(2), m_p q(F(x))+Γ^(2), m_p q(F(y))-2Γ^(2), m_p q(F(x+y/2))| ≤
                      |Γ^(2), m_p q(F(x))+Γ^(2), m_p q(F(y))-
                2Γ^(2), m_p q( 1/2 F(x)+1/2 F(y))
                |+c'|x-y|^s-1
                      ≤ c|F(x)-F(y)|^s-1+c'|x-y|^s-1≤ C|x-y|^s-1.
                
                which proves the claim (<ref>).
                
                By <cit.>, the space C^s-1_*(M)  with s>1 is an algebra, that is, the pointwise
                multiplication satisfies C^s-1_*(M) · C^s-1_*(M) ⊂ C^s-1_*(M). Thus
                it follows from (<ref>) that
                ^2 F^m/ x^i  x^j∈ C^s-1_*(B_1). This shows that
                F_g ∈ C^s+1_*(M;M) for g∈ G.
                
                Differentiating further (<ref>)
                with respect to variables x^j and
                repeating the above considerations, we see that if h∈ C^s_*(M) with s ∈ (2, 3] then F_g ∈ C^s+1_*(M,M). Iterating this construction, we obtain the following result.
                
                
                
                
                Let (M, h), h ∈ C^s_*(M), s>1, be a compact 
                Riemannian manifold.
                Let G be a Lie  group of transformations acting on M as isometries, i.e., for
                any g ∈ G,  the corresponding action F_g:M → M is an isometry.
                Then, 
                for each g∈ G  the map F_g is in C^s+1_*(M;M), and the norm
                 of F_g  in C^s+1_*(M;M) is
                uniformly bounded with respect to g ∈ G.
                
                
                
                We now turn to <cit.>. The corresponding result in <cit.> which we need
                is as follows (see Theorem on p. 212, sec. 5.2, <cit.>).
                
                 (Montgomery-Zippin)
                Let M be a differentiable manifold of class C^k, k ∈_+.
                Let G be a Lie group of transformations acting on M so that F_g(·) ∈ C^k(M;M) uniformly in g∈ G, where F_g is the action of g∈ G. 
                Define F:G× M→ M by F(g,x)=F_g(x).
                Then, 
                in   the local real-analytic coordinates 
                of G and the C^k-smooth coordinates of M,  we have
                
                    F ∈ C^k(G × M;M).
                
                Note that as G is a Lie group it has an analytic structure.
                
                
                Our goal is to extend Theorem <ref> to  spaces C^s_*.
                
                
                
                
                 (Generalization of  Montgomery-Zippin's theorem)
                
                 Let s>1 and (M, h) be a compact Riemannian manifold
                with  the C^s_*-smooth coordinates.
                Let G be a Lie  group of transformations acting on M
                
                
                
                
                
                so that F_g∈ C^s_*(M;M) uniformly in g∈ G.
                Then, in  the local real-analytic coordinates 
                of G and the C^s_*-smooth coordinates of M, 
                
                    F ∈ C^s_*(G × M;M).
                
                Let 1< s'< s”<s.
                Let us show that the map
                ℱ:G → C^s'_*(M;M),   g ↦ F_g(·)
                
                is continuous.  Essentially, this follows from the facts that s'<s and our assumption that
                 F_g∈ C^s_*(M;M) uniformly in g∈ G, and that
                F(g, x) is uniformly continuous with respect to (g, x)∈ G× M. 
                Let us explain this in details.
                First, we consider the case when 1<s ≤ 2 and  a scalar, uniformly continuous function f:G× M→.
                Denoting f_g(x)=f(g, x) and assuming
                that g↦ f_g is continuous map G→ C^s'_*(M),
                we have, for h,ρ∈, |h|<1, |ρ|≤ 1,
                      |(f_g-f_g')(x+h)/ x^i - (f_g-f_g')(x)/ x^i
                      -1/ρ((f_g-f_g')(x+h+ρ e_i)-(f_g-f_g')(x+h)
                          
                +(f_g-f_g')(x+ρ e_i)-(f_g-f_g')(x)  )| 
                ≤
                
                C(s”) ρ^s”,
                
                where C(s”) is uniform with respect to g, g'∈ G.
                
                
                
                On the other hand,
                1/|h|^s'| (f_g-f_g')/ x^i(x+h)-  (f_g-f_g')/ x^i(x) |
                ≤ C(s”) |h|^s”-s'.
                
                Thus, for any >0,  we can find h_0>0 such that the left hand side of (<ref>) is less than /2 for |h| <h_0.
                Moreover, we can find ρ_0>0 such that the right hand side of (<ref>) is less than /2 for
                |ρ| ≤ρ_0. As f:G× M→ is continuous and, therefore uniformly continuous,
                any g∈ G has a neighborhood V_(g) in  G, such that for g'∈ V_(g) and x∈ M,
                
                    |(f_g-f_g')(x)| < 1/8ρ_0  h_0^s'.
                
                Combining the above considerations we see that if g' ∈ V_(g), then
                
                    f_g-f_g'_C^s'_*(M) < .
                
                Applying the above for the coordinate representation of F:G× M→ M, we obtain (<ref>) in the case when 1<s'<s≤ 2. Analyzing the higher derivatives similarly,
                we obtain (<ref>) for all s>s'>1.
                
                
                
                
                Next, let g^, =1, 2,…, p, be real-analytic coordinates on
                G near the identity element id, for which the coordinates
                of the identity element are (0,0,…,0).
                Let e_=_g_|_g=id, and 
                t e_,t ∈, denote the elements of the one-parameter subgroup in G generated by e_.
                 Let
                 x^j,j=1,2,…,n, be the local
                coordinates in an open set B⊂ M and  x_0∈ B. Near (id,x_0)∈ G× M, we represent F in these coordinates
                as F(g,x)=(F^m(g,x))_m=1^n.
                As noted before, we will use  Latin indices i,k,l for 
                coordinates on M and Greek indices ,β,γ on coordinates on G.
                In particular, we will use this to indicate derivatives with respect to
                g^ and x^j. 
                
                
                
                Consider now the identity given in <cit.>,
                 
                F^m(ρ e_, x)-F^m(0, x)=∑_j=1^n (∫_0^1 _j F^m(tρ e_, x) dt )
                _ F^j(0, x),
                
                where ρ>0  is sufficiently small, and as noted before, 
                
                    _j F^m(g,x)= F^m/ x^j(g,x)  and  _ F^j(g,x)= F^j/ g^(g,x).
                 
                Since  
                
                    _j F^m(0, x)=δ^m_j,   x ∈ M,
                
                it follows from the uniform continuity of ∇_x F(g,x) with respect to (g,x) that the matrix
                [∫_0^1 _j F^m(tρ e_, x)  dt ]_j, m=1^n
                
                is invertible for sufficiently small ρ (note that ρ>0  can be chosen uniformly  with respect to x∈ M). 
                Denoting the inverse matrix (<ref>) by
                ϕ^m_j(ρ e_, x), we obtain from (<ref>) the identity
                _ F^m(0, x)= ϕ^m_j(ρ e_, x)  [ F^j(ρ e_, x)-F^j(0, x) ],
                
                when ρ >0 is sufficiently small. In the following, we can choose ρ_0>0
                so that (<ref>) is valid for all 0<ρ<ρ_0  and  x ∈ M. 
                Using our assumption that
                 F_g∈ C^s_*(M;M), uniformly in g∈ G, we see that the matrix
                (<ref>) is in C^s-1_*(M) and thus its inverse matrix satisfies ϕ^m_j(ρ_0 e_, ·) ∈ C^s-1_*(M).
                As F_g ∈ C^s_*(M;M) uniformly with respect to g∈ G, formula (<ref>) implies
                that
                _ F^m(0, x)∈ C^s-1_*(M).
                
                
                Our next goal is to show that 
                
                    _ F^m(g, x)= F^m(g, x)/ g^∈ C^s-1_*(G),
                
                uniformly with respect to x ∈ M.
                To this end, let g_0∈ G and x∈ M be fixed, and consider an element g∈ G which varies 
                in a neighborhood of g_0^-1.
                Let us consider function g̃(g):=g_0  g^-1, g̃:G→ G and
                ỹ(g)=F_g(x), ỹ:G→ M.
                The group G has real-analytic local coordinates (g^)_=1^p=
                (g^1,…,g^p) near g_0. 
                As the group operation (g',g”)↦ g'  (g”)^-1 is real-analytic, the local coordinates of g̃(g),
                denoted by g̃^(g^1,… ,g^p), are also real-analytic. 
                Then,
                
                    F(g̃(g), ỹ(g))= F(g_0  g^-1, F_g(x))=F(g_0  g^-1  g, x)=F(g_0, x).
                
                Thus, F(g̃(g), ỹ(g)) is independent of g, and using the
                chain rule, we obtain
                
                0   =   / g^(F(g̃(g), ỹ(g)))
                   =   _β F^m(g̃(g), ỹ(g)) g̃^β(g)/ g^ 
                +
                _j F^m(g̃(g), ỹ(g)) ỹ^j(g)/ g^.
                
                
                
                
                
                
                Let us next use that fact that for all z∈ M  we have  _j F^m (g,z)|_g=id=δ^m_j. 
                As derivatives of the map (g,z)↦_j F^m (g, z)
                 is continuous, the matrix (_j F^m (g, z))_j,m=1^n is invertible for (g,z) near (id,z).
                 Let us  denote this inverse matrix by
                by Φ^m_j(g,z). Composing this function with g̃(g) and ỹ(g), we see that  the matrix
                (_j F^m(g̃(g), y(g)))_j,m=1^n has an inverse
                Φ^m_j(g̃(g), y(g)) when g is sufficiently near to g_0.
                Then we obtain from (<ref>),
                _ F^k (g, x)   =   _ỹ^m(g, x)
                   =   -Φ^m_k(g̃(g),F(g, x)) ·_β F^k(g̃(g), F(g, x)) ·_g̃^β(g).
                
                Now we evaluate (<ref>) at g=g_0. As Φ^m_k(g̃(g_0),F(g_0, x))=δ^m_k, and in the local
                coordinates of G we have
                g̃(g_0)=0, we obtain 
                _ F^m (g_0, x)
                =-
                _β F^m (0, F(g_0, x)) ·_g̃^β(g_0).
                
                Note that g_0 and x above were arbitrary, even though we considered those as fixed parameters.
                Next we change our point of view and will consider those as variables. For this
                end, we use variable g∈ G instead of g_0 so that (<ref>) reads as
                _ F^m (g, x)
                =-
                _β F^m (0, F(g, x)) ·_g̃^β(g).
                
                Here g↦_g̃^β(g)  is real-analytic and by formula (<ref>), 
                the map z↦_β F^m (0, z) is in C^s-1_*(M).
                
                
                
                Let us next show that 
                _β F^m (0, F(g, x)) ∈ C^s-1_*(G × M).
                
                Once this is shown, the fact that
                _g̃^β(g)  is real-analytic  and the equation
                (<ref>) imply that  _ F^m(g, x) ∈ C^s-1_*(G × M).
                
                
                
                To show (<ref>),
                we start by considering the case when 1<s < 2. Then, using (<ref>)
                and that fact that  F ∈ C^1(G × M) by Theorem <ref>, we obtain in local coordinates
                      |_β F^m(0, F(g+h,x))- _β F^m(0, F(g,x)) | 
                
                   ≤    C|F(g+h, x)-F(g, x)|^s-1≤ C' |h|^s-1,
                
                where C, C'>0 are uniform with respect to (g, x)∈ G× M. 
                Hence (<ref>) is valid, i.e., _ F^m(g, x) ∈ C^s-1_*(G × M).
                Combining this with (<ref>) we see that
                for any x∈ M,  the function g↦_ F^m(g, x)  
                is in C^s-1_*(G) and its norm in  C^s-1_*(G)  is
                bounded by a constant which is independent of x∈ M.
                
                
                By assumption,
                functions x↦_j F^m(g, x)
                are in C^s-1_*(M), uniformly in g∈ G. 
                Thus, by  (<ref>), i.e., <cit.>, the map 
                F:G× M→ M  is C^s_*-smooth with respect to
                both g and x, that is,
                
                F(g,x) ∈ C^s_*(G × M).
                
                
                
                
                Next, let us consider the case when  s=2.
                To apply (<ref>), we need to estimate 
                
                    T_β^k :=_β F^k(0, F(g, x))
                
                and in particular its finite difference
                
                    Δ^2_h (T_β^k)=_β F^k(0, F(g+h, x))
                     +_β F^k (0, F(g-h, x)) -
                     2_β F^k (0, F(g, x)).
                
                 Using (<ref>), we observe that for h=(h^)_=1^p,
                
                    |F^k(g+h, x)- F^k(g, x)- _ F^k(g, x) · h^| ≤ C_s |h|^s',   for any  
                    s' <2.
                
                On the other hand, by (<ref>), we have
                
                    _β F^k (0, F(g, x)) ∈ C^t(G × M), 
                      for any  
                    t <1.
                
                Let us combine these two  observations in the 
                case when s'=3/2 and t=2/3. Then
                Δ^2_h (T_β^k )
                   =   _β F^k(0, F(g, x)+ (_ F^k)(g, x)  h^+O(|h|^s'))
                      + _β F^k (0, F(g, x)- (_ F^k)(g, x)   h^+O(|h|^s')) 
                 
                      -
                 2_β F^k (0, F(g, x)) 
                 
                   =   O(  (|h|^s')^ t)= O(|h|).
                
                Applying (<ref>) with s=2 to estimate the finite
                differences of T_β^k in the  x^j directions
                and  (<ref>)
                to estimate the finite
                differences of T_β^k  in the  g^ directions
                 in the formula (<ref>) for the Zygmund norm,
                 we obtain  T_β^k=_β F^k(g, x) ∈ C^1_*(G × M). 
                Moreover, by our assumption,
                functions x↦_j F^m(g, x)
                are in C^1_*(M) uniformly in g∈ G. 
                Thus, by applying (<ref>) again,  we see that F(g, x) ∈ C^2_*(G × M).  
                
                
                Next, we consider the case when 2<s ≤ 3.  By differentiating formula (<ref>) 
                with respect to g^β, we obtain
                __β F^k(g, x)   =   -/ g^β( _γ F^k (0, F(g, x))  _g̃^γ(g))
                
                    =   (_γ_j F^k (0, F(g, x)) _β F^j(g, x) _g̃^γ (g)
                
                      +_γ F^k (0, F(g, x))  ^2 g̃^γ/ g^ g^β(g).
                
                We see easily that all terms in the right hand side of 
                equation (<ref>), except maybe the term    _γ_j F^k (0, F(g, x)),
                are in C^1_*(G × M). To analyze this term, we
                observe that by (<ref>) we have _γ_j F^k (0, y) ∈ C^s-2_*(M). 
                Thus considerations similar to those leading to 
                the inequality (<ref>) show that
                
                    _γ_j F^k (0, F(g, x)) ∈ C^s-2_*(G × M).
                
                Hence by (<ref>),
                __β F^k(g, x)∈ C^s-2_*(G),
                
                uniformly with respect to x∈ M.
                By assumption,
                F_g(·) ∈ C^s_*(M;M) uniformly in g∈ G,
                so that
                _i_j F^k(g, x)∈ C^s-2_*(M),
                
                uniformly with respect to g∈ G.
                
                
                Applying (<ref>) with  (<ref>) and (<ref>),
                we see that F ∈ C^s_*(G × M).
                Iterating the procedures above we obtain (<ref>) for all s>1.
                
                
                Combining Lemma <ref> and Proposition <ref>, we obtain the following result.
                
                
                Let M be a compact, C^s+1_*-smooth differentiable manifold, s>1, and h ∈ C^s_*(M) be a Riemannian
                metric on M.
                Let G be a Lie group  acting on (M,h) for which the actions F_g: M→ M of the elements g∈ G are isometries. Define F:G× M→ M by F(g,x)=F_g(x).
                Then F ∈ C^s+1_*(G × M;M).
                
                
                § DEFINITION OF ORBIFOLDS
                
                
                Orbifolds were first introduced under the name "V-manifolds" in <cit.> and later studied in <cit.> under its current name.
                Here we review the definition in <cit.>.
                
                
                An (d-dimensional) orbifold chart on a topological space X is the triple (Ũ,G,φ) consisting of a connected open set Ũ⊂ℝ^d, a finite group G of diffeomorphisms of Ũ, and a continuous map φ from Ũ to an open set U⊂ X inducing a homeomorphism from the quotient space Ũ/G to U.
                
                Let (Ũ,G,φ) and (Ũ',G',φ') be two orbifold charts for U and U' with U⊂ U'.
                By an smooth embedding ι:(Ũ,G,φ)→ (Ũ',G',φ') we mean a smooth embedding ι:Ũ→Ũ' such that φ=φ'∘ι.
                
                
                An orbifold of dimension d consists of a Hausdorff topogical space X and a family 𝒜 of d-dimensional orbifold charts for open subsets in X satisfying the following conditions.
                
                
                (1) Every point x∈ X is contained in at least one orbifold chart in 𝒜 (i.e., there exists an open neighborhood U containing x such that there is an orbifold chart for U in 𝒜).
                If x is contained in two orbifold charts in 𝒜 for U_1 and U_2, then there exists an open set U_3⊂ U_1∩ U_2 containing x such that there is an orbifold chart for U_3 in 𝒜.
                
                
                (2) If (Ũ,G,φ) and (Ũ',G',φ') for U and U' are two orbifold charts in 𝒜 such that U⊂ U', then there exists a smooth embedding from (Ũ,G,φ) to (Ũ',G',φ').
                
                
                A point x of an orbifold is called singular, if there exists an orbifold chart (Ũ,G,φ) for some open set U containing x such that the isotropy group of φ^-1(x)∈Ũ in G is nontrivial. Points that are not singular are called regular, and the set of regular points has the structure of an open manifold.
                A closed manifold is an orbifold with the group G in orbifold charts being the trivial group.
                The quotient space ℝ^2/ℤ_m, m≥ 2, is a cone with cone angle 2π/m, where ℤ_m acts on ℝ^2 by a rotation of angle 2π/m around the origin.
                A manifold with boundary can be given an orbifold structure: every point on the boundary has a neighborhood homeomorphic to an open set in the half-plane ℝ^d/ℤ_2, where ℤ_2 acts on ℝ^d by reflection along the hyperplane that models the boundary.
                If M is a manifold and G is a group acting properly discontinuously on M, then the quotient M/G has an orbifold structure.
                An orbifold is said to be good if it is the quotient of a manifold by an action of a discrete group.
                
                
                A Riemannian metric on an orbifold is an assignment to each orbifold chart (Ũ,G,φ) of a G-invariant Riemannian metric g_Ũ on Ũ satisfying that each embedding in Definition <ref>(2) is an isometry.
                A Riemannian orbifold is a smooth orbifold equipped with a Riemannian metric.
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                
                99
                
                A90 M. Anderson, Convergence and rigidity of manifolds under Ricci curvature bounds, Invent. Math.  102 (1990), 429–445.
                
                AKKLT M. Anderson, A. Katsuda, Y. Kurylev, M. Lassas, M. Taylor, Boundary regularity for the Ricci equation, geometric convergence, and Gel'fand's inverse boundary problem, Invent. Math. 158 (2004), 261–321.
                
                book3
                T. Appelquist, A. Chodos, P. Freund, Modern Kaluza-Klein Theories, Frontiers in Physics, Addison Wesley, 1987.
                
                AG
                V. Arnaiz, C. Guillarmou, Stability estimates in inverse problems for the Schrödinger and wave equations with trapping, Rev. Mat. Iberoam. 39 (2023), 495–538.
                
                Ba
                S. Bando, Real analyticity of solutions of Hamilton's equation, Math. Z. 195 (1987), 93–97.
                
                BK M. Belishev, Y. Kurylev, To the reconstruction of a Riemannian manifold via its spectral data (BC-method), Comm. PDE. 17 (1992), 767–804.
                
                
                
                BD
                M. Bellassoued, D. Dos Santos Ferreira, Stability estimates for the anisotropic wave equation from the Dirichlet-to-Neumann map, Inverse Probl. Imag. 5 (2011), 745–773.
                
                BMR
                J. Bemelmans, Min-Oo, E. Ruh, Smoothing Riemannian metrics, Math. Z. 188 (1984), 69–74.
                
                
                Bhagavan 
                N. V. Bhagavan, C.-E. Ha,
                
                
                Essentials of Medical Biochemistry,  Academic Press, 2011.
                BKL1 R. Bosi, Y. Kurylev, M. Lassas, Stability of the unique continuation for the wave operator via Tataru inequality: the local case, J. Anal. Math. 134 (2018), 157–199.
                
                BKL2 R. Bosi, Y. Kurylev, M. Lassas, Stability of the unique continuation for the wave operator via Tataru inequality and applications, J. Diff. Eq. 260 (2016), 6451–6492.
                 
                BKL3
                R. Bosi, Y. Kurylev, M. Lassas, Reconstruction and stability in Gel'fand's inverse interior spectral problem, Anal. PDE. 15 (2022), 273–326.
                
                BBI D. Burago, Y. Burago, S. Ivanov, A course in metric geometry, Graduate Studies in Math. 33, AMS, 2001.
                
                BILL
                D. Burago, S. Ivanov, M. Lassas, J. Lu, Quantitative stability of Gel'fand's inverse boundary problem, to appear in Anal. PDE., arXiv:2012.04435v4.
                
                CH E. Calabi, P. Hartman, On the smoothness of isometries, Duke Math. J.  37 (1970),  741–750.
                
                C2
                J. Cheeger, Finiteness theorems for Riemannian manifolds, Amer. J. Math. 92 (1970), 61–74.
                
                
                CFG
                J. Cheeger, K. Fukaya, M. Gromov, Nilpotent structures and invariant metrics on collapsed manifolds, J. Amer. Math. Soc. 5 (1992), 327–372.
                
                
                ChBook
                Y. Choquet-Bruhat, General relativity and the Einstein equations, Oxford Mathematical Monographs, Oxford University Press, 2009.
                
                
                Coifman2
                R. Coifman, S. Lafon, Diffusion maps, Appl. Comput. Harmon. Anal. 21 (2006), no. 1, 5–30.
                
                Coifman1 
                R. Coifman, S. Lafon, A. Lee, M. Maggioni, B. Nadler, F. Warner, S. Zucker, Geometric diffusions as a tool for harmonic analysis and structure definition of data: Diffusion maps, PNAS 102 (2005), 7426–7431. 
                
                
                DLLO
                M. de Hoop, M. Lassas, J. Lu, L. Oksanen, Quantitative unique continuation for the elasticity system with application to the kinematic inverse rupture problem, Comm. PDE. 48 (2023), 286–314.
                
                DCGW E. Dryden, C. Gordon, S. Greenwald, D. Webb, Asymptotic 
                expansion of the heat kernel for orbifolds, Michigan Math. J.  56 (2008), 205–238.
                
                
                
                COLT
                C. Fefferman, S. Ivanov, Y. Kurylev, M. Lassas, H. Naranayan, Fitting a putative manifold to noisy data,
                Proceedings of the 31st Conference On Learning Theory (COLT),
                Proceedings of Machine Learning Research 75 (2018), 1–32. 
                
                FIKLN C. Fefferman, S. Ivanov, Y. Kurylev, M. Lassas, H. Naranayan,
                 Reconstruction and interpolation of manifolds I: The geometric Whitney problem, Found. Comp. Math.
                20 (2020), 1035–1133.
                
                
                FILH0
                C. Fefferman, S. Ivanov, M. Lassas, H. Naranayan,
                Reconstruction of a Riemannian manifold from noisy intrinsic distances, SIAM J. Math. Data Sci. 2 (2020), no. 3, 770–808.
                
                
                
                FILLN C. Fefferman, S. Ivanov, M. Lassas, J. Lu, H. Naranayan,
                 Reconstruction and interpolation of manifolds II: Inverse problems for Riemannian manifolds with partial distance data, to appear in Amer. J. Math., arXiv:2111.14528.
                
                
                FMN
                C. Fefferman, S. Mitter,  N. Narayanan, Testing the manifold hypothesis, J. Amer. Math. Soc. 29 (2016), no. 4, 983–1049.
                
                
                Cryo4  
                J. Frank, A. Ourmazd, Continuous changes in structure mapped by manifold embedding of single-particle data in cryo-EM, Methods 100 (2016), 61–67.
                
                
                
                book2
                T. Frankel, The geometry of physics: an introduction, Cambridge University Press, Cambridge, 1997.
                
                
                F86
                K. Fukaya, On a compactification of the set of Riemannian manifolds with bounded curvature and diameters, Curvature and Topology of Riemannian manifolds, Volume 1201, 1986, pp. 89–107.
                
                F87
                K. Fukaya, Collapsing of Riemannian manifolds and eigenvalues of Laplace operator, Invent. Math. 87 (1987), 517–547.
                
                
                F88
                K. Fukaya, A boundary of the set of the Riemannian manifolds with bounded curvatures and diameters, J. Differ. Geom. 28 (1988), 1–21.
                
                
                F90
                K. Fukaya, Hausdorff convergence of Riemannian manifolds and its applications, Recent Topics in Differential and Analytic Geometry, Adv. Stud. Pure Math. 18 (1990), 143–238.
                
                
                FY:fundgp
                K. Fukaya, T. Yamaguchi, The fundamental groups of almost nonnegatively curved  manifolds, Ann. of Math.  136 (1992), 253–333.
                
                
                G54
                I. Gel’fand, Some aspects of functional analysis and algebra, Proc. ICM. 1 (1954), 253–277.
                
                GW:lipschitz
                R. Greene, H. Wu, Lipschitz convergence of Riemannian manifolds, Pacific J. Math.  131  (1988), 119–141.
                
                
                GLP
                M. Gromov, Structures métriques pour les variétés riemanniennes, Edited by J. Lafontaine and P. Pansu, Textes Mathématiques  1, CEDIC, 1981.
                
                
                HH
                E. Hebey, M. Herzlich, Harmonic coordinates, harmonic radius and convergence of Riemannian manifolds, Rend. Mat. Appl. 17 (1997), 569–605.
                
                
                HV E. Hebey, M. Vaugon, The best constant problem in the Sobolev embedding theorem for complete Riemannian manifolds, Duke Math. J. 79 (1995), 235–279.
                
                HK
                E. Heintze, H. Karcher, A general comparison theorem with applications
                to volume estimate for submanifolds, Ann. Sci. Ecole
                Norm. Sup. series, t.11 (1978), 451–470.
                
                Hor3 L. Hörmander, Lectures on nonlinear hyperbolic differential equations, Mathem. et Applic. 26, Springer-Verlag, 1997.
                
                Jones-Marron 
                M. Jones, J. Marron, S. Sheather, A brief survey of
                bandwidth selection for density estimation, Journal of the American Statistical Association 91 (1996), no. 433, 401–407.
                
                
                
                JK
                J. Jost, H. Karcher, Geometrische methoden zur Gewinnung von a-priori- Schrander für harmonitsche Abbildungen, Manuscripta
                Math. 40 (1982), 27–77.
                
                K93
                A. Kasue, Measured Hausdorff convergence of Riemannian manifolds with Laplace operators, Osaka J. Math. 30 (1993), 613–651.
                
                KKL A. Katchalov, Y. Kurylev, M. Lassas, Inverse boundary spectral problems, Monographs and Surveys in Pure and Applied Mathematics 123, Chapman Hall/CRC-press, 2001.
                
                KrKL
                K. Krupchyk, Y. Kurylev, M. Lassas, Inverse spectral problems on a closed manifold. J. Math. Pures Appl.  90 (2008), 42–59.
                
                KLLY Y. Kurylev, M. Lassas, J. Lu, T. Yamaguchi, Inverse spectral problems for collapsing manifolds I: uniqueness and stability, arXiv:1209.5875v3.
                
                
                
                
                Meila  M. Meila, H. Zhang, Manifold Learning: What, How, and Why, to appear in
                Annual Review of Statistics and Its Application, Volume 11, March 2024.
                
                
                MZ D. Montgomery, L. Zippin, Topological transformation groups, R.E.Kruger Publ., New York, 1974.
                
                
                Pets:conv
                S. Peters, Convergence of Riemannian manifolds, Compositio Math.   62  (1987),  no. 1, 3–16. 
                
                
                PP P. Petersen, Riemannian geometry, 2nd edition, Springer, 2006.
                
                Pet
                A. Petrunin, Parallel transportation for Alexandrov space with
                curvature bounded below, Geom. Funct. Anal.  8  (1998),  no. 1, 123–148.
                
                
                Satake
                I. Satake, The Gauss-Bonnet theorem for V-manifolds, J. Math. Soc. Japan 9 (1957), 464–492.
                
                
                Singer 
                A. Singer, Mathematics for Cryo-Electron Microscopy,
                Proceedings of ICM 2018, 3995–4014, 2019.
                
                
                
                SU P. Stefanov, G. Uhlmann, Stability estimates for the hyperbolic Dirichlet-to-Neumann map in anisotropic media, J. Funct. Anal. 154 (1998), 330–357.
                
                SU2
                P. Stefanov, G. Uhlmann, Stable determination of generic simple metrics from the hyperbolic Dirichlet-to-Neumann map, Int. Math. Res. Not. 17 (2005), 1047–1061.
                
                
                book1 A. Strominger, Kaluza-Klein compatifications, supersymmetry, and
                Calabi-Yau spaces, Quantum fields and strings: A course for mathematicians, Vol 1 and 2, by Deligne P. et al (eds.), AMS, Providence, RI, 1999.
                
                ISOMAP  
                J. Tenenbaum, V. de Silva, J. Langford, A global geometric framework for nonlinear dimensionality reduction, Science 290 (2000), 2319–2323.
                
                
                Thurston W. Thurston, The geometry and topology of three-manifolds, Lecture notes, Princeton University, 1980.
                
                
                Tri H. Triebel, Interpolation theory, function spaces, differential operators, second edition, Johann Ambrosius Barth, Heidelberg, 1995.
                
                Wald R. Wald, General relativity, University of Chicago press, 1984.
                
                Wang-Marron
                X. Wang, J. Marron, A scale-based approach to finding effective dimensionality in manifold learning,
                Electron. J. Stat. 2 (2008), 127–148.
