Counterroating incommensurate magnetic order and strong quantum fluctuations in the honeycomb layers of \rm NaNi_{2}BiO_{6}

Counterroating incommensurate magnetic order and strong quantum fluctuations in the honeycomb layers of

A. Scheie Institute for Quantum Matter and Department of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218    K. Ross Department of Physics, Colorado State University, Fort Collins, CO 80523 Quantum Materials Program, Canadian Institute for Advanced Research (CIFAR), Toronto, Ontario M5G 1Z8, Canada    E. Seibel Department of Chemistry, Princeton University, Princeton, NJ 08544    J. A. Rodriguez-Rivera NIST Center for Neutron Research, National Institute of Standards and Technology, Gaithersburg, MD 20899 Department of Materials Sciences, University of Maryland, College Park, MD 20742    J. A. Tang Department of Chemistry, Johns Hopkins University, Baltimore, MD 21218    Yi Li Institute for Quantum Matter and Department of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218    R. J. Cava Department of Chemistry, Princeton University, Princeton, NJ 08544    C. Broholm Institute for Quantum Matter and Department of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218 NIST Center for Neutron Research, National Institute of Standards and Technology, Gaithersburg, MD 20899 Department of Materials Science and Engineering, Johns Hopkins University, Baltimore, MD 21218
August 23, 2019
Abstract

We report the magnetic structure of the honeycomb antiferromagnet . We find magnetic order with moments along the axis below and then in the honeycomb plane below with a counterrotating pattern and an incommensurate ordering wave vector . Neutron scattering shows half of the expected static magnetic moment remains fluctuating at low temperatures, suggesting proximity to a quantum disordered phase. The ordering wave vector, in-plane magnetic correlations, and superexchange pathways are all consistent with a bond-dependent 120 compass model exchange in .

I Introduction

The discovery of the exactly solvable Kitaev model with a spin liquid ground state Kitaev (2006) has attracted much attention to the realization and consequences of anisotropic bond dependent exchange interactions on the honeycomb lattice. In the last decade many materials have been discovered which exhibit Kitaev Hamiltonians such as Plumb et al. (2014); Sears et al. (2015); Banerjee et al. (2016, 2017) and the iridates Ganesh et al. (2011); Biffin et al. (2014); Hwan Chun et al. (2015). Although the Kitaev model has received the most attention, there are other bond-dependent magnetic exchanges that have yet to be explored experimentally Nussinov and van den Brink (2015). An example is the honeycomb 120 compass model. This Hamiltonian has been studied by theorists for more than a decade Nussinov et al. (2004); Wu (2008); Zou et al. (2016), but experimental examples are lacking. Here we present a first experimental realization of the magnetic 120 compass model exchange on a honeycomb lattice. The resulting magnetism is incommensurate and characterized by strong quantum fluctuations. In addition to the unique magnetic order, these results demonstrate that 3d transition ions such as Ni exhibit bond-dependent exchange.

Figure 1: Crystal structure of from ref. Seibel et al. (2014).

Recently, Siebel et al. discovered and reported which features magnetic Ni ions on a honeycomb lattice Seibel et al. (2014) as shown in Figure 1. The space group is , with lattice parameters and measured at 2 K. Magnetic susceptibility and thermogravimetric analysis indicate that , leading to a ratio of about 2/3 Ni and 1/3 Ni ( and respectively in octahedral environments) due to oxygen vacancies. A Curie-Weiss fit to high temperature susceptibility data yields a Weiss temperature of K and an effective moment of 2.26(1) /Ni Seibel et al. (2014), which is consistent with the effective moment of 2.1 calculated for a 2:1 ratio of Ni/Ni. Zero-field heat capacity measurements, shown in Fig. 2, reveal two peaks, indicating second-order phase transitions at and . The strong field dependence of these peaks show these transitions are magnetic in nature. Here we report the magnetic structure and properties of based on heat capacity, electron spin resonance, and neutron scattering. We argue that the counterrotating incommensurate magnetic order and strong quantum fluctuations that we have discovered result from a first realization of the 120 compass model within the honeycomb lattices of .

Ii Experimental Results

Figure 2: Low temperature heat capacity of . (a) Plot of measured heat capacity. The inset shows two transitions at and . (b) Entropy obtained from integrating heat capacity (extrapolated to zero using a fit).

We measured the heat capacity of from K to K using a Quantum Design PPMS NIS , shown in Fig. 2. We estimated the entropy by computing (extrapolating to at using a cubic -dependence).

Figure 3: ESR data from , measured between K and K. The inset showns an example of a two Lorentzian derivative curve fit to the data. Figure 4 shows the data extracted from these fits.
Figure 4: Data extracted from two Lorentzian derivative fits to ESR data in Fig. 3. (a) Temperature dependence of the effective factor, where the red triangles represent the sharp feature (peak 1), and the blue circles represent the broad feature (peak 2). (b) Total integrated intensity of the two features. (c) Width of the resonance features.

We collected the X-band electron spin resonance (ESR) data shown in Fig. 3 on 200 mg of loose powder using a Bruker EMX spectrometer NIS . The powder was sealed in a quartz tube filled with argon gas to avoid contact with air. Magnetic field scans at temperatures between K and K were performed at 9.440 GHz, with and without sample so we can display and analyze difference data that reflect ESR from the sample. We analyzed the ESR data by fitting to two Lorentzian derivative curves: one for the Ni resonance and one for Ni, with the results shown in Fig. 4. There is a very small resonance feature right at , but we did not consider it in our analysis. The lack of temperature dependence and tiny integrated intensity (0.005(1)% of the broad resonance) suggests this feature is from contaminants in the sample chamber.

Finally, we performed a neutron scattering experiment on using MACS at the NCNR with g loose powder of anhydrous loaded in a sealed aluminum can under 1 atm helium at room temperature. The monochromator was set to doubly focusing with a pre-monochromator aperture of 360 mm x 360 mm. The data are shown in Fig. 5. We measured the momentum (Q) dependence of elastic (, ) and inelastic (, , ) scattering for temperatures between K and K. We also measured the full excitation spectrum at K (below both heat capacity peaks), K (in between the heat capacity peaks), and at K (above both heat capacity peaks). We converted the ratio of detector to monitor count rates to absolute values of the partial differential scattering cross section

by normalizing to the (001) nuclear Bragg peak in accord with Ref. Xu et al. (2013). Here cm, is the g-factor for Ni, is the magnetic form factor and is the dynamic correlation function. Empty can measurements were subtracted from the data presented in Fig. 5(a) and Fig. 5(c)-(e) with a self-shielding factor of 0.93. The horizontal line of diminished intensity at in panels (c)-(e) is is a systematic error in normalization associated with removal of the incident beam beryllium filter for (). This causes a slight offset in intensity for a small range of near the filter edge.

Iii Analysis

iii.1 Heat Capacity and Entropy

Heat capacity reveals much less entropy recovered across the phase transitions than one would expect for complete magnetic order. The expected entropy recovered for this compound is because we have a presumably quenched mixture of and . The total entropy recovered between K and K is only 46% of this total spin entropy [see Fig. 2(b)]. Even without subtracting the lattice contribution to heat capacity, we find far less than . Unfortunately no non-magnetic analogue compound is available, so we are unable to determine how much additional magnetic entropy is recovered at higher temperatures. (See Supplementary Materials for details.) Be that as it may, we can conclude that the transitions at and recover only a small portion of the magnetic entropy, which means either that significant magnetic correlations persist above K or that a significant amount of residual entropy exists within the ordered magnetic phase.

iii.2 Electron Spin Resonance

As Figures 3 and 4 show, the high temperature ESR spectrum has a sharp resonance at and a broad resonance at . The double resonance connotes separate signals from Ni and Ni ions: the sharp resonance probably from the Ni Kramers doublet, and the broad resonance probably from the Ni. ions typically have broader signals for a powder sample due to single-ion anisotropies Abragam and Bleaney (1970). Upon cooling, the sharp resonance weakens and shifts to lower field (higher effective -factor) while the broad resonance grows stronger and shifts to higher fields (lower effective -factor). The overall signal intensity follows a Curie-Weiss law [Fig. 4(b)] consistent with typical transition ion behavior Abragam and Bleaney (1970). A fit to the ESR intensity data beyond K yields K, in agreement with magnetic susceptibility measurements.

Naively, one might expect the sharp ESR signal to carry most of the intensity because 2/3 of the ions are Ni. However, the reverse is true: the broad feature carries 85-90% of the signal for K K, and more at low temperatures. This indicates exchange broadening where Ni sites broaden and shift the ESR signal from Ni neighbors. Such effects have been observed for Ni impurities in other compounds Barra et al. (1999); Glazkov et al. (2001); Chappel et al. (2000). Quantitative analysis supports this explanation: assuming a random distribution of missing O sites (with missing site O per octahedron) which cause the formation of Ni then of ions either are Ni or border on a Ni site. (This number is higher if Ni sites are correlated in any way.) This means we expect over 80% of the ESR signal to be exchange coupled to Ni impurities, which is consistent with our finding that most of the ESR signal is dramatically exchange broadened. Meanwhile, the exchange broadening of Ni surrounded by Ni neighbors is less dramatic because of the lack of single ion anisotropy.

This interpretation is consistent with the temperature dependence of the overall ESR signal: upon cooling, magnetic correlations develop and exchange coupled spin clusters around Ni impurities grow in size. As this happens, spectral weight is transferred to the broad resonance as fewer Ni can remain unaffected by the single ion anisotropy associated with Ni impurities.

This defect interpretation is also consistent with the temperature dependence of the average effective -factor. ESR studies of Ni impurities in Stoyanova et al. (1994); Barra et al. (1999) and Glazkov et al. (2001) also show the ESR resonance -factor decreasing with temperature. Such behavior can arise from anisotropic exchange terms Krug von Nidda et al. (2002), though the specific mechanism is not known for .

It is worth noting that the effective -factor for the sharp Ni resonance dramatically increases as temperature decreases below K. The reasons for this are not clear, but it could indicate the onset of short-range order Nagata and Tazuke (1972), or changes in the orbital occupation due to a low-lying orbital state.

In the end, we draw three conclusions from the ESR data. First, Ni and Ni ions are present in a 2:1 mixture with inter-site exchange interactions revealing no clustering of the two different valence states. Second, the Ni ions impact the low temperature magnetism by nucleating anisotropic correlated spin clusters. Third, changes (possibly due to magnetic correlations) in the Ni signal appear at K. This is well above the ordering temperature and consistent with the reduced spin entropy of the low phase transitions.

iii.3 Neutron Scattering

The temperature-dependent elastic neutron scattering data in Fig. 5(b) show new Bragg peaks appearing at low tmeperatures. The onset temperature matches and from heat capacity, indicating that these anomalies mark the onset of magnetic order. The inelastic scattering data in panel (a) show an increase in paramagnetic diffuse scattering above the transitions, and in particular near the 0.81 Å magnetic peak. The integrated intensity of this inelastic peak, shown in Fig. 8(d), is highest at K, and then gradually diminishes upon warming.

Figure 5: Neutron scattering of . (a) Inelastic temperature scan at , (). (b) Elastic temperature scan (, ), revealing magnetic Bragg peaks that emerge below the transition temperatures. The strong feature at 1.1 Å is remnants of a subtracted nuclear Bragg peak. (c-e) Inelastic neutron scattering data at K, K, and K.

While it may look like the intensity of the peak in inelastic scattering near the 1.1 Å nuclear Bragg peak in Fig. 5(a) is enhanced above the transition, Gaussian fits to each temperature show the integrated peak intensity is constant. The increase is actually in the Q-independent diffuse background intensity.

The fixed temperature full-spectrum scans in Fig. 5(c)-(e) provide a more complete picture of the excitation spectra. The data in Fig. 5(e) resembles powder-averaged inelastic scattering from spin waves with a bandwidth meV, which is less than the bandwidth estimated from the Curie-Weiss temperature: meV if (or meV if ). (A derivation of this equation is given in the Supplemental Materials.) This bandwidth suppression may be associated with geometrical frustration which precludes simultaneous minimization of the interactions that contribute to the Weiss temperature. In any case, the spin-wave-like excitations and the appearance of low temperature Bragg peaks show the transitions around K are to a long-range ordered magnetic state.

The higher temperature scattering at K in Fig. 5(d) shows that spin correlations persist at temperatures well above the upper phase transition. This is consistent with expectations for a frustrated quasi-two-dimensional magnet is consonant with the possibility of magnetic entropy recovered above , and with the strong temperature-dependence of the -factor in the ESR signal above .

The dynamic magnetic moment can be computed from the inelastic spectral weight per formula unit using

integrated from 0.3 meV to 2.5 meV and from 0.5 Å to 1.9 Å, where detailed balance from thermal equilibrium has been employed. We find /Ni ion at 1.8 K, 3.8(8) /Ni at 5 K, and 4.3(9) /Ni at 10 K. (Comparison to total moment estimates is made below.) These values ought to be taken cautiously because inelastic spectral weight from phonons was not excluded from the scattering that contributed to the integrals. That being said, the phonon scattering at K and at low is relatively weak (phonon intensity varies as ), making in partiuclar the result at K reliable.

Magnetic Structure:

Figure 6: Ordering wave vector of . (a) Plot of elastic neutron scattering at 2K and 10K, showing the appearance of additional Bragg peaks at low temperatures. The five strongest temperature-dependent Bragg peaks are indicated with pink vertical bars. (b) Plot of theoretically predicted magnetic peaks (green lines) from with , where varies along the axis. The horizontal yellow line shows which correctly indexes the observed Bragg peaks. The vertical pink dashed lines show smaller Bragg peaks also indexed by .

Using the elastic scattering data from , we were able to determine the magnetic structure below the phase transitions. The first step was to identify the wave vector characterizing the magnetic order. We compared the wave vectors of the five strongest temperature-dependent Bragg peaks to those calculated from . Here are nuclear Bragg peaks and is a symmetry-allowed ordering wave vector in the space group Kovalev (1965). The error bars in experimental peak locations (represented visually by the widths of the vertical bands in Fig. 6) were determined from the range of Gaussian fit peak widths for elastic data at temperatures below 4K. Visual comparisons, as in Fig. 6(b), allowed us to identify the correct magnetic wave vector . The only ordering wave vector that can account for the five strongest magnetic Bragg peaks is . No other symmetry-allowed propagation vectors besides comes close. As Fig. 6 shows, this ordering wave vector also correctly indexes weaker magnetic Bragg peaks at 1.49 Å and 1.85 Å. This wave vector means the magnetic unit cell encompasses three nuclear unit cells in the plane, and has a characteristic wave length of Å along the axis.

IRs component Ni1 Ni2 (5 K) BVs (1.8 K) BVs (1.8 K) BVs
Real (1.5 0 0) (0 -1.5 0) 14.1 9.6 0.202
Imaginary ( 0) ( 0)
Real (1.5 0 0) (0 1.5 0) 14.1 9.7 0.183
Imaginary ( - 0) (- 0)
Real (1.5 0 0) (0 -1.5 0) 5.8 0.0
Imaginary ( 0) (- 0)
Real (0 0 3) (0 0 -3) 0.314 0.366 0.337
Imaginary (0 0 0) (0 0 0)
Table 1: Irreducible Representations and associated basis vectors (BVs) for space group and propagation vector . The are from a Rietveld refinement to the 1.8K data using the FullProf suite. Refinements to the 1.8 K data cannot distinguish between and , but symmetry considerations forbid and from causing the transition at .
Figure 7: Magnetic refinement of . (a) Refinement of K elastic neutron scattering data with high temperature data (K-K) subtracted. The slight increase with Q in the background difference intensity level away from magnetic Bragg peaks can result from the change in the Debye Waller factor with T. This produces a dependence of the difference intensity for low . The resulting magnetic structure, shown in panel (c), has an in-plane component to the spins (). (b) Refinement of K neutron data with high temperature data subtracted. The magnetic structure, shown in panel (d), has all spins aligned along the axis (). Panels (e) and (f) show the temperature dependence of refined irrep weights of in combination with and , respectively. The red vertical lines indicate and . In both cases, (-axis magnetism) is associated with , and and (in-plane magnetism) with . Error bars indicate one standard deviation.

The next step in determining the magnetic structure was fitting the neutron scattering intensity data to symmetry allowed structures with the given magnetic wave vector via Rietveld refinement. We used group-theoretical analysis to generate the irreducible representations ("irreps") of the little space group, which are shown in Table S.IV. These irreps were computed by hand via the method outlined by Ref. Izyumov and Naish (1979) (see Supplemental Materials for these calculations), and were cross-checked with the program SARAh Wills (2000). There were originally four basis vectors in the two-dimensional irrep treating the two magnetic ions in the unit cell separately. The basis vectors of were combined so as to preserve the equivalency of the two Ni sites. This site-equivalency is necessary to permit a second transition at (see Supplemental Materials for details). We refined the elastic scattering data at K (below the first transition) and at K (below the second transition) using the Fullprof software package Rodríguez-Carvajal (1993), with high temperature data (K - K) subtracted to isolate the temperature-dependent Bragg peaks. The space groups and their respective best fit values are listed in Table S.IV, and the refinements are shown in Fig. 7.

While chemical analysis and ESR clearly show that Ni defects are present, they are specifically not included in the refined model. If they indeed occupy random sites in the lattice, the associated magnetic disorder is "Fourier-filtered" out of the magnetic neutron Bragg diffraction pattern and shows up as a coherent diffuse background. Such diffuse scattering can be difficult to detect and distinguish from residual inelastic scattering in powder diffraction data. Accordingly, we carried out the refinements assuming only one type of magnetic ion, and the resulting model fits the data quite well.

In refining the magnetic structure at K, we used just one irrep at a time because only one second-order phase transition has been passed at K. yielded the best fit. For the K data we fit to combinations of (the K irrep) with and and found both combinations fit the K data almost equally well (right two columns in Table S.IV). To test this two-stage order, we repeated the refinements allowing multiple irreps at all temperatures. As Fig. 7(e)-(f) show, the relative weights of and refine to zero above , meaning that only is present for .

The refined magnetic structure between and [Fig. 7(d)] has all spins aligned along the axis, with the moment size modulated versus displacements within the basal plane and along the axis. This means that some spins fluctuate more than others within the ordered phase. In the magnetic structure below [Fig. 7(c)] every spin gains a counterrotating plane component (where the two Ni spins in the unit cell rotate opposite directions versus displacement) while the amplitude of the -axis component continues to increase upon cooling. Thus we conclude that is associated with ordering the -component of spins while the in-plane spin components only order for .

Although neutron diffraction cannot distinguish in-plane spin structures based on and , symmetry analysis identifies the one based on as the correct low temperature structure. This is because the addition of would not reduce the symmetry of the system, and could not result in a phase transition at . Meanwhile, breaks a mirror-plane plane that is present in the structure so its appearance must be associated with a phase transition (see Supplemental Materials for details). Therefore, we can identify as the proper in-plane magnetic structure. has ferromagnetic in-plane bond-dependent correlations (see Fig. 9). Although the magnetic structure breaks inversion symmetry, the counterrotation precludes a definite handedness as seen in spiral incommensurate ferroelectrics Lawes et al. (2005), so we do not expect ferroelectricity in this compound.

The peak widths in the refined model in Fig. 7 were defined by the nuclear peak refinement (see Supplementary Information), but the magnetic Bragg peaks are slightly wider than the peaks from the refined model. This indicates the magnetic correlation length is less than the correlation length of the nuclear structure. We can quantify this by fitting the 0.81 Å peak with a convolution of a Gaussian (with peak width defined by the nuclear phase) and a Lorentzian profile, where the inverse of the Lorentzian HWHM is the magnetic correlation length. Using this method, we calculate a magnetic correlation length of Å. (See Supplemental Material for details.) It is noteworthy that the spin correlations apparently extend well beyond the correlation length anticipated for oxygen vacancies and Ni2+/Ni3+ disorder.

At K, the refined ordered moments have a fixed in-plane magnitude while their axis component is spatially modulated [see Fig. 7(c)]. The overall size of the ordered moments range from 1.43 /Ni to 0.32 /Ni, with a mean value of /Ni. These values are taken from refinements which allow the magnetic peak width to be larger than the nuclear peak width so that all the elastic magnetic diffraction is accounted for. Adding this to the fluctuating moment from the inelastic sum-rule analysis above, we find that the total magnetic neutron scattering corresponds to a mean squared moment of , which is within error bars of inferred from high- susceptibility data through Curie-Weiss analysis. We also find that 79(4)% of the magnetism remains dynamic within the ordered phase (K).

Theoretically, the neutron spectral weight from elastic magnetic scattering is proportional to and the total magnetic scattering is proportional to Squires (1978), so that the ratio for elastic vs total magnetic spectral weight for a fully static spin configuration is . So theoretically, in fully-ordered % of the moment should be elastic. Our measured ratio % means that we observe only 54(8)% of the maximum static moment in the ordered phase. This evidences strong magnetic fluctuations in the ground state. We note that the observed elastic-inelastic spectral weight ratio does match the calculated ratio if order resides on the only and the is completely dynamic: %. However, powder neutron diffraction cannot distinguish this situation from collective quantum fluctuations.

Iv Discussion

has a larger magnetic unit cell and a more complex magnetic ground state than related Ni honeycomb compounds Seibel et al. (2013); Zvereva et al. (2015). As we shall now show, the in-plane spin structure is particularly interesting, and points to bond dependent magnetic interactions.

Two-step longitudinal to transverse polarized incommensurate magnetic ordering has been seen in other materials with easy axis anisotropy and competing interactions such as Kenzelmann et al. (2005) and Kenzelmann et al. (2006) in a mechanism described theoretically by Nagamiya Nagamiya (1967). Association of the transitions in with this mechanism is supported by reports of a Ni easy-axis anisotropy along in the honeycomb compounds and Zvereva et al. (2015), which have similar Ni ligand environments to . With an easy-axis anisotropy, one would expect low energy structure in the spin-wave spectrum at energy transfer of K, or meV. However, our neutron experiment does not resolve the spectrum below meV, so we could not detect such structure.

One puzzling aspect of the magnetic order is the temperature-dependent elastic scattering [see Fig. 8(a)-8(c)], which does not follow the typical single-exponent order-parameter curve for a second order transition. The temperature dependent magnetic Bragg diffraction intensity is similar for the magnetic Bragg peaks: it increases linearly as temperature decreases between and , and then flattens off and decreases slightly at the lowest temperatures. This low-temperature decrease in elastic intensity is accompanied by an enhancement of inelastic fluctuations, revealed by the small upturn in Fig. 8(d). This indicates a weakening of the counter-rotating spin order as might occur near a transition to a different phase. We leave this feature to be explored in future studies.

The observed ordering wave vector is very unusual for honeycomb compounds; in fact unprecedented to our knowledge. The in-plane wave vector is difficult to stabilize, and suggests a highly frustrated set of exchange interactions. honeycomb order is found in phase diagrams of isotropic exchanges only in the "spiral phase" when exactly Li et al. (2012). We consider this possibility unlikely because (i) it is stabilized in a vanishingly small region of parameter space, and (ii) the spiral phase in-plane structure is corotating, and does not match the counterrotating in-plane structure. As we will explain shortly, the most natural explanation for the structure is a bond-dependent exchange interaction.

The incommensurate modulation along the axis requires competing interactions along . High -axis frustration is likely in because the exchange pathways for the first, second, and third nearest inter-plane neighbors are Ni-O-Na-O-Ni, so they all may have similar strengths. However, mean-field level Luttinger-Tisza analysis Lyons and Kaplan (1960) with the nearest five neighbor exchanges failed to identify a set of inter-plane exchange parameters with a finite extent in parameter space which stabilize incommensurate -axis order. (See Supplemental Materials for details.) This suggests the incommensurate order is stabilized by quantum effects or an interplay with non-trivial in-plane interactions. We also note that in the hydrate version of (with molecules in-between the planes) the transition temperature from heat capacity (K) barely changes at all in comparison to the anhydrate Seibel et al. (2014), which indicates that inter-plane interactions are not very significant and the ordered magnetism is quasi-2D. According to the Mermin-Wagner theorem, magnetic order in such a 2D system requires anisotropic in-plane interactions.

Figure 8: Temperature dependence of magnetic peaks, extracted from Gaussian fits. (a)-(c) Order parameter curves for three elastic magnetic peaks. (d) and (e) respectively show the area and FWHM of the E=0.5 meV inelastic peak at 0.815 Å. (d) has the same units as (a-c). Red vertical lines indicate and from heat capacity. Error bars indicate one standard deviation.

Perhaps the most intriguing aspect of the magnetic order is the counterrotating in-plane structure, shown in Fig. 9. This structure is unusual because the the mean field component of isotropic exchange interactions average to zero for such structures. Specifically, for nearest neighbor, next-nearest neighbor, and all further neighbor spin pairs forming a (1/3, 1/3) counterrotating spin state on the honeycomb lattice. This can be proved as follows: Fig. 9 shows the angles between nearest-neighbor spins are always 180  60  and 300. Thus, for nearest neighbor exchange on any site, . Extending this analysis to further neighbors is straightforward and yields the same result. (This result holds for other layers where the spins are rotated about the axis as shown in the lightly-shaded structures in Fig. 9.) This means that this magnetic structure cannot be stabilized by isotropic exchange interactions at the mean-field level. This condition holds for each bond even if the three-fold axis is broken and the three bond directions have different interaction strengths, as in Zvereva et al. (2015), because for each of the three distinct bond directions considered as groups. Confirming this conclusion is the fact that this structure is not found in theoretical phase diagrams for isotropic exchange interactions on the ideal honeycomb lattice Albuquerque et al. (2011); Li et al. (2012); Clark et al. (2011). Dzyaloshinskii-Moriya (DM) interactions are forbidden on the Honeycomb lattice because the midpoint between magnetic ions is a point of inversion Moriya (1960). This leaves two possibilities: either Ni disorder influences the magnetic interactions in such a way as to stabilize this structure (DM interactions are allowed on bonds with oxygen vacancies), or there must be more exotic anisotropic interactions at play.

Figure 9: In-plane magnetic structure of below K. The lightly-shaded structures in the background show how the spins rotate in the plane as one goes down the axis.

Certain anisotropic exchange interactions are possible through bond-dependent orbital interactions. In , the magnetic structure is consistent with a type of exchange called a 120 compass model exchange Nussinov and van den Brink (2015); Mostovoy and Khomskii (2002) (here referred to as the 120 exchange). This interaction is analogous to the the Kitaev bond-dependent interaction, but the Ising-like exchange directions are coplanar and 120 apart. This model is ordinarily studied in the context of orbital order, but it is symmetry-allowed for magnetic interactions as well. The general exchange Hamiltonian can be written

where are the three distinct nearest neighbor bonds, and are projections of an equivalent spin operator (usually denoted ) along three directions in the plane: , , and Nussinov and van den Brink (2015). In the general 120 exchange model, the Ising-like exchange directions can point along any direction in the plane (so long as they stay 120 apart), but in this case the correlations match a 120 Hamiltonian with exchange direction perpendicular to the bonds. The net effect in the honeycomb lattice is that each bond energy only depends on the spin components perpendicular to the bond and in the honeycomb plane.

This type of symmetry-allowed interaction stabilizes the peculiar in-plane structure: for every spin pair, the in-plane spin components that are perpendicular to the bond direction are parallel (see Fig. 9). Thus, a ferromagnetic 120 exchange interaction stabilizes the structure.

Theoretical studies have examined the 120 exchange on the honeycomb lattice, and generally find that the ground state is in-plane order stabilized via order-by-disorder Wu (2008); Zhao and Liu (2008); Nasu et al. (2008). Intriguingly, spin wave analysis of the antiferromagnetic 120 exchange on the honeycomb lattice with exchange perpendicular to bonds predicts the exact counterrotating in-plane structure that accounts for our diffraction data Wu (2008) (with the difference that ref. Wu (2008) assumed antiferromagnetic 120 exchange whereas we have ferromagnetic 120 exchange). The theory on this, however, is not settled; others have predicted a dimer-like ground state for 120 exchange on the honeycomb lattice Zou et al. (2016). The ground-state dispute notwithstanding, the in-plane order and the evidence for bond-dependent correlations indicates a 120 exchange model is relevant for . In addition, the large fluctuating moment observed even within the ordered phase indicates proximity to quantum criticality or a quantum-spin-liquid phase.

Theoretical studies of the honeycomb Heisenberg-Kitaev model also predict a in-plane order dubbed the "120 phase" Rau et al. (2014). Intriguingly, the Kitaev 120 phase spin-structure is that of the basis vector, which fits the neutron diffraction data equally well as . If there were just one phase transition we would be required to affirm the Heisenberg-Kitaev 120 structure, but the existence of a double phase transition disallow this spin structure, which leads us to favor the 120 exchange model over the Heisenberg-Kitaev model as the explanation.

A similar incommensurate counterrotating magnetic order was observed in honeycomb with ordering wave vector Williams et al. (2016), or expressed in the reciprocal lattice. In this case, the spin structure is attributed to a Kitaev-like Hamiltonian with different couplings on the vertical and zig-zag bond directions Biffin et al. (2014); Williams et al. (2016); Kimchi and Coldea (2016). Although and share a 1/3 counterrotating structure, there are important differences between them. First, has two magnetic phase transitions and has one. Second, the counterrotating structures are different and the spin structure is inconsistent with theoretical predictions from the 120 model. Third, does not have three-fold rotation symmetry about its magnetic sites, and its structure requires either or an additional Ising term on the bonds to stabilize the counterrotating order Williams et al. (2016). Meanwhile, can be explained by a Hamiltonian which preserves the three-fold axis.

Figure 10: Ni-O-Ni exchange pathways of . The Ni-O-Ni bond angle is 97.9(4)  leading to an anisotropic exchange perpendicular to the ligand-ion plane (shown by the red, green, and blue arrows) plus an isotropic term.

The microscopic origin of the 120 exchange in cannot be determined from the data reported here, but we present one possibility: In the case of large spin orbit coupling and a 90 ion-ligand-ion bond, the nearest neighbor exchange is an Ising-like anisotropic exchange in the direction perpendicular to the plane formed by the superexchange pathway Jackeli and Khaliullin (2009). This effect emerges also for ions with intermediate spin-orbit coupling, such as Ru in Kim et al. (2015); Banerjee et al. (2016), due to direct overlap of orbitals Kim et al. (2015). In , the situation is imperfect with a 97.9(4) Ni-O-Ni bond (shown in Fig. 10) and a Ni spin-orbit coupling constant meV Abragam and Bleaney (1970) (26% of meV Banerjee et al. (2016), and 9.5% of mev Jackeli and Khaliullin (2009)). While these features would weaken the anisotropic exchange, the anisotropy mechanisms (particularly the orbital overlap) are still present in Ni—which means that the nearest neighbor exchange can be written ( and denote the directions in and perpendicular to the the Ni-O-Ni plane) and . The resulting exchange anisotropies, shown in Fig. 10, do not lie in the plane, but are rotated 38.5 out of the plane. So, rather than 120  the the anisotropy directions actually are 94.3(5) apart ( in the nomenclature of ref. Zou et al. (2016)), making this exchange very close to the celebrated Kitaev model where . According to recent theory work Zou et al. (2016), this lies within the regime where the system behaves according to 120 exchange, but extremely close to the critical point of the Kitaev spin liquid regime (). This proximity to the Kitaev spin-liquid quantum critical point may explain the large fluctuating moment within the ordered phase.

We note that because the normal directions to planes in Fig. 10 have a component along these interactions could also be responsible for the -polarized phase for . An extension of the Nagamiya phase diagram to the 120 interaction on the honeycomb lattice is needed to ascertain this.

One caveat to all this is that we have ignored the chemical disorder from oxygen-site deficiency. The ESR data shows that Ni sites impact more than 90% of the magnetic sublattice at low temperature so this disorder may play a significant role in the magnetism of . Nevertheless, we consider it unlikely that random disorder could alone induce the spatially coherent counter-rotating spin order that we observe. Meanwhile, the bond-dependent counterrotating structure that we observe for is favored by the 120  compass model at the mean field level.

V Conclusion

We have acquired and analyzed heat capacity, ESR, and neutron scattering data on , and have shown this compound orders magnetically in two stages: axis spin components order below and the plane components order below . ESR data indicate Ni impurities that are projected to exists at the 30% level through chemical analysis influence the low-temperature magnetic properties. Half of the expected magnetic moment remains fluctuating in the low temperature magnetic phase. Specific heat, ESR, and neutron scattering also indicate substantial short ranged magnetic correlations in the paramagnetic phase.

The observed magnetic structure has unusual counterrotating in-plane correlations, which are not favored by isotropic interactions. These correlations and the special ligand environment suggest a bond-dependent 120 compass model exchange. The in-plane structure predicted through a spin wave analysis of the 120 exchange on the honeycomb lattice Wu (2008) is the only symmetry allowed structure that is consistent with our diffraction data. This constitutes the only experimental evidence that we are aware of for magnetic 120 exchange on the honeycomb lattice.

These results are significant firstly because new types of magnetic exchange offer the chance to study new physics. Bond dependent interactions beyond the Kitaev model have not been well-studied, but they give rise to a wealth of model systems to explore. Secondly, bond-dependence in Ni is surprising; conventional wisdom says its weaker spin-orbit coupling would make bond-dependent effects too weak to impact magnetism Jackeli and Khaliullin (2009). But in , the effects seem to be dominant, which raises the possibility for discovering Kitaev-like spin-liquid phases in 3d transition metal oxides with edge sharing six-fold coordination.

Acknowledgments

This work was supported through the Institute for Quantum Matter at Johns Hopkins University, by the U.S. Department of Energy, Division of Basic Energy Sciences, Grant DE-FG02-08ER46544. AS and CB were supported through the Gordon and Betty Moore foundation under the EPIQS program GBMF4532. Access to MACS was provided by the Center for High Resolution Neutron Scattering, a partnership between the National Institute of Standards and Technology and the National Science Foundation under Agreement No. DMR-1508249. AS acknowledges helpful discussions with Kemp Plumb.

References

Supplemental Material

s.i Magnetic Entropy from Phonon Subtraction

No nonmagnetic analogue to is currently available to measure the phonon specific heat and isolate the magnetic contribution to heat capacity in . Therefore, we attempted to estimate the magnetic entropy by subtracting a phonon background calculated using the Debye equation for heat capacity Ashcroft and Mermin (1976). Here and were fitted using the ten highest temperature data points, which gave values of per unit cell and K. The results are shown in Fig. S1, and indicate that between K and K the entropy only reaches 65% of the expected entropy [see Fig. S1(b)].

In the Dulong-Petit limit should be 5 (the number of atoms per Ni). Our fitted value is 2.90(6). This discrepancy is a sign that the Debye estimate for heat capacity is unrealistic. Therefore, we do not have much confidence in the entropy computed from this background subtraction, and leave the presence of high temperature magnetic entropy as an open question.

s.ii Nuclear Refinements

In addition to the neutron experiment on the MACS spectrometer, we used the same powder sample to collect neutron diffraction data using the BT1 powder diffractometer at the NCNR. We used 18.9 meV neutrons with 60’ collimation before the monochromator and 20’ after the sample, measuring for 8.5 hours at 1.5 K, 6 hours at 4.8 K, and 6 hours at 28 K. These measurements cover a much larger range than the MACS measurements with better Q-resolution, and allows a more complete determination of the nuclear structure.

Before refining the magnetic structure, we refined the nuclear structure using both the MACS and BT1 neutron diffraction data sets. The refinements are shown in Fig. S2. Both these data sets were taken below , and thus the refinements include the magnetic phase. The refinement in panel (a) includes a nuclear phase with Bragg peaks located as indicated by the upper vertical green lines, and a magnetic phase indicated by the lower row of vertical green lines. The refinement to the BT1 data set in panel (b) includes the nuclear phase (topmost vertical green lines), an additional NiO powder phase with 1.5% of the refined intensity of the nuclear phase (second row of green lines), aluminum peaks from the sample can (third row of green lines), and the magnetic phase (fourth row of green lines). The BT1 data do not show the magnetic structure as clearly as the MACS data—only the peak at 0.81 Å is visible at 1.5 K—but the BT1 data include many more nuclear peaks and provide a better view of the nuclear structure.

Figure S1: Low temperature heat capacity of . (a) Plot of measured heat capacity (black) and calculated phonon background (green). The blue data show the measured heat capacity minus the calculated background , which we interpret as the magnetic contribution to heat capacity. The inset shows two transitions at and . (b) Entropy obtained from integrating heat capacity (extrapolated to zero using a fit).
atom type label S.O.F.
Na Na1 1/3 2/3 1/2 1
Ni Ni1 1/3 2/3 0 1
Bi Bi1 0 0 0 0.912
Bi Bi2 0 0 0.114 0.080
O O1 0.344(4) 1.0 0.180(1) 0.944
Table S.I: Refined nuclear positions for , .
Figure S2: Nuclear refinement of . (a) Refinement of unsubtracted MACS data at K, which was used to define the peak widths and intensities of the magnetic peaks in the magnetic refinement. The vertical green lines show the peak locations from the nuclear and magnetic phases. Features for arise from instrumental direct beam backgrounds. (b) Refinement of BT1 data, showing fits to much higher peaks. This refinement includes the nuclear phase, a NiO powder phase, aluminum from the sample can, and the magnetic phase.

The refined nuclear model, given in Table S.I, displayed in Fig. 1 of the main text, and described in detail in Ref Seibel et al. (2014), provides a reasonable fit. However, some small peaks are not accounted for, most noticeably a weak Bragg peak at Å. It is unclear what causes these deviations, whether there exists a nuclear supercell associated with oxygen vacancies or an additional phase in the sample. Be that as it may, none of the unindexed peaks are temperature-dependent, which means that the magnetic signal from temperature subtraction is reliably from alone. This signal could be indexed by a single ordering wave vector and fit to a consistent model based on the proposed chemical structure.

s.iii Relating Neutron Bandwidth to Curie Temperature

Based on a spin Hamiltonian where represents bond energies, the Curie temperature (the temperature at which spontaneous magnetization occurs in the mean field approximation) is

(S.1)

Ashcroft and Mermin (1976). This is the same Curie temperature which appears in the Curie-Weiss law .

Meanwhile, the expression for a spin wave dispersion for collinear antiferromagnetic order is (see for example Lovesey eq. 9.245 Lovesey (1984)), where is the coordination number, is the exchange interaction, is the single ion anisotropy, and where lists the nearest neighbors. is maximal when is minimal, which in the honeycomb lattice goes to zero when , being the nearest neighbor distance. Thus, is maximal at

(S.2)

By combining eq. S.1 with S.2, we can solve for and relate the curie temperature to the neutron scattering bandwidth:

(S.3)

Since this result was derived for a very simple AFM we should only expect the general form of it to hold true more generally. Frustrated interactions will contribute to , but not towards increasing the spin wave bandwidth. Thus we can expect a reduced bandwidth relative to S.3 for frustrated magnets.

s.iv Magnetic Refinements

Figure S3: Comparison between the magnetic refinements using (grey line) and (orange line). The differences between the two structures are too subtle to differentiate using the neutron diffraction data from our experiments.

s.iv.1 Irrep Decomposition

Here we summarize our analysis to generate basis vectors of the space group with the ordering vector .

Space group (also written in Schoenflie notation) has 12 point symmetry operations. Half of them preserve up to a reciprocal lattice vector, which is an ordering vector of type in Kovalev’s notation, yielding a group of the propagation wave vector with the following point operations:

where the unit vectors (100), (010), and (001) are along the , , and axes respectively. Generating the permutation, axial, and magnetic representations yields a character table in Table S.II. In single valued representations, there are three irreducible representations listed in Kovalev’s tables Kovalev (1965), shown in Table S.III.

3 0 0 -1 -1 -1
2 2 2 0 0 0
6 0 0 0 0 0
Table S.II: Character table for with propagation vector .
1 1 1 1 1 1
1 1 1 -1 -1 -1
Table S.III: Irreducible representations for with propagation vector . .

Given their dimensionality, and have one basis vector and has two basis vectors. To find the basis vectors, we project onto the test functions: , , . Using the projection equation Izyumov and Naish (1979)

(S.4)

we have, throwing away all the zero pairs of basis vectors, the set of basis vectors listed in Table S.IV.

As is immediately clear, this procedure yields more than two basis vectors for which therefore cannot be orthogonal. Thus, we must combine them into two pairs of linear combintions. Two of the basis vectors describe one triangular Ni sub-lattice site and two describe the other equivalent Ni lattice, the two together forming the honeycomb structure. The sets of basis vectors are identical except for a sign change, but they describe the two lattice sites separately. Linear combinations that link the two sublattices take the form (where ), since the diffraction pattern is independent of . The existence of another phase transition, however, requires that site equivalency be enforced. A value of would result in symmetry, but any value other than would result in symmetry (only the identity operation)—which has no symmetry elements beyond translations. This precludes a second order phase transition to a lower-symmetry state that does not modify the magnetic wave vector. We do have an additional phase transition at , so a point group symmetry must remain for . Therefore, we neglect the possibility of spontaneous sublattice symmetry breaking and set . This results in the basis vectors listed in Table I of the main text. The in-plane 120exchange bond-dependent correlations are described by the one-dimensional irrep which is not subject to these considerations.

IRs component Ni1 Ni2
Real (1.5 0 0) (0 -1.5 0)
Imaginary ( 0) ( 0)
Real (1.5 0 0) (0 1.5 0)
Imaginary ( - 0) (- 0)
Real (0 0 0) (0 -1.5 0)
Imaginary (0 0 0) (- 0)
Real (0 0 3) (0 0 0)
Imaginary (0 0 0) (0 0 0)
Real (1.5 0 0) (0 0 0)
Imaginary ( 0) (0 0 0)
Real (0 0 0) (0 0 -3)
Imaginary (0 0 0) (0 0 0)
Table S.IV: Original irreducible representations and associated basis vectors for space group and ordering wave vector .

s.iv.2 In-plane Structure and Symmetry

As noted in the text, there are two combinations of irreducible representations that fit the low temperature phase of : (with antiferromagnetic 120exchange in-plane correlations) and (with ferromagnetic 120exchange in-plane correlations). The two different predicted diffraction patterns are shown in Fig. S3, with peak widths fit to the magnetic peaks. There are subtle differences between the patterns, but the differences are so small that we are unable to distinguish between them with neutron diffraction. Therefore, we look to symmetry considerations to determine the correct ground state.

Second order phase transitions are directly associated with symmetry breaking. This means that for in-plane spin ordering to account for the phase transition that we observe at , the in-plane spin structure must break a symmetry operation of the intermediate temperature phase. This allows us to identify the in-plane spin order for .

The structure of the intermediate temperature () phase, when site-equivalency is enforced, has only two valid symmetry operations from the group of the propagation vector : (identity) and (reflection about [110]). Meanwhile, preserves all symmetries (, , , , , and ), preserves all the 3-fold rotation symmetries but no mirror planes (, , and ), and the in-plane structure preserves only the identity and the (110) mirror plane ( and ). A structure would have and symmetry, resulting in no broken symmetries. A structure would have only symmetry, resulting in a broken symmetry. A structure would have and symmetry, resulting in no broken symmetries. The only in-plane structure that breaks a symmetry of the intermediate temperature phase is . This means the magnetic structure of must form a reducible representation of based on irreps . The corresponding spin structure is depicted in Fig. 7(c) of the main text.

s.iv.3 Correlation Length

The peak widths of the nuclear Bragg peaks in Fig. S2 are smaller than the magnetic Bragg peaks widths in the temperature-subtracted data. This indicates the magnetic correlation length is less than that of the underlying crystal structure. To determine the magnetic correlation length, we fit the the strongest magnetic Bragg peak (at 0.81 Å) with a convolution of a Gaussian (with peak width defined by the nuclear phase) and a Lorentzian profile, as shown in Fig. S4. The inverse of the Lorentzian HWHM is the magnetic correlation length, which has a best fit value Å, indicated by the minimum in reduced of the convoluted profile fit in Fig. S4.

Figure S4: Profile fits to the 0.81 Å magnetic Bragg peak, with the data in red, the original Gaussian fit with the with the width defined by the nuclear peaks in black, and the best fit convoluted profile in blue. (a) shows fits to MACS data, and (b) shows fits to BT1 data. The insets show reduced of a convoluted profile fit vs magnetic correlation length.

The 0.81 Å peak is at peak, where is the magnetic propagation vector . The remainder of the magnetic Bragg peaks are much weaker and the fits are consequently less reliable, but the results are consistent: fitting the MACS data for 1.22 Å [], 1.61 Å [], and 2.12 Å [] peaks simultaneously yielded a correlation length of Å, which agrees with the fit of the 0.81 Å peak to within uncertainty. These peaks are 27, 85, and 86  from the axis respectively, which means the first is associated mostly with axis correlations and the last two are associated with in-plane correlations. If we treat the peak (mostly along the axis) separate from and (mostly in-plane), we find a correlation length of Å  for and a correlation length of Å  for the and peaks. This indicates a correlation length three times smaller along the -axis than in the plane, consistent with quasi-2D magnetic order.

We carried out the same analysis on the 0.81 Å peak from the BT1 data in Fig. S4(b), and found a correlation length of Å. We consider this value to be more reliable than the MACS data because the BT1 nuclear peak width is defined by many peaks and is well constrained, but the nuclear peak width for the MACS data is defined only by three peaks (see Fig. 7) and is underconstrained.

s.v Luttinger Tisza Analysis

As noted in the text, the magnetic ordering wave vector of is unusual on the honeycomb lattice, both because of the axis incommensurability and the in-plane order. To explore whether such a magnetic ordering wave vector can be stabilized by Heisenberg interactions at the mean-field level, we used Luttinger-Tisza theory Lyons and Kaplan (1960). While this methods has its limitations and does not account emergent interactions resulting from thermal or quantum fluctuations, it does give a basic picture of what orders are readily stabilized.

We began with the Fourier transform of the Heisenberg Hamiltonian for helical order:

(S.5)

where and sum over sites within the paramagnetic unit cell and

(S.6)

With two nickel sites per unit cell, the Hamiltonian can be written as a matrix whose smallest eigenvalue is the minimum energy. Although is complex, when summing over all atoms in the unit cell the eigenvalues of the matrix are always real. With this equation, one can find the ordering wave vector which minimizes for a given set of —i.e., we identify the the magnetic wave vector stabilized by a given set of interactions.

To search for a set of exchange constants which stabilize ( signifies incommensurability), we systematically defined a series of exchange constants and found the minimizing . We used a L-BFGS-B minimization routine Morales and Nocedal (2011), always with as the starting .

Figure S5: Exchange constants considered in the Luttinger Tisza analysis of . and are in-plane interactions while , , and are out-of-plane interactions.

In our analysis we considered five exchange interactions: two in-plane, and three out-of plane interactions (Fig. S5). The super-exchange paths of all three out of plane interactions involve the same number of atoms: Ni-O-Na-O-Ni (Fig. S5), which indicates their strength could be comparable. Including both Ni sites in the Hamiltonian ensured that is always real. We set (nearest neighbor exchange), and let the other interactions vary from -1 to 1, and from -2 to 2.

As a rough cross-check, we compared our Luttinger-Tisza results to more sophisticated calculations of the honeycomb phase diagram Albuquerque et al. (2011); Li et al. (2012). In our calculations the transition from (0,0) to (1/2, 1/2) occurs at when all other interactions are zero. In refs. Albuquerque et al. (2011); Li et al. (2012), the transition is closer to , though there is an intermediate disordered phase in between that does not appear in the Luttinger-Tisza analysis.

Figure S6: Out-of-plane exchange pathways in . The first (green), second (purple), and third (blue) nearest inter-honeycomb-plane exchange pathways all follow a Ni-O-Na-O-Ni pathway, which means they could be of comparable strength.

A selection of results of this analysis are shown in Figs. S7 and S8. Most of the exchange parameter space considered stabilizes commensurate order, but never in-plane order. The boundaries between phases [for example between and ] sometimes show incommensurability, but "zooming in" and increasing the resolution of the parameter search shows the incommensurate regions exist only on the boundaries. No finite regions of parameter space stabilize order. The failure to account for the observed magnetic order suggests that the exchange Hamiltonian is not of the conventional isotropic form and this concords with the analysis presented in the main text that the counterrotating state is not favored by conventional bond independent exchange interactions.

Figure S7: Ordering wave vectors stabilized by isotropic exchange interactions , , , and (defined in Fig. S5) relative to , as calculated by Luttinger Tisza theory. Each panel shows a range of and values for specific values of and . The ordering vector is indicated by the colorscale to the right. refers to an incommensurate order, which only appears on the boundaries between other phases.
Figure S8: 3D plot of ordering wave vectors stabilized by isotropic exchange interactions , , and relative to with , as calculated by Luttinger Tisza theory. (a) and (b) show two views of the same data set. The ordering vector is indicated by the colorscale to the right. Note that incommensuarate order (ordering vectors with ) only appears on the boundaries between other phases.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
""
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
   
Add comment
Cancel
Loading ...
305088
This is a comment super asjknd jkasnjk adsnkj
Upvote
Downvote
""
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters
Submit
Cancel