loops — which arise due to the non-trivial holonomies of the torus. For an Abelian gauge theory the Polyakov loops are defined as $$\Phi_x(y) = \int_0^{L_x} dx \ A_x(x,y) - \varphi_x(y), \
\Phi_y(x) = \int_0^{L_y} dy \ A_y(x,y) - \varphi_y(x).$$ In our case, they are given by $$\Phi_x(y) = B L_x y - \frac{\theta_x}{e}, \
\Phi_y(x) = B L_y x - \frac{\theta_y}{e}.$$

In order to respect gauge invariance of the theory on the torus, under shifts the wave function must also be gauge transformed accordingly $$\begin{aligned}
\label{boundary}
&&\Psi(x + L_x,y) = \exp\left(i e \varphi_x(y)\right) \Psi(x,y) =
\exp\left(i \theta_x - i e B L_x y\right) \Psi(x,y), 
\nonumber \\
&&\Psi(x,y + L_y) = \exp\left(i e \varphi_y(x)\right) \Psi(x,y) =
\exp\left(i \theta_y\right) \Psi(x,y).\end{aligned}$$ The angles $\theta_x$ and $\theta_y$ parametrize a family of self-adjoint extensions of the Hamiltonian on the torus. Applying the boundary conditions from above in two different orders one obtains $$\begin{aligned}
\Psi(x + L_x,y + L_y)&=& 
\exp\left(i \theta_x - i e B L_x (y + L_y)\right) 
\Psi(x,y + L_y) \nonumber \\
&=&\exp\left(i \theta_x + i \theta_y - 
i e B L_x (y + L_y)\right) \Psi(x,y), \nonumber \\
\Psi(x + L_x,y + L_y)&=&\exp\left(i \theta_y\right) \Psi(x + L_x,y)
\nonumber \\
&=&\exp\left(i \theta_x + i 