See Refs. [23, 24] for recent reviews and a more complete list of references for the modern numerical conformal bootstrap.

Here, $P_p^2$ contains the 2-RDMs that encode all $p$-particle correlations in an $N$-particle system. Thus, we still have $\mathrm{Tr}[{}^{2}D]=N(N-1)/2,\forall {}^{2}D\in P_p^2$. It should really be labeled by $P_{p,N}^2$, but we drop the $N$ for simplicity.

In the way we defined $(2,p)$ constraints, this statement is trivial. However, we can define $(2,p)$ constraints to be those constructed from operators which are polynomials whose order is exactly $p$ rather than at most $p$. In this case, it is not obvious that $(2,p)$ constraints contain $(2,p^{\prime })$ constraints for $p^{\prime }<p$, but the statement still turns out to be true.

For integrals involving $V_{\mathit{m}q}$, which is always accompanied by the Gaussian factor, we can always assume they are replaced by a properly periodized version. This periodization does not make a difference for sufficiently large systems due to the rapidly decaying Gaussian.

We always assume our states satisfy all symmetries of the problem. For symmetry-broken states, we can always consider a mixed state averaged over all possible symmetry sectors.

See Appendix pp2 for more details.

Since we always use an odd $\times$ odd grid, we cannot do exactly half filling, but for sufficiently large system sizes, we can get arbitrarily close.

The peak in $S_q$ obtained in our approach is around $1.8k_F$.

For bootstrap, we have $E_{32/63}/V_s=0.808$ and $E_{31/63}/V_s=0.77$, while for VMC, we have $E_{32/64}/V_s=0.805$ and $E_{31/62}/V_s=0.801$.

In fact, when the filling is $1/5$ and the Hamiltonian is $V_1+V_3$, we still get a set of degenerate ground states after the bootstrap with $(2,2)+T2$. The interior point method will always pick up a state that is very close to the Laughlin state.

In this section, we neglect the sparsity of the matrices which further reduces the complexity.

This can also be confirmed from the bootstrap results with no mirror symmetry imposed which show automatically the mirror symmetry.

The particle number matters a lot in the optimization but does not affect the computational complexity since it does not change the dimensionality of the constraint matrices or the optimization variables.

We checked that the oscillating structure factor at $\nu =\frac{1}{3}$ has a peak at $q_{*}\approx 1.59{\ell }_B^{-1}$. A fictitious CFL at $\nu =\frac{1}{3}$ will have $\pi k_F^2=2\pi {\ell }^{-2}\nu$ such that $2k_F\approx 1.62$. A triangular lattice Wigner crystal, with translation breaking wave vector determined by $(\sqrt{3}/2)q_{\mathrm{wc}}^2=2\pi {\ell }_B^{-2}$, is another possible competing state. Here, the expected peak is $q_{\mathrm{wc}}\approx 1.56{\ell }_B^{-1}$, which is also close to the bootstrap peak.

If the ground state breaks such symmetry, we can always form a mixed state density matrix which averages over all symmetry sectors. Note that our approach does not distinguish pure and mixed density matrices.