Combinatorial exploration of quantum spin liquid candidates in the herbertsmithite material family

TL;DR

This study systematically explores structurally related kagome lattice compounds using first-principles calculations to identify promising quantum spin liquid candidates among herbertsmithite-like materials.

Contribution

It introduces a combinatorial first-principles approach to evaluate and predict new quantum spin liquid candidates within the herbertsmithite material family.

Findings

Identified several compounds with favorable structural and thermodynamic properties for quantum spin liquids.

Assessed synthesis feasibility through convex hull calculations.

Evaluated disorder and distortion likelihood in candidate materials.

Abstract

Geometric frustration of magnetic ions can lead to a quantum spin liquid ground state where long range magnetic order is avoided despite strong exchange interactions. The physical realization of quantum spin liquids comprises a major unresolved area of contemporary materials science. One prominent magnetically-frustrated structure is the kagome lattice. The naturally occurring minerals herbertsmithite [ZnCu$_3$(OH)$_6$Cl$_2$] and Zn-substituted barlowite [ZnCu$_3$(OH)$_6$BrF] both feature perfect kagome layers of spin-$1/2$ copper ions and display experimental signatures consistent with a quantum spin liquid state at low temperatures. To investigate other possible candidates within this material family, we perform a systematic first-principles combinatorial exploration of structurally related compounds [$A$Cu$_3$(OH)$_6B_2$ and $A$Cu$_3$(OH)$_6BC$] by substituting non-magnetic divalent cations ($A$) and halide anions ($B$, $C$). After optimizing such structures using density functional theory, we compare various structural and thermodynamic parameters to determine which compounds are most likely to favor a quantum spin liquid state. Convex hull calculations using binary compounds are performed to determine feasibility of synthesis. We also estimate the likelihood of interlayer substitutional disorder and spontaneous distortions of the kagome layers. After considering all of these factors as a whole, we select several promising candidate materials that we believe deserve further attention.

Tables

Table 1. Table 1: Properties of the most promising QSL candidate materials as compared to the reference materials. The references (herbertsmithite and Zn-barlowite) are highlighted in gray, and the final candidates (with no instabilities throughout the Brillouin zone) are marked with asterisks.

Compound	$f_0$ (THz)	$E_{\mathrm{hull}}$ (meV/atom)	$E_d^f$ (eV)	$\theta$ (deg)	$d_{\mathrm{inter}}$ (Å)	$d_{\mathrm{in}}$ (Å)
BaCu3(OH)6I2	0.41	42.6	2.42	128.0	6.09	3.53
CaCu3(OH)6Br2	0.50	30.7	0.87	125.7	5.19	3.53
CaCu3(OH)6Cl2	0.70	44.8	0.57	125.8	5.06	3.51
MgCu3(OH)6Br2 *	2.23	36.0	0.36	125.2	4.65	3.57
ZnCu3(OH)6Cl2	2.63	41.2	0.13	125.0	4.58	3.53
CaCu3(OH)6IBr	0.77	31.6	0.74	127.7	5.20	3.58
CaCu3(OH)6ICl *	0.94	19.2	0.72	125.4	5.17	3.54
MgCu3(OH)6ClF *	1.09	39.6	0.39	118.1	4.60	3.38
MgCu3(OH)6BrCl	0.35	26.9	0.30	126.1	4.61	3.56
ZnCu3(OH)6BrF	1.41	38.6	0.10	118.0	4.69	3.39
ZnCu3(OH)6ClF	0.89	43.1	0.07	118.5	4.64	3.38

Equations

$$E_d^f=E[\mathrm{defect}]-E[\mathrm{bulk}]+({\mu }_A-{\mu }_{\mathrm{Cu}})=\Delta E_s+\Delta \mu ,$$

$${\mu }_A+3{\mu }_{\mathrm{Cu}}+6{\mu }_{\mathrm{OH}}+2{\mu }_B>E[A{\mathrm{Cu}}_3(\mathrm{OH}{)}_6{B}_2].$$

$${\mu }_A+2{\mu }_B<E[A{B}_2],$$

$${\mu }_{\mathrm{Cu}}+2{\mu }_B<E[\mathrm{Cu}{B}_2],$$

$${\mu }_A+2{\mu }_{\mathrm{OH}}<E[A(\mathrm{OH}{)}_2].$$

Full text

I Introduction

In a quantum spin liquid (QSL), frustrated antiferromagnetic exchange interactions prevent localized spins from ordering at low temperatures, instead forming a fluid-like phase. The large degeneracy of this state can give rise to novel phenomena such as fractionalized quasiparticles, emergent gauge fields, and long-range entanglement [1, 2, 3, 4]. The kagome lattice of corner-sharing triangles is known to have high geometric frustration and is capable of hosting such a phase. A leading QSL material candidate possessing this structure is herbertsmithite [ZnCu3(OH)6Cl2], which contains perfect kagome layers of spin-$1/2$ copper cations separated by non-magnetic Zn and Cl ions [5, 6], as shown in Fig. 1(a,c). Indeed, although herbertsmithite has strong antiferromagnetic exchange interactions, no magnetic phase transition is observed down to sub-kelvin temperatures [7, 8, 9, 10, 11], and an array of experimental and theoretical work favors a possible QSL scenario [12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 19, 24].

Despite its many promising features, herbertsmithite is prone to cation substitutional disorder, where Cu may occupy interlayer sites and Zn may occupy intralayer kagome sites [7, 14, 25]. The precise amount of this disorder is debated. Several studies suggest that while there is minimal substitution of Zn on the kagome layers, the interlayer sites can be occupied by up to 15% Cu [26, 27, 12, 28], resulting in a decidedly off-stoichiometric compound. These interlayer “orphan” spin-$1/2$ Cu2+ defects are highly problematic for the QSL state, causing weak ferromagnetic interactions between kagome layers and distorting the surrounding matrix of magnetic ions [13]. Zn-substituted barlowite (Zn-barlowite), a structurally related compound and another potential QSL candidate [29, 30], is thought to have a much lower interlayer disorder concentration, largely due to the greater chemical distinction between the interlayer and intralayer sites, as shown in Fig. 1(b,d) [31, 32]. Experiments indicate that in Zn-barlowite, off-center interlayer $C_{2v}$ sites can contain up to  5% Cu defects. Like herbertsmithite, however, Zn-barlowite does not order magnetically, even with these large concentrations of magnetic defects [33, 34, 35]. While progress on this class of materials is encouraging, it is nevertheless desirable to further minimize orphan Cu spins to realize a clean QSL ground state.

Synthesizing compounds structurally similar to herbertsmithite and Zn-barlowite is a promising route to discover new QSL candidates. For example, Mg-substituted herbertsmithite, MgxCu4-x(OH)6Cl2 (tondiite), has been successfully synthesized and shows no magnetic phase transition down to 1.8 K [36, 37, 38], and a Cd analog [CdCu3(OH)6Cl2] shows no magnetic ordering down to 2 K, although it exhibits significant distortions of the kagome planes [39]. Synthesis of the bromide analog of herbertsmithite [ZnCu3(OH)6Br2] was attempted but unsuccessful [40]. A Zn-barlowite related structure, Zn-claringbullite [ZnCu3(OH)6ClF], shows no obvious magnetic transition down to 2 K, but a perfectly stoichiometric compound was not achieved [41]. While the Mg analog of barlowite cannot be synthesized due to the insolubility of MgF2 in water, the bromide analog was attempted [MgCu3(OH)6Br2], but did not have the Zn-barlowite structure and ordered antiferromagnetically at 5.4 K [42].

Clearly, more work is needed to search for and identify viable candidates in this material family. Only a few computational studies exist exploring cation substitution in barlowite [31, 32], and a complete exploration of the structural families of herbertsmithite and Zn-barlowite using computational methods has not been performed. In this paper, we use ab initio calculations to systematically explore compounds within the herbertsmithite and Zn-barlowite families. We compare the thermodynamic stability, structural properties, and tendency towards disorder. After considering all these criteria together, we select promising QSL candidates that merit further experimental and theoretical examination.

II Computational Procedure

We carry out a systematic exploration of the structural relatives of herbertsmithite [$A$Cu3(OH)${}_{6}B_{2}$] and Zn-barlowite [$A$Cu3(OH)${}_{6}BC$] by substituting closed-shell (spinless) $2+$ cations ($A$ = Ba, Be, Ca, Cd, Ge, Hg, Mg, Pb, Sn, Sr, Zn) and halide anions ($B,C$ = Br, Cl, F, I). We investigate all 44 possible herbertsmithite relatives. While there are 176 possible Zn-barlowite relatives, we eliminate compounds where $B=C$ because the herbertsmithite structure always has lower energy in these cases. We also do not consider compounds in which the less electronegative anion occupies the $C$ site [the site occupied by F in Fig. 1(b,d)]. All hydrogen bonds are oriented towards the $C$ site, so the more electronegative ion will always occupy this position to minimize energy. Thus, a total of 66 relatives in the Zn-barlowite family were selected for consideration.

We perform high-throughput calculations where the structural optimization of each candidate is followed by a static calculation to extract the ground-state energy and to compute phonon frequencies at the $\Gamma$ point to confirm structural stability. In addition to confirming the stability of the relaxed structures, we perform convex hull calculations to determine if synthesis of the candidate compounds is thermodynamically feasible. For the most promising materials, we also calculate defect formation energies and full phonon dispersions throughout the first Brillouin zone to verify stability at $k$-points away from the zone center.

All structures were calculated by allowing the lattice parameters, cell volume, and atomic positions to fully relax using density functional theory (DFT) as implemented in the Vienna ab initio simulation package (vasp) [43, 44, 45]. We used the supplied projector augmented wave potentials [46] within the generalized gradient approximation and Perdew-Burke-Ernzerhof scheme [47]. Electronic wave functions were expanded in a plane wave basis set with an energy cutoff of 800 eV, and reciprocal space was sampled using an $8\times 8\times 8$ $k$-point mesh for herbertsmithite-related structures and an $8\times 8\times 5$ $k$-point mesh for Zn-barlowite-related structures. A $\Gamma$-centered mesh is necessary due to the hexagonal symmetry of Zn-barlowite. The spacing between $k$-points was $\sim$0.15 Å*-1* for both structural families, and this spacing was also used for calculating the energies of binary compounds used in the convex hull analysis. All structures were relaxed until forces on the atoms were less than 1 meV/Å. Calculations were non-spin-polarized. Input files for all calculations can be found in the Supplemental Material [48].

III Results and Discussion

III.1 Phonon Calculations

Phonon calculations at the $\Gamma$ point for the fully-relaxed structures were performed in vasp within the finite differences approximation to confirm structural stability. As expected, many structures have unstable phonon modes. Fig. 2(a,b) shows the frequency of the lowest energy optical phonon mode, $f_{0}$, for all compounds. In all subsequent plots, the unstable compounds (with $f_{0}<0$) are marked with an ‘X’ to distinguish them from structurally stable and potentially viable candidates. Cations are shown on the vertical axis and anions on the horizontal axis, in order of increasing ionic radius from bottom to top and left to right, respectively. The reference compound, either herbertsmithite or Zn-barlowite, is shown in white and marked with an asterisk. Compounds with parameter values more favorable than the reference compound are shown with warm colors, and values less favorable are shown with cool colors. For example, a higher frequency of the lowest energy optical mode indicates higher dynamical stability, so higher frequencies are shown with warm colors. Most compounds containing group IV elements (Ge, Sn, Pb) tend to be unstable, with the exception of GeCu3(OH)6F2 and PbCu3(OH)6F2. Compounds containing larger cations are generally unstable, as well as Zn-barlowite relatives containing Be.

III.2 Convex Hull Calculations

The convex hull of a compound is useful for determining if synthesis is thermodynamically feasible, usually through a comparison of the compound’s formation energy to the sum of the energies of all other possible combinations of crystal structures that could be created from the same set of elements in the same ratios. Due to the prohibitive size of the phase space for our candidate materials, we perform a simplified procedure. Instead of considering all possible crystal structures, we consider only simple binary ionic compounds [e.g. $A$(OH)2, $AB_{2}$], which are most likely to yield the lowest convex hull energies (see Supplemental Material [48]). Starting structures for these binary compounds were obtained from the Materials Project [49] and then re-relaxed with our settings.

Insulators with energies less than $\sim$50 meV above the convex hull tend to be stable [50]. We therefore use an energy cutoff of 50 meV/atom as our criteria for thermodynamic stability when identifying candidate materials. The calculated energy above the hull for each compound is shown in Fig. 2(c,d). Energies higher than the reference compound are considered unfavorable and are represented with cool colors, while energies lower than the reference compound are favorable and represented with warm colors. Again, the reference compounds are shown in white and marked with an asterisk, and compounds with structural instabilities (as determined by phonon calculations) are marked with an ‘X’. There does not appear to be a clear connection between convex hull energy and structural stability or ion size.

III.3 Comparing Structural Parameters

In addition to structural and thermodynamic stability, we use Cu-O-Cu bond angles and spacings between kagome layers as additional metrics to rank the candidate compounds. A Cu-O-Cu bond angle approaching 180*∘* leads to a large antiferromagnetic superexchange interaction while minimizing undesirable Dzyaloshinskii–Moriya interactions. Larger bond angles are therefore highly desirable. Additionally, a greater separation between the kagome layers isolates the two-dimensional magnetic subsystems and suppresses unwanted coupling between planes. In Fig. 3, these two structural properties are displayed for all candidate compounds. Squares corresponding to specific compounds are colored and marked according to the same system described for Fig. 2, where bond angles and interplane distances larger (smaller) than the reference compounds are favorable (unfavorable) and represented with warm (cool) colors, and structurally unstable compounds continue to be marked with an ‘X’. Compounds with larger cation and anion radii generally lead to larger bond angles and interplane distances, but also tend to be structurally unstable. Compounds containing group IV elements are unstable and tend to have smaller bond angles.

In Fig. 4, we investigate the effects of ion size on the physical properties of the candidate compounds in more detail. In Fig. 4(a), the Cu-O-Cu bond angle is plotted versus anion radius for the structurally stable materials. The anion size plotted on the horizontal axis for Zn-barlowite relatives refers to the $C$-site anion that occupies the same position as F in the reference compound [ZnCu3(OH)6BrF] because it has the largest influence on bond angle. For all materials, bond angle increases with increasing anion size, and for a given anion, the bond angle also increases with increasing cation size. Figure 4(b) shows the kagome plane spacing versus cation radius for stable compounds, with separate traces for each anion. As expected, a larger cation radius leads to greater distance between the kagome layers. For a given cation, interplane distance also increases with increasing anion size. In Fig. 4(c), we find that while the $C$-site anion has the greatest effect on the Cu-O-Cu bond angle, larger bond angles are obtained when the $B$-site anion is similar in size to the $C$-site anion.

We examine the effect of ion size on the lattice parameters of stable compounds in Fig. 4(d). The $c$-axis length primarily increases with cation size while the $a$-axis length primarily increases with anion size, although anion size has a much weaker affect on the $a$-axis than cation size does on the $c$-axis. The frequency of the lowest optical phonon mode ($f_{0}$) is plotted against $c$-axis length in Fig. 4(e) for both stable (filled markers) and unstable (empty markers) structures. Of all the structural parameters, the $c$-axis length has the highest correlation with $f_{0}$. For herbertsmithite relatives, as the $c$-axis increases, $f_{0}$ decreases, meaning compounds tend to be less dynamically stable. Compounds containing group IV ions (Ge, Sn, Pb) are plotted in darker shades for both structural families because nearly all compounds containing these elements are unstable. Of the compounds not containing group IV ions, $c$-axis lengths that are very small or very large lead to structural instabilities. Compounds containing cations from groups IIA and IIB which are close in size to Zn tend to be most stable. Fig. 4(f) shows Cu-O-Cu bond angle versus $a$-axis length. We find that a larger $a$-axis leads to a larger bond angle, which agrees with the results in Fig. 4(a), where bond angle is positively correlated with anion radius, and Fig. 4(d), which shows the positive correlation between anion size and the length of the $a$-axis. It should be noted that many unstable compounds containing group IV elements have much smaller bond angles than most other candidates.

We also explored correlations between Cu-O-Cu bond angle, interplane distance, and in-plane Cu-Cu bond length. These plots can be found in the Supplemental Material [48]. The Cu-O-Cu bond angle has a weak positive correlation with interplane distance. There is also a positive correlation between in-plane Cu-Cu distance and Cu-O-Cu bond angle, as both are influenced by the length of the $a$-axis, which increases with increasing anion size. There is no obvious correlation between the interplane kagome distance and the in-plane Cu-Cu bond length, as the interplane distance depends mostly on cation size, and in-plane bond length depends on anion size. Overall, for both structural families, compounds with cations of intermediate size (Mg, Zn, Cd, and Hg) are most stable. Compounds containing group IV elements (Ge, Sn, Pb) are mostly unstable. Larger anions and cations lead to favorable structural properties, such as larger bond angles and interplane distances, but may also lead to distortions of the kagome layers or other structural instabilities.

III.4 Defect Formation Energy

Herbertsmithite and Zn-barlowite are both susceptible to cation disorder. In herbertsmithite, the Jahn-Teller active $d^{9}$ Cu2+ ion occupies the tetragonally elongated site in the center of the CuO4Cl2 octahedra. The $d^{10}$ Zn2+ ions are not Jahn-Teller active, and occupy the higher-symmetry trigonally compressed octahedral sites between the kagome layers. Due to the electronic configurations of the ions and distinct coordination environments, it is not favorable for Zn to occupy the in-plane sites within the kagome layer. However, herbertsmithite is the $x=1$ end member of the Zn-paratacamite family [ZnxCu4-x(OH)6Cl2], and there is a preference for some Cu to exist on the interlayer site instead of full occupation with Zn alone [7]. The equilibrium occupation of the interlayer site by Cu has been estimated to be as large as 15% in herbertsmithite [26, 27].

In Zn-barlowite, the interlayer site has a trigonal prismatic geometry, making it even less favorable for the Jahn-Teller active Cu2+ ion. As a result, the interlayer Cu occupation is only $\sim$5% in Zn-barlowite [33], confirming early computational predictions [31, 32]. Site-specific x-ray diffraction measurements have shown that there are two distinct interlayer sites in Zn-barlowite: an off-center $C_{2v}$ site and a central $D_{3h}$ site. The interlayer Cu defects occupy the $C_{2v}$ sites. It should be noted that even for large concentrations of magnetic impurities on the interlayer site, Zn-barlowite does not show signs of magnetic ordering, indicating that the possible QSL phase is somewhat robust against interlayer magnetic impurities [33].

An ideal QSL candidate will have only non-magnetic ions on the interlayer sites, and therefore must have a high energy cost for interlayer Cu substitution. We calculated the formation energy of such defects in a select number of our most promising candidates (those structurally stable, with $E_{\mathrm{hull}}<50$ meV/atom, and with bond angles and interplane distances larger than the reference compounds). Since nearly all experimental and computational studies indicate that there is negligible substitution of non-magnetic ions within the kagome layers, we consider only interlayer defects. The general expression for the formation energy of a charge-neutral substitutional defect is

[TABLE]

where $\Delta E_{s}$ is the difference in energy between a structure with a single defect and the pristine bulk structure and $\Delta\mu$ is the chemical potential difference of $A$ and Cu. To calculate $E[\mathrm{defect}]$, we construct defect structures from $2\times 2\times 2$ supercells of herbertsmithite relatives and $2\times 2\times 1$ supercells of Zn-barlowite relatives, with a single Cu substitution. A depiction of our defect configuration can be found in the Supplemental Material [48]. We relax the atomic positions of the defect structures and subtract the energy of the original defect-free structure to obtain $\Delta E_{s}$.

The chemical formulas for the defect-containing and defect-free configurations are not equivalent, so the chemical potential difference $\Delta\mu=\mu_{A}-\mu_{\mathrm{Cu}}$ must be considered. Interlayer defects are primarily created during the initial growth of the material. During synthesis of $A$Cu3(OH)${}_{6}B_{2}$, the chemical potentials of the constituent elements must satisfy the inequality

[TABLE]

Individual chemical potentials must all be less than zero ($\mu_{A}<0$, $\mu_{B}<0$, $\mu_{\mathrm{OH}}<0$, and $\mu_{\mathrm{Cu}}<0$). Additionally, the formation of unwanted side products must be avoided, imposing the additional inequalities

[TABLE]

[TABLE]

[TABLE]

Similar inequality constraints exist for $A$Cu3(OH)${}_{6}BC$. A higher defect formation energy is preferable to minimize disorder. To maximize $E_{d}^{f}$, we must maximize the chemical potential difference $\Delta\mu$ subject to the above inequality constraints. The defect formation energies calculated with these optimal values of $\Delta\mu$ are given in Table 1. All candidate compounds investigated had a higher energy cost for interlayer defects than herbertsmithite and Zn-barlowite except ZnCu3(OH)6ClF (Zn-substituted claringbullite).

Two previous computational studies investigated doping selectivity in barlowite [31, 32]. In both cases, the authors investigated the likelihood of substituting various non-magnetic ions into the interlayer and intralayer sites of barlowite, in contrast to the present work where we examine the energy cost of a Cu defect on an interlayer site in fully-substituted $A$-barlowite ($A$ = Zn, Mg, Ca). Despite differences in the methodology used to construct defect structures and calculate the chemical potential differences, our findings are generally consistent with those studies, which suggested Zn and Mg to be the most favorable ions for synthesizing barlowite-related compounds. More details on our defect formation energy calculations can be found in the Supplemental Material [48].

III.5 Selecting Promising Candidates

After eliminating all compounds with structural instabilities at the $\Gamma$ point, formation energies greater than 50 meV/atom above the convex hull, and Cu-O-Cu bond angles smaller than the reference compounds, 9 candidate materials remained. For these candidates, we calculated the defect formation energy $E_{d}^{f}$. To determine a final ranking, we used the following criteria:

Structural stability ($f_{0}>0$) 2. 2.

Convex hull energy (E${}_{\mathrm{hull}}<50$ meV/atom) 3. 3.

Defect energy cost ($E_{d}^{f}[\mathrm{candidate}]>E_{d}^{f}[\mathrm{ref}]$) 4. 4.

Cu-O-Cu bond angle ($\theta>\theta^{\mathrm{ref}}$)

All compounds satisfying these criteria are listed with their associated properties in Table 1. Complete data sets for all 44 herbertsmithite relatives and 66 Zn-barlowite relatives can be found in the Supplemental Material [48].

We also verified structural stability by calculating the full phonon dispersion throughout the entire Brillouin zone using the finite displacement method within the phonopy code [51]. Such calculations can identify structural instabilities associated with an enlargement of the unit cell. Dispersion curves were calculated for all candidates in Table 1. However, only one compound in the herbertsmithite family and two compounds in the Zn-barlowite family were found to be stable throughout the entire Brillouin zone. The dispersion curves of these compounds are shown in Fig. 5, while dispersions for all compounds in Table 1 can be found in the Supplemental Material [48]. Surprisingly, while Zn-claringbullite [ZnCu3(OH)6ClF] is known to have perfect kagome layers at room temperature [41], our ground state dispersion shows instabilities at the $M$ and $K$ points (see Supplemental Material [48]). The instabilities we observe in DFT may be avoided by thermal fluctuations at room temperature, which could explain the discrepancy between our calculations and the experimental results. Two other Zn-barlowite-related candidate compounds listed in Table 1, CaCu3(OH)6IBr and MgCu3(OH)6BrF, showed similar instabilities, and therefore may also be stable at room temperature (see Supplemental Material [48]).

Our calculations identify MgCu3(OH)6Br2 as a potential candidate within the herbertsmithite family, as well as CaCu3(OH)6ICl and MgCu3(OH)6ClF in the Zn-barlowite family. However, some practical considerations related to synthesis may require further investigation. For instance, the Mg analog of Zn-barlowite [MgCu3(OH)6BrF] has not been synthesized due to the insolubility of MgF2 in water. While synthesis of Zn-barlowite using NH4F yields a structurally equivalent compound, crystals obtained using this method show a similar magnetic transition to barlowite, suggesting possible differences in defect structures between the two synthesis methods [52]. The insolubility of MgF2 may therefore present difficulty in synthesizing our candidate MgCu3(OH)6ClF [41]. Synthesis of MgCu3(OH)6Br2 has been attempted, but the desired product was a Zn-barlowite analog [42]. The synthesis method, which followed the typical hydrothermal procedure, resulted in a compound with $P\bar{3}m1$ symmetry, which may mean that the herbertsmithite $R\bar{3}m$ structure is not favored in this reaction. It is possible that other synthesis methods could yield different results. To our knowledge, no experimental studies have been performed on the Ca analog of either herbertsmithite or Zn-barlowite, nor any related compounds containing I.

IV Conclusion

In summary, we performed a systematic combinatorial exploration of herbertsmithite and Zn-barlowite material relatives and identified those with properties that may enhance the likelihood of an ideal QSL ground state. We found several promising candidates—MgCu3(OH)6Br2, CaCu3(OH)6ICl, and MgCu3(OH)6ClF—that are structurally stable, thermodynamically feasible to synthesize, have high energy costs for interlayer defects, and whose structural properties may result in antiferromagnetic superexchange interactions stronger than herbertsmithite or Zn-barlowite. These compounds, if they can be synthesized, may prove to be better QSL candidates than their well-studied counterparts.

Acknowledgments

We would like to thank Siavash Karbasizadeh for helpful discussions. This work was supported by the Air Force Office of Scientific Research under AFOSR award no. FA9550-21-1-0337. C.A. acknowledges support from the UCSB Quantum Foundry Internship Program, which is funded by the National Science Foundation (NSF) through Enabling Quantum Leap: Convergent Accelerated Discovery Foundries for Quantum Materials Science, Engineering, and Information (Q-AMASE-i): Quantum Foundry at UC Santa Barbara (DMR-1906325). Use was made of computational facilities purchased with funds from the NSF (CNS-1725797) and administered by the Center for Scientific Computing (CSC). The CSC is supported by the California NanoSystems Institute and the Research Science and Engineering Center (MRSEC; NSF DMR-1720256) at UC Santa Barbara.