The principle of learning sign rules by neural networks in qubit lattice models

A qubit is commonly used to represent a quantum state $\rvert n\rangle$ in various fields of condensed matter physics. Examples include a spin-$1/2$ in quantum magnets Vasiliev et al. (2018), a single fermion state in ultra-cold atomic systems Gardiner and Zoller (2017), a two-level atom in quantum cavities Meher and Sivakumar (2022), and so on Makhlin et al. (2001); Luo (2008); Kjaergaard et al. (2020). The binary value $n$ corresponds to an empty or occupied fermion level, or a spin-$1/2$ polarizing $\uparrow$ or $\downarrow$ in the z-axis. For a lattice with $L$ qubits, the basis can be expressed as $\rvert\mathbf{n}\rangle=\otimes^{L}_{l=1}\rvert n_{l}\rangle$, where $\rvert n_{l}\rangle$ represents the local basis at site-$l$, and the quantum indices $n_{l}=0$, $1$ form a vector $\mathbf{n}=(n_{1},\ \dots,\ n_{L})$.

Without loss of generality, let us consider a spin-$1/2$ as an example. A spin operator $\hat{\mathbf{S}}_{l}=(\hat{S}^{x}_{l}$, $\hat{S}^{y}_{l}$, $\hat{S}^{z}_{l})$ defined at site-$l$ has three components in the $x$, $y$ and $z$-axes, respectively. For any state, there are only two free real variables out of two complex coefficients in front of the basis $\rvert\sigma_{l}\rangle$ of the $\hat{S}^{z}_{l}$-representation. The index $\sigma_{l}=\uparrow$, $\downarrow$, corresponding to values of $\pm 1/2$. These variables are governed by a pair of site-dependent angles $\theta_{l}$ and $\phi_{l}$ appearing in a spin-coherent state Penc and Läuchli (2010)

$\displaystyle\rvert\Omega_{l}\rangle=c^{\uparrow}_{l}\left\rvert\uparrow\right\rangle+c^{\downarrow}_{l}\left\rvert\downarrow\right\rangle\ ,$

(1)

where the coefficients are given by

$\displaystyle c^{\uparrow}_{l}=\cos(\theta_{l}/2)\quad\text{and}\quad c^{\downarrow}_{l}=\sin(\theta_{l}/2)e^{i\phi_{l}}\ .$

(2)

As convention, $\theta_{l}\in[0$, $\pi]$ and $\phi_{l}\in[0$, $2\pi)$. In such a state, $\mathbf{S}_{l}=\langle\hat{\mathbf{S}}_{l}\rangle=\mathbf{\Omega}_{l}/2$ behaves as half of the unit vector in three-dimensional coordinates, that is,

$\displaystyle\mathbf{\Omega}_{l}=(\sin\theta_{l}\cos\phi_{l},\ \sin\theta_{l}\sin\phi_{l},\ \cos\theta_{l})\ .$

(3)

The phase factor for the basis $\rvert\sigma_{l}\rangle$ with a non-vanishing amplitude only depends on $\phi_{l}$, since both $\cos(\theta_{l}/2)$ and $\sin(\theta_{l}/2)$ are positive. Besides, a phase angle $h_{l}$ can modulate the phase factor in front of the spin-coherent state (1), i.e., $e^{ih_{l}}\rvert\Omega_{l}\rangle$.

In the GWMF theory Gutzwiller (1963, 1965), the wave function of the ground state $\rvert\psi^{\text{gw}}\rangle$ is a product of bases for $L$ spin-$1/2$s, i.e.,

$\displaystyle\rvert\psi^{\text{gw}}\rangle=\bigotimes^{L}_{l=1}\rvert\Omega_{l}\rangle=\sum_{\{\bm{\sigma}\}}\bar{a}^{\text{gw}}_{\bm{\sigma}}p^{\text{gw}}_{\bm{\sigma}}\rvert\bm{\sigma}\rangle\ ,$

(4)

where

$\displaystyle\bar{a}^{\text{gw}}_{\bm{\sigma}}=\prod_{l}\lvert c^{\sigma_{l}}_{l}\rvert\quad\text{and}\quad p^{\text{gw}}_{\bm{\sigma}}=\exp\left[i\left(\bm{\phi}\cdot\mathbf{n}+\sum_{l=1}^{L}h_{l}\right)\right]\quad$

(5)

represent the positive amplitude and the phase factor, respectively. Here, we define the angle vector as $\bm{\phi}=(\phi_{1},\ \cdots,\ \phi_{L})$, and the index vector $\bm{\sigma}=(\sigma_{1},\ \cdots,\ \sigma_{L})$ satisfies the corresponding relation $\mathbf{n}=1/2-\bm{\sigma}$. For a quantum model featuring a real Hamiltonian discussed in this work, the complex conjugate of the GWMF wave function $\rvert\psi^{\text{gw}}\rangle^{*}$ also indicates a state sharing the ground-state energy. Thus, the modified GWMF wave function $e^{ih_{0}}|\psi^{\text{gw}}\rangle=|\psi^{\text{gw}}_{r}\rangle+i|\psi^{\text{gw}}_{i}\rangle$ for the ground state can be decomposed into the real part $|\psi^{\text{gw}}_{r}\rangle$ and imaginary part $|\psi^{\text{gw}}_{i}\rangle$, given a specific global phase angle $h_{0}$. Both parts are real-valued and have an extra orthogonalization relation of $\langle\psi^{\text{gw}}_{r}|\psi^{\text{gw}}_{i}\rangle=\langle\psi^{\text{gw}}_{i}|\psi^{\text{gw}}_{r}\rangle=0$. Concretely, two parts are expressed as

$\displaystyle\begin{split}\rvert\psi^{\text{gw}}_{r}\rangle&=\sum_{\{\bm{\sigma}\}}\bar{a}^{\text{gw}}_{\bm{\sigma}}\cos(\bm{\phi}\cdot\mathbf{n}+\tilde{h})\rvert\bm{\sigma}\rangle\ ,\\ \rvert\psi^{\text{gw}}_{i}\rangle&=\sum_{\{\bm{\sigma}\}}\bar{a}^{\text{gw}}_{\bm{\sigma}}\sin(\bm{\phi}\cdot\mathbf{n}+\tilde{h})\rvert\bm{\sigma}\rangle\ ,\end{split}$

(6)

where the phase angle $\tilde{h}=\sum^{L}_{l=1}h_{l}+h_{0}$ is given. Regardless of whether the imaginary part $|\psi^{\text{gw}}_{i}\rangle$ is null or two parts share the same energy, we can always obtain a real ground-state wave function. It is worth noticing that the mean fields in the GWMF theory prefer selecting one of the degenerate manifolds if they exist, which artificially breaks the corresponding symmetry. So, the above-mentioned rotation, adjusted by the phase angle $h_{0}$, would introduce an energy split between the real part $|\psi^{\text{gw}}_{r}\rangle$ and the imaginary part $|\psi^{\text{gw}}_{i}\rangle$. In our theory, we ignore this effect and suppose their degeneracy still survives.

To assume that the real part $\rvert\psi^{\text{gw}}_{r}\rangle=\sum_{\{\bm{\sigma}\}}c^{\text{gw}}_{\bm{\sigma}}\rvert\bm{\sigma}\rangle$ can be expanded in the representation of bases $\rvert\bm{\sigma}\rangle$, the real expansion coefficient $c^{\text{gw}}_{\bm{\sigma}}=\bar{a}^{\text{gw}}_{\bm{\sigma}}\cos(\bm{\phi}\cdot\mathbf{n}+\tilde{h})=s^{\text{gw}}_{\bm{\sigma}}a^{\text{gw}}_{\bm{\sigma}}$ comprises of an amplitude $a^{\text{gw}}_{\bm{\sigma}}=\lvert c^{\text{gw}}_{\bm{\sigma}}\rvert$ and a sign $s^{\text{gw}}_{\bm{\sigma}}=\text{Sgn}(c^{\text{gw}}_{\bm{\sigma}})$ following a rule:

$\displaystyle s^{\text{gw}}_{\bm{\sigma}}=\text{Sgn}[\cos(\bm{\phi}\cdot\mathbf{n}+\tilde{h})]\ .$

(7)

The rule is called the leading-order version, removing short-range fluctuations completely. The phase angle $\tilde{h}$ is determined by other necessarily preserved global symmetries for a specified eigenstate, e.g., translational and inversion symmetries. And “Sgn” denotes the standard sign function. When examining the sign rule for the nonzero imaginary part $\rvert\psi^{\text{gw}}_{i}\rangle$, an extra $\pi/2$ needs to be added to the phase angle $\tilde{h}$, which equivalently replaces the cosine function with the sine function in Eq. (7). In the alternative notation-$\mathbf{n}$ of utilizing bases $\rvert\mathbf{n}\rangle$, the sign rule remains the same, i.e., $s^{\text{gw}}_{\mathbf{n}}\equiv s^{\text{gw}}_{\bm{\sigma}}$. For the convenience of the following discussions, we only use the notation-$\mathbf{n}$.

In the GWMF scenario, spins in the ordered state are typically visualized as classical vectors that follow a regular profile $\{\mathbf{\Omega}_{l}\}$ in space. Our derivation shows that the leading-order sign rule, which depends on angles $\phi_{l}$, is closely related to the spin-order profile $\{\mathbf{\Omega}_{l}\}$. The above conclusion is still valid for general qubit lattice models. In Sec. II, we will demonstrate that the leading-order sign rule can, in principle, be learned by shn-FNN.

The feed-forward neural network (FNN) is a powerful tool for approximating continuous functions Cybenko (1989); Hornik et al. (1989) and sorting samples by discrete values of characters Goodfellow et al. (2016). For instance, it can be applied to the classification of the double-valued sign $s_{\mathbf{n}}$ for arbitrary basis $\rvert\mathbf{n}\rangle$, when the expansion coefficient $c_{\mathbf{n}}$ in the ground-state wave function $\rvert\psi\rangle=\sum_{\mathbf{n}}c_{\mathbf{n}}\rvert\mathbf{n}\rangle$ for a real Hamiltonian consists of an amplitude $a_{\mathbf{n}}$ and a sign $s_{\mathbf{n}}$. However, as FNN becomes deeper, its complexity grows, making it difficult to understand the sign rule and its connection to meaningful physics insights. To address this issue, we introduce shn-FNN, similar to previous shallow FNNs He et al. (2020); Yuan and Weng (2021) but distinct from recently developed operations in a compact latent space Iten et al. (2020); Wang et al. (2021).

The shn-FNN consists of an input, hidden, and output layer, as illustrated in Fig. 1. The configuration $\mathbf{n}$ is assigned to the input layer of shn-FNN by simply setting a $L$-dimensional vector $\bm{y}_{I}=\mathbf{n}$. The hidden layer contains one neuron that produces a one-dimensional vector $y_{H}$. The output layer consists of two neurons that form a one-hot vector $\bm{y}_{O}=(y_{O,1},\ y_{O,2})$. These three layers are connected by two weight vectors $\bm{w}$ and $\bm{w}_{O}$.

The activation function cosine is empirically chosen for the hidden layer so that

$\displaystyle y_{H}=\cos(\bm{w}\cdot\mathbf{n})\ .$

(8)

The vector $\bm{y}_{O}$ in the output layer is determined by applying the softmax function, i.e.,

$\displaystyle\bm{y}_{O}=\text{softmax}(y_{H}\bm{w}_{O})\ ,$

(9)

where the weight vector $\bm{w}_{O}=(w_{O,1}$, $w_{O,2})$ is fixed. The function softmax executes normalization by an exponential function to obtain probabilities, which is usually used in classification tasks Goodfellow et al. (2016). Specifically, the function softmax gives two neurons $y_{O,1}$ and $y_{O,2}$ in the output layer, given by

$\displaystyle\begin{split}y_{O,1}&=\frac{e^{w_{O,1}y_{H}}}{e^{w_{O,1}y_{H}}+e^{w_{O,2}y_{H}}}\ ,\\ y_{O,2}&=\frac{e^{w_{O,2}y_{H}}}{e^{w_{O,1}y_{H}}+e^{w_{O,2}y_{H}}}\ .\end{split}$

(10)

In this work, we choose $\bm{w}_{O}=(1,-1)$, so that $e^{w_{O,1}y_{H}}$ and $e^{w_{O,2}y_{H}}$ are distinguishable. Thus, the desired sign $s^{\text{fnn}}_{\mathbf{n}}$ is determined as follows:

$\displaystyle s^{\text{fnn}}_{\mathbf{n}}=\left\{\begin{array}[]{lc}+1&\quad\text{for}\quad y_{O,1}>y_{O,2}\\ \\ -1&\quad\text{for}\quad y_{O,1}\leq y_{O,2}\end{array}\right.\ ,$

(14)

which equivalently indicates the sign of $y_{H}$, i.e., $s^{\text{fnn}}_{\mathbf{n}}\equiv\text{Sgn}(y_{H})$ in principle. Without one-hot representation, it has been proven that FNN performs worse since a categorization task turns into a regression task compulsively Goodfellow et al. (2016).

In summary, the shn-FNN representation for the sign rule corresponds to a function

$\displaystyle s^{\text{fnn}}_{\mathbf{n}}=\text{Sgn}\left[\cos(\bm{w}\cdot\mathbf{n})\right]\ .$

(15)

Compared to the leading-order sign rule (7) that is derived from the GWMF theory, we usually get

$\displaystyle w_{l}=\phi_{l}+\tilde{h}/N\quad\text{and}\quad N=\sum^{L}_{l=1}n_{l}$

(16)

for models with constant $N$. Thus, Eq. (15) is also called the leading-order sign rule.

Data sets. We use the exact diagonalization (ED) method to obtain the ground-state wave function. In the wave function, the sign $s_{\mathbf{n}}$ and the configuration $\mathbf{n}$ for a specified basis $\rvert\mathbf{n}\rangle$ constitute a sample in the data set $\mathbf{T}$. Each sign $s_{\mathbf{n}}$ is encoded as a one-hot vector $\bm{y}^{(\mathbf{n})}=(y^{(\mathbf{n})}_{1}$, $y^{(\mathbf{n})}_{2})$, which can only take two valid combinations:

$\displaystyle\bm{y}^{(\mathbf{n})}=\left\{\begin{array}[]{lc}(1,\ 0)&\quad\text{for}\quad s_{\mathbf{n}}>0\\ \\ (0,\ 1)&\quad\text{for}\quad s_{\mathbf{n}}<0\end{array}\right.\ .$

(20)

After arranging the samples in descending order of amplitude $a_{\mathbf{n}}$, we discard those with $a_{\mathbf{n}}<10^{-15}$ to avoid any artificial effects caused by the limited numeric precision. The remaining $\mathcal{N}_{s}$ samples in the data set $\mathbf{T}$ are divided into a training set $\mathbf{T}_{\text{train}}$ and a testing set $\mathbf{T}_{\text{test}}$ in a $4:1$ ratio Goodfellow et al. (2016). Thus, the number of samples in the testing set $\mathbf{T}_{\text{test}}$ is given by $\mathcal{N}_{\text{test}}=\mathcal{N}_{s}/5$.

Training. During the training scheme-$\rm{I}$, we employ the back-propagation (BP) algorithm Rumelhart et al. (1986) to optimize the variables in the weight vector $\bm{w}$ while adaptively adjusting the learning rate using the Adam algorithm P and Ba (2014). The process aims to minimize the cross entropy, defined as

$\displaystyle{\cal S}_{\times}=-\sum_{\{\mathbf{n}\}}\left(y^{(\mathbf{n})}_{1}\ln y^{(\mathbf{n})}_{O,1}+y^{(\mathbf{n})}_{2}\ln y^{(\mathbf{n})}_{O,2}\right)\ ,$

(21)

which sums over samples in the entire training set. Here, the one-hot vector $\bm{y}^{(\mathbf{n})}_{O}=(y^{(\mathbf{n})}_{O,1}$, $y^{(\mathbf{n})}_{O,2})$ is the output of shn-FNN as the vector $\mathbf{n}$ is input.

We employ the mini-batch method based on the stochastic gradient descent (SGD) Goodfellow et al. (2016); Wilson and Martinez (2003) to reduce the huge computational costs. Instead of using the entire training set $\mathbf{T}_{\text{train}}$ directly, we randomly select $\mathcal{N}_{\text{step}}=100$ samples to calculate the gradients of the weights $w_{l}$ at each training step. In such a case, Eq. (21) only sums over selected $\mathcal{N}_{\text{step}}$ samples. This method performs well in accuracy and speed, and the random selection helps prevent the “over-fitting” problem. In our program, we implement shn-FNN and the Adam optimization using the ML library “TensorFlow” Abadi et al. (2016).

To evaluate the performance of shn-FNN, we introduce the accuracy rate (AR), which is calculated as the ratio of the number of successfully classified samples $\mathcal{N}^{c}_{s}$ to the total number of samples $\mathcal{N}_{s}$ in the data set $\mathbf{T}$, i.e., $\text{AR}=\mathcal{N}^{c}_{s}/\mathcal{N}_{s}$. To monitor the optimization process, we define two additional accuracy rates. At each training step, we calculate the number of samples successfully classified $\mathcal{N}^{c}_{\text{step}}$. Then, we evaluate the optimization by computing the accuracy rate $\text{AR}_{\text{step}}=\mathcal{N}^{c}_{\text{step}}/\mathcal{N}_{\text{step}}$. To provide a comprehensive evaluation, we utilize the testing set $\mathbf{T}_{\text{test}}$ and define the accuracy rate $\text{AR}_{\text{test}}=\mathcal{N}^{c}_{\text{test}}/\mathcal{N}_{\text{test}}$, where $\mathcal{N}^{c}_{\text{test}}$ represents the number of correctly classified samples in the testing set $\mathbf{T}_{\text{test}}$.

After each training step, we assess the convergence criterion $\epsilon$, which measures the absolute value of the difference between the accuracy rates $\text{AR}_{\text{test}}$ obtained in the current step and the previous step. We halt the training process once the convergent criterion $\epsilon$ falls below a threshold $\epsilon_{0}=10^{-3}$.

To exemplify how shn-FNN learns the sign rule, we study the ground state of a generalized Ising ring described by the Hamiltonian

$\displaystyle\hat{H}_{\text{Ising}}(J,\ \mathbf{\Omega})=J\sum^{L}_{l=1}(\hat{\mathbf{S}}_{l}\cdot\mathbf{\Omega})(\hat{\mathbf{S}}_{l+1}\cdot\mathbf{\Omega})\ ,$

(22)

where $J$ represents the strength of the coupling, and $\mathbf{\Omega}$ determines the orientation of the spin polarization.

For the case of ferromagnetic coupling ($J<0$) and spin polarization along the x-axis ($\mathbf{\Omega}=\hat{x}$), we obtain the ground-state wave function for $L=16$ spins, and then prepare the data sets $\mathbf{T}$, $\mathbf{T}_{\text{train}}$ and $\mathbf{T}_{\text{test}}$, as stated above. We initialize the weight vectors $\bm{w}$ and $\bm{w}_{O}$ in shn-FNN and start the training process. As illustrated in Fig. 2, the cross entropy $\mathcal{S}_{\times}$ rapidly decreases around the $170$-th training step. After approximately $200$ training steps, the cross entropy $\mathcal{S}_{\times}$ tends towards stability, and the three accuracy rates $\text{AR}_{\text{step}}$, $\text{AR}_{\text{train}}$ and AR consistently reach a maximum value $1$ or $100\%$. We terminate the training process when the convergence criterion $\epsilon<\epsilon_{0}$ is met. In addition to scheme-$\rm{I}$, we will introduce scheme-$\rm{II}$ in Sec. IV.2.1, tailored for frustrated spin models and the Fermi-Hubbard ring later.

Using shn-FNN, we analyze the leading-order sign rules for various ordered ground states in qubit lattice models, including non-frustrated spin models in Sec. IV.1, frustrated spin models in Sec. IV.2, and interacting fermions in Sec. IV.3.

We have successfully developed a Gutzwiller mean-field theory of sign rules for the ordered ground states in qubit lattice models, which perfectly matches the sign predicted by a shallow FNN with a single hidden neuron, called shn-FNN. By utilizing this principle, we provide a consistent explanation for the excellent performance of activation functions in the neural network and offer a vivid interpretation of the sign rule represented by FNN.

We systematically test our theory on various spin models and the Fermi-Hubbard ring. For non-frustrated spin-$1/2$ models, such as a generalized Ising ring, (twisted) XY rings, and an antiferromagnetic Heisenberg ring, the sign rules for ground states with magnetic orders can be fully captured by shn-FNN, where the accurate rate of the prediction can archive $1$ exactly. However, in the case of frustrated models where interactions compete, the complexity of sign rules for ground states is significantly enhanced, reducing prediction accuracy. Nonetheless, the leading-order sign rules obtained by optimizing shn-FNN still provide a visual scenario of orders in spins, with the characteristic weight vector closely related to pitch angles. In the Fermi-Hubbard ring, we can obtain a unified sign rule by selecting suitable bases.

GWMF may not be suitable for 1D models since quantum fluctuations tend to destroy long-range orders. However, our current work presents a fresh perspective by demonstrating that GWMF can effectively capture the leading-order sign rule in the wave function, where fluctuations in amplitudes are erased. Our theory is a simple starting point by removing short-range details in ordered states. It would be intriguing to explore the information encoded in high-order microscopic processes instead of focusing solely on the leading-order ones. Of course, the theory for general lattice models also deserves profound studies in the future.

We thank Tao Li, Rui Wang, Ji-Lu He, Wei Su, and Wei Pan for the grateful discussion. S. H. acknowledges funding from the Ministry of Science and Technology of China (Grant No. 2022YFA1402700) and the National Science Foundation of China (Grants No. 12174020). Z. P. Y. acknowledges funding from the National Science Foundation of China (Grants No. 12074041). S. H. and K. X. further acknowledge support from Grant NSAF-U2230402. The computations were performed on the Tianhe-2JK at the Beijing Computational Science Research Center (CSRC) and the high-performance computing cluster of Beijing Normal University in Zhuhai.