La terrible rareza
$\def\R\{\mathbb{R}}\def\pd#1#2{\frac{\partial #1}{\partial #2}}$ Para ver por qué la fórmula dada tiene sentido, creo que es mejor empezar con el caso euclidiano para tener una idea de lo que está pasando. (Se podría hacer exactamente lo mismo con el espacio de Minkowski, pero como aquí voy a la intuición me voy a ceñir a lo más conocido).
Por lo tanto, suponemos por ahora que $U = M = \R^n$ con la métrica estándar, por lo que podemos escribir explícitamente la función de distancia al cuadrado: es el escalar $\sigma : \R^n \times \R^n \to \R$ definido por $$\sigma(x,y) = \frac12 |x-y|^2 = \frac12 \sum_i (x^i - y^i)^2.$$
Como tenemos esta fórmula explícita, sabemos inmediatamente que $\sigma$ es diferenciable y podemos escribir sus derivadas parciales $$\pd{\sigma}{x^i} = x^i-y^i,\;\pd{\sigma}{y^i}=y^i-x^i.$$
Para conectar esto con la fórmula general que das, tenemos que cambiar nuestra perspectiva y pensar en términos de derivadas direccionales en su lugar. Dado que $\sigma$ se define en $M \times M$ el espacio de direcciones que podemos diferenciar en un punto $(x,y) \in M\times M$ es el espacio tangente $$T_{(x,y)}(M\times M) = T_x M \times T_y M,$$ que para el caso euclidiano es de nuevo simplemente $\R^n \times \R^n.$ Si tomamos $(v,w)$ en este espacio (por lo que estamos variando simultáneamente $x$ en la dirección $v$ y $y$ en la dirección $w$ ) entonces la derivada direccional correspondiente es (a partir de la fórmula habitual en términos de derivadas parciales) $$D_{(x,y)}\sigma(v,w) = \pd\sigma{x^i}v^i+\pd\sigma{y^i}w^i= (x^i - y^i)(v^i - w^i) = (x-y)\cdot(v-w).$$ Podemos recuperar la fórmula que dio eligiendo $w=0$ (es decir, sólo variando el punto $x$ ) y observando que si requerimos una parametrización de velocidad unitaria en $[a,b]$ el vector tangente unitario es $(x-y)/(b-a)$ . Por intuición, piense en lo que esto implica sobre la función de distancia puntual $d_y(x) = \sqrt{2\sigma(x,y)}$ La regla de la cadena hace que esto se convierta en $$D_xd_y (v) = \frac{x-y}{|x-y|}\cdot v,$$ es decir, la proyección de $v$ en la línea a través de $x$ y $y$ . Si haces un dibujo deberías poder tener alguna intuición geométrica para esto. El hecho de que la fórmula de $D\sigma$ depende de la longitud de $(b-a)$ es sólo una consecuencia de la cuadratura: al moverse $x$ más allá de $y$ todas las pendientes de la parábola $\sigma$ son cada vez más empinadas.
La prueba general
La clave es conocer la fórmula de la primera variación. Esta es la versión para el funcional de energía $E(\gamma) = \frac12 \int |\dot \gamma|^2,$ que funciona para las variedades semirriemannianas generales:
[O'Neill Cap 10 Prop 39] Si $\gamma : (-\epsilon,\epsilon) \times [a,b] \to M$ es una familia suave de curvas $\gamma_s(t) = \gamma(s,t)$ entonces $$\frac{d}{ds}\Big|_{s=0}E(\gamma_s) = g(V,\dot \gamma_0)\Big|_a^b - \int_a^b g(V, \ddot \gamma_0)$$ donde $V = \partial_s \gamma|_{s=0}$ es el campo de variación.
Ahora, por su suposición de convexidad sabemos que las geodésicas que unen pares de puntos son únicas; así que si producimos una geodésica $\gamma_s$ unirse a $x_s$ a $y_s$ podemos concluir inmediatamente $\sigma(x_s,y_s) = E(\gamma_s).$ Así, para determinar $D_{(x,y)}\sigma(v,w)$ podemos simplemente calcular la primera variación de una familia de geodésicas $^1$ $\gamma_s$ tal que $\gamma_0$ es la geodésica desde $x$ a $y$ , $V(a) = v$ y $V(b) = w$ . Resulta que esta información es todo lo que necesitamos para calcular la derivada - utilizando la fórmula que he citado anteriormente obtenemos (observando $\ddot \gamma_0 = 0$ ya que es una geodésica) $$D_{(x,y)}\sigma(v,w) = g(w, \dot \gamma_0(b)) - g(v, \dot \gamma_0(a)).$$ En el caso $w=0$ , eligiendo $v = \partial_\beta$ podemos escribir esto en su notación de índice como $$\nabla_\beta \sigma = -g_{\beta \alpha} \dot \gamma_0^\alpha(a).$$ Hemos perdido la escala $b-a$ aquí, que creo que se debe a la confusión (en general incorrecta) de la longitud al cuadrado $\frac12 \left(\int_a^b |\dot \gamma|\right)^2$ y la energía $\frac12 \int_a^b |\dot \gamma|^2$ . Para las geodésicas (que tienen constante $|\dot \gamma|$ ) estos difieren exactamente por la constante multiplicativa $b-a$ . Como alternativa, puede fijar $[a,b] = [0,1]$ y todo está bien.
- La forma que conozco para demostrar que tal familia geodésica existe requiere el requisito adicional de que $\gamma_0$ no tiene puntos conjugados - lee sobre los campos de Jacobi si quieres conocer los detalles. No estoy seguro de si esta suposición es necesaria para obtener la suavidad de la función de distancia.