Let $M,N$ be smooth Riemannian manifolds with boundary (In particular, we assume the boundaries are smooth).
Suppose we have a map $\phi:M \to N$ which satisfies the following properties:
$$(1) \, \, \phi:M \to N \, \, \text{is a bijection}$$
$$ (2) \, \, \phi(\operatorname{int}M)=\operatorname{int}N,\phi(\partial M)=\partial N $$
$$ (3) \, \, \phi:M \to N \, \,\text{is a metric isometry}$$
By the Myers-steenrod theorem, applied to $\phi|_{\operatorname{int}M} :\operatorname{int} M \to \operatorname{int}N $, $\phi$ is a diffeomorphism between $\operatorname{int} M , \operatorname{int}N$.
Applying the theorem for $\phi|_{\partial M}:{\partial M} \to {\partial N}$ (the boundaries are manifolds without boundary), we also get that $\phi|_{\partial M}:{\partial M} \to {\partial N}$ is smooth (it's a diffeomorphism).
Question: Is $\phi$ necessarily smooth as a map $M \to N$?
When looking at the proof of Myers-steenrod theorem (for example here), part of the problem seems to be that initial conditions do not determine a unique geodesic, if the starting point is on the boundary.
The basic idea of the proof is to express $\phi$ in exponential coordinates, and then show that this exponential representation is linear, hence smooth. However, constructing this representation relies on the uniqueness of geodesics.
I suspect there might be a counter example where singularity can occur at the boundary, but I could not fine one.