Origin of integral of field strength tensor in path-ordered exponential in gauge field theory When studying some gauge theories approach to problems in Mechanics, I've found the following integral
$$P\exp\left[\oint A \ dt\right]=1+\dfrac{1}{2}\oint_{\partial D}\sum_{\mu,\nu}F_{\mu\nu}\gamma^{\mu}(t) \dot{\gamma}^{\nu}(t)dt,$$
where $A$ is the gauge potential, and $F$ is the field strength tensor (i.e. the pull-back of the curvature two form by a certain choice of gauge map). This integral appeared, in the articles (like this one, on page 564) I've seem, on the computation of the path-ordered exponential, but I couldn't understand where it comes from. It seems, on this formula, that we are integrating over a path, but $F$ is a $2$-form, so it should be integrated over a $2$-chain.
On the article there is one derivation, but I really didn't understand what they did, it doesn't seem very rigorous. Also, when I've studied principal fiber bundles and connections on those bundles, I didn't see this integral. I've also searched on some math books and didn't find it.
So, where this integral comes from, what it rigorously means and how it relates to the path-ordered exponential?
 A: It seems that you ask for the following thing. Suppose one has the Wilson operator
$$U\left(b,a\right)=P\exp\left[\dfrac{\mathbf{i}g}{\hbar}\int_{0}^{1}dz^{\mu}\left(s\right)A_{\mu}\left(z\left(s\right)\right)\right]$$
with $P$ the path-ordered operator, $A$ the gauge-potential (connexion 1-form) and $z\left(s\right)$ the path along which the integral is defined, from $s=0$ to $s=1$. Then the infinitesimal variation of $U$ with respect to its end-points reads
$$\delta U\left(b,a\right)=\dfrac{\mathbf{i}g}{\hbar}A\left(b\right)U\left(b,a\right)db-\dfrac{\mathbf{i}g}{\hbar}U\left(b,a\right)A\left(a\right)da\\-\dfrac{\mathbf{i}g}{\hbar}\int_{0}^{1}\left[U\left(b,z\left(s\right)\right)F_{\mu\nu}\left(z\left(s\right)\right)U\left(z\left(s\right),x\right)\right]\dfrac{dz^{\mu}}{ds}\left(\dfrac{dz^{\nu}}{db^{\lambda}}db^{\lambda}+\dfrac{dz^{\nu}}{da^{\lambda}}dx^{\lambda}\right)ds$$
where $F_{\mu\nu}=\partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}-\dfrac{\mathbf{i}g}{\hbar}\left[A_{\mu},A_{\nu}\right]$ the gauge-field strength. $g$ is a constant throughout, and $\hbar$ is for decoration here. 
The above formula is demonstrated in 


*

*J. Winter, Covariant Extension of the Wigner Transformation To Non-Abelian Yang-Mills Symmetries for a Vlasov Equation Approach To the Quark-Gluon Plasma Le J. Phys. Colloq. 45, C6.53 (1984).

*H.-T. Elze, M. Gyulassy, and D. Vasak, Transport equations for the QCD quark Wigner operator Nucl. Phys. B 276, 706 (1986).
see also 


*

*H. Weigert and U. Heinz, Kinetic equations for the quark-gluon plasma and their semiclassical expansion Zeitschrift Für Phys. C Part. Fields 50, 195 (1991).


I also have some notes about the derivation in a simpler (quasi-heuristic) way, but I'm not sure if it's exactly what you're looking for. It seems the formula I gave above contains yours (it's difficult to see since you didn't even provide an equality...
