"Second" Lie derivative?
Many of us are familiar with the standard definition of the Lie derivative of some smooth function $\varphi \in \Omega^{0}(X)$ as $$ \frac{d}{dt}(f_{t}^{*}\varphi) \big|_{t=0} =: L_{V}\varphi, $$ where $f_{t}$ is the flow generated by the vector field $V$. My question is given this, can we compute $$ \frac{d^{2}}{dt^{2}}(f_{t}^{*}\varphi) \big|_{t=0} $$ in terms of $L_{V}\varphi$? My naive hunch is that this would simply be $L_{V}(L_{V}\varphi)$, but I'm cautious about this and haven't found a convincing way of proving it. Is this true?
Moreover, I wonder if whatever we compute this "second derivative" to be, should it hold not only for $\varphi \in \Omega^{0}(X)$ but for any $\omega \in \Omega^{k}(X)$, $0 \leq k \leq \mathrm{dim}(X)$? That would be even nicer.
Thanks in advance to anyone who thinks a bit about this and shares any ideas.
We have \begin{align}L_V(L_V\phi):&= \frac{d}{ds}\bigg\vert_{s=0}f^*_s\bigg(\frac{d}{dt}\bigg\vert_{t=0}f^*_t\phi\bigg)\\&=\frac{d}{ds}\bigg\vert_{s=0}\bigg(\frac{d}{dt}\bigg\vert_{t=0}f^*_sf^*_t\phi\bigg)\\ &=\frac{\partial^2}{\partial s\partial t}\bigg\vert_{s=t=0}\bigg(f^*_sf^*_t\phi\bigg)\\ &=\frac{\partial^2}{\partial s\partial t}\bigg\vert_{s=t=0}\bigg(f^*_{s+t}\phi\bigg)\\ &=\frac{1}{2}\frac{d^2}{dt^2}\bigg\vert_{t=0}\bigg(f^*_{2t}\phi\bigg)=\bigg(\frac{d^2}{dt^2}\bigg\vert_{t=0}f^*_{t}\phi\bigg) \end{align}
So we can calculate this second derivative by Cartan's magic formula: $L_V\phi=d\iota_V\phi+\iota_Vd\phi$. Applying this formula twice gives: \begin{align} L_V(L_V\phi)&=d\iota_V(d\iota_V\phi+\iota_Vd\phi)+\iota_Vd(d\iota_V\phi+\iota_Vd\phi)\\ &=d\iota_Vd\iota_V\phi+d(\iota_V\iota_V\phi)+\iota_Vd^2\iota_V\phi+\iota_Vd\iota_Vd\phi\\ &=d\iota_Vd\iota_V\phi+\iota_Vd\iota_Vd\phi \end{align} Since Cartan's magic formula applies to all $k$-forms, we can apply this to any form.
Edit: Since the second to last equality given in my proof is not obvious I will elaborate. (after having deleted my post due to a few errors I felt conscious about, I have decided to rewrite my edit for more clarity and less objectionable notation).
Take $\alpha\in \Omega^k(M)$ and $v_1,v_2,\cdots, v_k\in T_pM$. By definition $$(L_V\alpha)(v_1,v_2,\cdots, v_k):=\frac{d}{dt}\bigg\vert_{t=0}(f^*_t\alpha)(v_1,v_2,\cdots, v_k)=\frac{d}{dt}\bigg\vert_{t=0}\alpha(Tf_{t}\cdot v_1, Tf_t\cdot v_2,\cdots, Tf_t\cdot v_k).$$ We can see that $f^*_t\alpha(v_1,v_2,\cdots, v_k)$ is a function from $(a,b)\to \mathbb{R}$, since $Tf_t\cdot v_1$ gives us a vector at $f_t(p)$ that we can apply $\alpha$ to. We can also apply $f^*_s$ to $f^*_t\alpha$ to get $f^*_{s}f^*_{t}\alpha=f^*_{s+t}\alpha$.
Since $(f^*_{t+s}\alpha)(v_1,v_2,\cdots, v_k):=F(s,t)$ is a smooth function from $(a,b)\times(a,b)\to \mathbb{R}$ we can apply the regular rules of multivariate calculus to it. Specifically the fact that for $F(s,t)=G(s+t)$, we have $\frac{\partial^2 F}{\partial s\partial t}(0,0)=G''(0)$.
This means that for each $v\in T_pM$, we have $$\frac{\partial^2}{\partial s\partial t}\bigg\vert_{s=t=0}(f^*_{s+t}\phi)(v_1,v_2,\cdots, v_k)=\frac{d^2}{dt^2}\bigg\vert_{t=0} (f^*_t\phi)(v_1,v_2,\cdots, v_k)$$