$\DeclareMathOperator\Exp{Exp}$A. Gavrilov has several nice works studying the double exponential map and its properties ([1] and references therein). 

Given a complete Riemannian manifold $(M,g)$ and a point $p \in M$, the double exponential map is defined as $\Exp_p : T_p M \times  T_p M \rightarrow M$ such that, if $v$ is of small norm, 
\begin{align}
    \Exp_p (v,u) := \Exp_{\Exp_p (v)} (u_{\Exp_p (v)}), 
\end{align} 
where $ u_{\Exp_p (v)} $ is the parallel transport of $u \in T_p M$ to $T_{\Exp_p (v)} M$ along the minimizing geodesic connecting $p$ and $\Exp_p (v)$. 

Given $p \in T_p M$, one can define a smooth map $h_p : T_p M \times  T_p M \rightarrow T_p M$ such that, for $v,u \in T_p M$ of small norm, it holds that
\begin{align}
    \Exp_p (h_p (v,u)) = \Exp_p (v,u). 
\end{align}

A. Gavrilov derived a Taylor series approximation for $h_p (v,u)$, which is 
\begin{align}
h_p (v,u) =&  v + u +  \frac{1}{6} R (u,v) v +  \frac{1}{3} R (u,v) u + \frac{1}{12} \nabla_v R (u,v) v \\
            &+ \frac{1}{24} \nabla_u R (u,v) v + \frac{5}{24} \nabla_v R (u,v) u
            +  \frac{1}{12} \nabla_u R (u,v) u + \cdots, \tag{1}\label{1}
\end{align} 
where $R$ is the Riemann curvature tensor. 
This formula is Eq.3 in [1]. 

I strongly suspect that the above result can be thought as a corollary of the celebrated Baker–Campbell–Hausdorff–Dynkin (BCHD) formula. To link these two results, one can perhaps construct "constant vector fields" so that the BCHD formula can be applied. For $p$, within a small neighborhood $p$ (written $U_p$), let $X$ be the vector field such that $X_q = v_q$ for all $q \in U_p$, and let $Y$ be the vector field such that $Y_q = u_q$ for all $q \in U_p$. For vector fields $X$ and $Y$ over $U_p$, we have, if $\exp(Z) = \exp(X) \exp(Y)$, then it holds that
\begin{align}
    Z = X + Y + \frac{1}{2} [X,Y] + \frac{1}{12} [ X ,[X,Y]] + \cdots. \tag{2}\label{2}
\end{align}
It seems \eqref{1} and \eqref{2} don't coincide. 

I probably missed something very basic here, or I'm not using the right vector fields, or perhaps the results of Gavrilov cannot be thought of as a corollary of the BCHD formula. Any pointers? 

**Additionally**:
I've been debating: on one hand, the Lie algebra for the vector fields over $M$ is infinite-dimensional and thus the BCHD series may not converge; on the other hand, the BCHD series should converge for  constant vector fields (of small magnitude).  

**References**:
[1] [A. Gavrilov (2007), The Double Exponential Map and Covariant Derivation, Siberian Mathematical Journal](https://www.academia.edu/7675502/The_double_exponential_map_and_covariant_derivation).