Consider a random geometric graph $G_M(n;r)$ whose vertex set consists of $n$ points chosen independently and uniformly from a Riemannian manifold $M$, with edges joining pairs of vertices whose distance in the metric $d_M$ is at most $r$. Let $Δ$ denote the expected average degree of the graph. As is the case for Erdős-Rényi graphs, there is a critical value $Δ_c$, depending only on the dimension of $M$, such that if $Δ> Δ_c$ then $G_M(n;r)$ has a giant component. We show that whenever $Δ> Δ_c$, the giant component of $G_M(n;r)$, equipped with the graph distance, converges to the underlying manifold $M$ in the Gromov-Hausdorff distance after rescaling by an appropriate deterministic factor. Our result holds for $Δ$ depending on $n$ as well, provided $Δ= o(n)$ and $Δ\geq Δ_c + \varepsilon$ for any fixed $\varepsilon > 0$. As a consequence, we show that for any pair of non-isometric compact Riemannian manifolds $M_1$ and $M_2$, there is a polynomial-time algorithm that distinguishes random geometric graphs on $M_1$ and $M_2$ throughout this regime of $Δ.$ In the thermodynamic regime -- i.e.\ when $Δ$ is constant -- our results appear to be new even in the classical cases where $M$ is a sphere or a torus. Our proof makes use of techniques from first-passage percolation which allow us to understand the long-range behavior of the graph distance on small, approximately Euclidean patches of $M$, together with global arguments that glue these local estimates into a global description.