quantum mechanics - Calculating the Berry curvature in case of degenerate levels (Non abelian Berry curvature): issue

The Berry phase accumulated on a path can be described by a matrix when we look at adiabatic time evolution with a Hamiltonian with degenerate energy levels. The Berry phase matrix is given by $$\gamma_{mn}= \int_\mathcal{C} \left\langle m(R) \right | \nabla_R \left| n(R) \right \rangle . d R.$$

here $R$ parametrizes the said path and $A_{mn}= \left\langle m(R) \right | \nabla_R \left| n(R) \right \rangle$. Now what I want to do is calculate the Berry Curvature, something that, if I assume my path above is closed and has three determining coordinates$R_1$, $R_2$ , $R_3$ is $\vec{F}$ such that

$$\mathbf{\gamma}=\int_\mathcal{S} \vec{F}.d\vec{s}$$ note that the $\gamma$ and $F$ here are matrices and we're integrating over the surface $\mathcal{S}$ enclosed by curve $\mathcal{C}$

What is stopping me from applying the stokes theorem to $\gamma_{mn}$ and getting $\vec{F}_{mn}=\nabla_R\times A_{mn}$?

It is said that the answer contains a matrix commutator $[A_i,A_j]_{mn}$ c because this berry phase is non abelian. but I seem to be missing something fundamental.

Edit: Note: this also corresponds to problem 2 Chapter 2 of Topological insulators and superconductors by Bernevig and Hughes

Answer

The short answer is that only if the Berry curvature is defined by: (in matrix notation): $$F_{\mu \nu} = \partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu} + [A_{\mu}, A_{\nu}]$$ it becomes gauge covariant, i.e., for a gauge transformation: $$A_{\mu} \rightarrow g^{-1}A_{\mu} g+g^{-1}\partial_{\mu}g$$ $g \in U(N)$ ($N$ is the degeneracy of the level), the curvature transforms as: $$F_{\mu \nu} \rightarrow g^{-1}F_{\mu \nu} g$$ The partial expression $F_{\mu \nu} = \partial_{\mu}A_{\nu}-\partial_{\nu}A_{\mu}$ will be covariant (actually invariant) only with respect to Abelian gauge transformations and not the whole $U(N)$ group.

The reason that the Berry curvature cannot be obtained as in the Abelian case by a direct application of the Stokes theorem is that the Stokes theorem does not exist(in the usual sense) in the non Abelian case becuse the holonomy needs to be path ordered:

$$\mathcal{Hol}(A) = \mathrm{P} e^{\int_C A_{\mu}dx^{\mu}}$$

However, it can be applied sequentially for an infinitesimal path. let us choose a square path in the $(xy)$ plane centered at $(x,y)$ with infinitesimal sides of length $\Delta x$ and $\Delta y$ as depicted in the picture.

Therefore

$$\begin{matrix} \mathcal{Hol}(A) \approx e^{\int_a^b A_{\mu}dx^{\mu}} e^{\int_b^c A_{\mu}dx^{\mu}} e^{\int_c^d A_{\mu}dx^{\mu}} e^{\int_d^a A_{\mu}dx^{\mu}} \\\approx (1+A_x\Delta x-\frac{1}{2}\partial_yA_x\Delta x \Delta y) (1+A_y\Delta y+\frac{1}{2}\partial_x A_y\Delta x \Delta y) (1-A_x\Delta x-\frac{1}{2}\partial_yA_x\Delta x \Delta y)(1-A_y\Delta y+\frac{1}{2}\partial_x A_y\Delta x \Delta y) \\ \approx 1+(\partial_x A_y - \partial_y A_x + [A_x, A_y])\Delta x \Delta y\approx (1+ F_{xy}\Delta x \Delta y) \end{matrix}$$

Please notice that the commutator term is created due to the need to exchange the order of the second and the third expressions in the product in order to cancel $A_x \Delta_x$.