Noncollinear magnetic order, in-plane anisotropy, and magnetoelectric coupling in a pyroelectric honeycomb antiferromagnet Ni$_{2}$Mo$_{3}$O$_{8}$

Ni$_{2}$Mo$_{3}$O$_{8}$ is a pyroelectric honeycomb antiferromagnet exhibiting peculiar changes of its electric polarization at magnetic transitions. Ni$_{2}$Mo$_{3}$O$_{8}$ stands out from the isostructural magnetic compounds, showing an anomalously low magnetic transition temperature and unique magnetic anisotropy. We determine the magnetic structure of Ni$_{2}$Mo$_{3}$O$_{8}$ utilizing high-resolution powder and single-crystal neutron diffraction. A noncollinear stripy antiferromagnetic order is found in the honeycomb planes. The magnetic space group is \textit{P$_C$na}2$_1$. The in-plane magnetic connection is of the stripy type both for the $ab$-plane and the $c$-axis spin components. This is a simpler connection than the one proposed previously. The ferromagnetic interlayer order of the $c$-axis spin components in our model is also distinct. The magnetic anisotropy of Ni$_{2}$Mo$_{3}$O$_{8}$ is characterized by orientation-dependent magnetic susceptibility measurements on a single crystal, consistent with neutron diffraction analysis. The local magnetoelectric tensor analysis using our magnetic models provides new insights into its magnetoelectric coupling and polarization. Thus, our results deliver essential information for understanding both the unusual magnetoelectric properties of Ni$_{2}$Mo$_{3}$O$_{8}$ and the prospects for observation of exotic nonreciprocal, Hall, and magnonic effects characteristic to this compound family.


I. INTRODUCTION
Multiferroic materials are compounds in which electric polarization and magnetic order coexist in the same phase.The magnetic and dielectric properties of multiferroics can often be manipulated by the conjugate fields.In particular, their magnetic order may be affected by an applied electric field, and the electric polarization by a magnetic field.Such effects are called magnetoelectric (ME).Compounds exhibiting large ME effects are highly sought because of fundamental and technological interest.
Pyroelectric magnets are a subset of polar multiferroics with fixed-direction polarization.Because such a built-in polarization often has a significant magnitude, it holds the potential for large magnetically-induced changes.Many pyroelectrics can be grown as monodomain crystals, and therefore, they do not require complex multifield poling procedures to induce a single-domain multiferroic state.These properties are important from both the technological and fundamental perspectives because they allow simplified manipulations, as well as definitive studies, of the ME effects with a potentially large magnitude.
A family of polar magnetoelectric A 2 Mo 3 O 8 (A = Fe, Mn, Co, Ni) [1,2] compounds with widely varying ME properties has recently attracted significant interest.Fe 2 Mo 3 O 8 and its derivative Fe 2−x Zn x Mo 3 O 8 have been studied the most extensively so far.They were found to exhibit such intriguing properties as a hidden ferrimagnetic order with a hybridized band-Mott gap [3,4], tunable magnetoelectricity [5], the giant thermal Hall effect [6], the optical diode effect [7], and axion-type coupling at optical ME resonance [8].This compound family exhibits many other exotic properties, including diagonal ME susceptibility in Mn 2 Mo 3 O 8 [9][10][11] and successive electric polarization transitions and nonreciprocal light propagation in Co 2 Mo 3 O 8 [12][13][14].
The A 2 Mo 3 O 8 structure is shown in Figs.1(a) and 1 (b).This structure is noncentrosymmetric and polar.It consists of honeycomb layers formed by alternating NiO 4 tetrahedra and NiO 6 octahedra.The formal spin of the Ni 2+ ions in the d 8 state is 1.These layers are separated by trimerized nonmagnetic Mo ions.The A=Fe, Mn, and Co compounds exhibit magnetic ordering temperatures in the range of 40 to 60 K [3,12].Their magnetic order within the honeycomb layers is of the simple collinear Néel type, in which the nearby spins are antiparallel and point along the c axis [15,16].
On the other hand, the magnetic properties of Ni 2 Mo 3 O 8 are unique among the compounds of the A 2 Mo 3 O 8 family.
Its magnetic ordering temperature (∼5.5 K) is anomalously low, and the magnetic anisotropy is of the easy-plane type, with magnetic moments lying predominantly in the ab plane [17].Recent powder diffraction measurements suggest that the magnetic order is noncollinear [18].There are several patterns of magnetic order typical of the honeycomb lattice.In addition to the Néel order described above, they are the stripy and zigzag orders.Both orders consist of ferromagnetic spin chains antiferromagnetically coupled to each other.The chains are of the linear and zigzag types in the stripy and zigzag order, respectively [see Figs.1(c) and 1(d)].The magnetic structure of Ni 2 Mo 3 O 8 proposed in Ref. [18] is very complex.The in-plane magnetic connection is of the stripy type for the ab-plane components of the spins, while the c-axis spin components exhibit the zigzag-type connection.Both the reduced magnetic transition temperature and the complex noncollinear structure indicate significant magnetic frustration.
Frustrated honeycomb magnetic lattices are of significant interest.For example, topological magnons have been predicted in the stripy and zigzag states [19].In the A 2 Mo 3 O 8 compounds, the potentially exotic magnetism is coupled to the electric polarization.In the specific case of Ni 2 Mo 3 O 8 , a rich and complex magnetoelectric behavior was observed [17].The polarization dependence on the applied magnetic field exhibited linear, parabolic, and plateau-like features, depending on the field direction and its magnitude.These unusual phenomena were discussed using symmetry-based local ME tensor analy- sis.It was proposed that both the spin current and p-d hybridization mechanisms could contribute to the electric polarization and the ME response.The magnetic structure is the necessary starting point for understanding the mechanism of the ME coupling and the possible observation of the topological magnons and exotic magnetic states in Ni 2 Mo 3 O 8 , as well as for the analysis of the potential nonreciprocity and exotic Hall effects, as found in the other compounds in the A 2 Mo 3 O 8 family.
In this work, we revise the magnetic structure of Ni 2 Mo 3 O 8 using both powder and single-crystal neutron diffraction.We find that the in-plane magnetic connection is of the stripy type for both the ab-plane and the c-axis spin components.This connection is referred to as stripy/stripy hereafter, in contrast to the stripy/zigzag connection found in the previous powder diffraction studies [18].The revised magnetic structure presents a significantly simpler magnetic connection than the one proposed previously.The ferromagnetic interlayer order of the c-axis spin components in our model is also distinct.neutron diffraction analysis were consistent, verifying the quality of our samples (see Appendix A for the details of sample growth).
The temperature dependences of the specific heat, dielectric constant, and pyroelectric currents of single crystals are measured using a physical property measurement system (PPMS; Quantum Design).Magnetic susceptibility was measured using either a superconducting quantum interference device (magnetic property measurement system, Quantum Design) or a vibrating sample magnetometer in a PPMS.For the orientation-dependent susceptibility measurements, the crystallographic axes of single crystals were predetermined via Laue x-ray diffraction (see Appendix A for the details of measurement methods).
Time-of-flight neutron diffraction measurements of polycrystalline and single-crystal samples were conducted on the WISH (Wide angle In a Single Histogram) [20] and SXD (Single Crystal Diffractometer) [21] instruments, respectively, at ISIS, United Kingdom.Nuclear and magnetic refinements were performed using the JANA2006 package [22].An absorption correction for the single-crystal data was employed analytically using a multifaceted crystal model [23] implemented in the SXD2001 software [21,24].The observed reflections with intensities I > 3.0 × σ(I) were used in the refinements for both the powder and single-crystal data; I and σ represent the intensity and standard deviation, respectively.Symmetry analysis is performed to determine the candidate magnetic space groups based on the observed magnetic wave vector q = (1/2, 0, 0) in the parent hexagonal notation using the Bilbao Crystallographic Server [25] (see the detailed information in Appendix A).H || ab and H || c measured in an applied magnetic field of 0.1 T. The large difference in the effective moments indicates significant magnetic anisotropy between the ab plane and c axis, consistent with previous reports [17,18].χ(T ) exhibits a sharp maximum at T N = 5. Figure 2(c) shows the temperature dependence of the magnetic specific heat divided by the temperature, C m /T .The magnetic specific heat is obtained by subtracting the lattice contribution from the nonmagnetic counterpart, Zn 2 Mo 3 O 8 .C m /T exhibits a gradual increase followed by a maximum at around T N as the temperature is lowered.The magnetic entropy S m was calculated by integrating C m /T with respect to temperature.Figure 2(d) shows that S m reaches the theoretically expected value Rln(2S + 1) = Rln3 at 30 K. S m decreases quite gradually and releases 53 % of its value through T N upon cooling.This indicates significant frustration in the magnetic system.The dielectric constant ϵ, electric polarization change ∆P (P (T ) -P (T = 30 K)), and pyroelectric current density J, shown in Figs.2(e) and 2(f), present sharp maxima (or a kink) at T N .The observed anomalies in the magnetic, thermodynamic, and electric results corroborate the increase in the polarization associated with the onset of the magnetic order, confirming strong ME coupling in Ni 2 Mo 3 O 8 .FIG. 6. Magnetic structural refinement fits for model 1 in Table I.The magnetic space group is PC na21.The inset shows the resolved (002), (111), and (200) peaks.The solid Gaussian curves are guides to the eye.The asterisk ( * ) marks the background peak; others are not marked for simplicity (see Fig. 5 for the rest).

IV. SINGLE-CRYSTAL X-RAY DIFFRACTION
Single-crystal X-ray diffraction was performed to determine the crystal structure at room temperature.The crystal exhibited high quality; more than 98 % of the detected Bragg reflections originated from a single hexagonal domain (1124 out of 1146 peaks).In our measurements, we have observed numerous weak reflections that could not be indexed within the published space group of Ni 2 Mo 3 O 8 (as well as other A 2 Mo 3 O 8 compounds), P 6 3 mc (No. 186).One possibility is that the symmetry of Ni 2 Mo 3 O 8 at room temperature is lower.Trigonal P 3m (No. 156) or P 3 (No.143) space groups, for example, would accommodate all the observed reflections.However, before making such conclusions, one must confirm that the forbidden reflections do not originate from multiple scattering.
For this purpose, we performed so-called azimuthal scans for a set of the forbidden and allowed peaks.The sample is rotated around the scattering vector Q in the azimuthal scan.The direction and the magnitude of Q are fixed, and the Bragg condition is therefore maintained.In such scans, the intensities of the forbidden peaks should show extreme variations because the conditions for the multiple scattering are broken and restored as the azimuthal angle changes.On the other hand, the intensity of the allowed peaks shows little variation (neglecting the absorption effects).
Figure 3 shows the evolution of the normalized intensities I i in a wide range of azimuthal angles Ψ i for seven representative reflections.To minimize the extrinsic effects (such as geometry-dependent x-ray absorption), the reflections are chosen to be nearby in momentum space.This set includes both types of forbidden reflections in the P 6 3 mc space group, (H H odd) and (0 0 odd), as well as allowed reflections for a cross-check, based on the reflection conditions for the general Wyckoff site of P 6 3 mc.Each raw intensity I i is first divided by the corresponding TABLE I. Magnetic refinement results of the powder neutron diffraction analysis for the magnetic space group PC na21 at 1.5 K. Our best results are listed as models 1 to 4. We also show our fit results using the published magnetic structures [18] (solutions 1 and 2).Note that solutions 1# and 2# are manually modified structures from solutions 1 and 2 to fit our data.GOF means the goodness of fit.Definitions of all reliability factors presented in this paper are given in Table XIV 3 clearly demonstrates that the intensities of the forbidden reflections vary widely, whereas the allowed reflections show little variation with the azimuthal angle.We therefore conclude that the forbidden reflections originate from multiple scattering [26] and that the space group of Ni 2 Mo 3 O 8 at room temperature is indeed P 6 3 mc, a polar group.We refined the crystallographic structure of Ni 2 Mo 3 O 8 using this space group and found results very similar to those obtained by neutron diffraction at 5.5 K, described in detail below.

V. POWDER NEUTRON DIFFRACTION
Neutron diffraction measurements were performed to determine the magnetic order of Ni 2 Mo 3 O 8 .We start with powder diffraction.The bulk magnetic characteristics of our high-quality polycrystalline samples can be found in the Appendix B (see Fig. 12).The nuclear structure was determined first.The data were collected at 5.5 K, just above the magnetic ordering temperature.The corresponding Rietveld refinement is shown in Fig. 4. It reveals the same noncentrosymmetric structure (P 6 3 mc) found at room temperature in our x-ray measurements, as well as in the literature [1,18], and shown in Figs.1(a) and 1(b).The refined lattice parameters and atomic positions are shown in Table III in Appendix C.This result is consistent with a smooth change in the bulk characteristics above the magnetic transition [Fig.2(f)].Tiny amounts of impurity phases of NiO, MoO 2 , and NiMoO 4 were detected and fitted with Le Bail fits (see Fig. 13 in Appendix C).The vertical bars corresponding to these phases are not shown in Fig. 4 because the mass fractions of these phases were very small (less than 2 %).
For the magnetic structural refinement, the data were collected at 1.5 K. Figure 5 shows both the 1.5 K and 5.5 K data for comparison.Several new peaks appear at 1.5 K, indicating the onset of the long-range magnetic order.These peaks are listed in Table VII.The magnetic order is indexed with propagation vector q = (1/2, 0, 0).That is, the high-temperature unit cell is doubled along the a axis, becoming orthorhombic.In this work, we use a nonstandard magnetic unit cell for the analysis of the magnetic structure for simplicity [shown as a black parallelogram in Fig. 7(a) below].Group analysis performed using the parent structure with the q = (1/2, 0, 0) magnetic propagation vector yields four possible maximal magnetic space groups: and P C mc2 1 (see Fig. 14 and Table IV in Appendix C).
Among the four candidate groups, P C mn2 1 and P C mc2 1 can be excluded because the (1 1 0) magnetic peaks are clearly observed at 1.5 K (see Fig. 5).Indeed, these groups allow only the magnetic moments along the b direction of the parent unit cell, as shown in Table IV in Appendix C. The (1 1 0) wave vector is parallel to the b axis, and therefore, no magnetic neutron diffraction signal is allowed [27] at this Q.Magnetic refinements using the P C ca2 1 group yielded poor results, as shown in detail in Appendix C.
Magnetic refinements using the P C na2 1 (No.33.154) magnetic space group, on the other hand, produced very good results.We therefore conclude that P C na2 1 is the correct magnetic group for Ni 2 Mo 3 O 8 .We tested various magnetic starting structures, distinguished by the relative sizes of the Ni 2+ moments on the tetrahedral and octahedral sites, the dominant in-or out-of-plane spin component, and the magnetic connection type (stripy/stripy or stripy/zigzag).The refinements produced four magnetic structures of similar fit quality.They are listed as models 1 through 4 in Table I.Table I contains the magnetic connection type, the refined magnetic moment vectors, and the standard reliability parameters.The obtained quality of the fits is illustrated in Fig. 6, which shows the refinement result for model 1.The inset demonstrates the high resolution of our measurements: the magnetic (1 1 1) peak is clearly resolved from the nearby (2 0 0) and (0 0 2) nuclear peaks.
The refined magnetic structures of models 1 and 2 are shown in Fig. 7.Both structures are noncollinear and of the fully stripy type; that is, they exhibit the stripy/stripy connection in the honeycomb planes.The ab spin components are quite similar in the two models.The main difference between models 1 and 2 is in the relative values of the c-axis spin components on the tetrahedral and octahedral sites (see the bottom panels in Fig. 7).Models 3 and 4 exhibit the same features, except for the ratio of the tetrahedral to the octahedral Ni 2+ magnetic moment values.This ratio is larger than 1 for models 1 and 2; it is smaller than 1 for models 3 and 4. In other words, the magnitudes of the octahedral and tetrahedral moments are interchanged as one goes from models 1 and 2 to models 3 and 4.
In addition to T = 1.5 K, neutron diffraction data were collected at other temperatures and analyzed.The magnetic peaks smoothly disappear at the magnetic transition temperature, as expected.The refined values of the Ni 2+ moments at the two different sites exhibit smooth monotonic curves and go to zero at the transition temperature.These data can be found in Figs.19

and 20 in Appendix C.
Importantly, all the refinements using the powder data converge to the stripy/stripy structure.This is qualitatively different from previously published results, which presented a complex stripy/zigzag connection in the honeycomb planes [18].We therefore undertook a thorough comparative analysis of our models and the published structures.We could not stabilize the stripy/zigzag solutions [18] in the refinement because they always converged to the stripy/stripy structures.We therefore fixed the magnetic connection and the corresponding moment components of solutions 1 and 2 [18] in the refinements.The corresponding refinement results are much worse than those of our stripy/stripy models in Table I.To find better solutions within the stripy/zigzag structure, we manually tweaked the magnetic moments and found two better candidates, labelled as Solutions 1# and 2# in Table I.As one can see from all the reliability factors, these fits are again consistently worse than those of models 1 to 4. This can be seen, in particular, from RFw(Mag), which is known as a good metric to compare how well different models match the same set of experimental data [28].
To illustrate the better fit quality of the stripy/stripy models directly, we compare the obtained fits in the selected ranges of the d spacing (2π/Q) in Fig. 8.The peaks shown in Fig. 8 are the strongest magnetic peaks representing the most reliable measurements and having a significant effect on the refinements (see the full range data in Fig. 17).We emphasize that all the stripy/stripy models in Table I exhibit better fits than those for all the stripy/zigzag solutions.The observed differences between the RF(Mag) and RFw(Mag) of the stripy/stripy and stripy/zigzag structures lie in the 0.11 − 0.8 range.Such differences are accepted as statistically significant in the community [29][30][31][32][33]. Based on the combined results of the powder diffraction data refinement, we conclude that the magnetic structure of Ni 2 Mo 3 O 8 is of the stripy/stripy type.It is, however, difficult to determine which structures are better among the four models from the powder diffraction data alone.

VI. SINGLE-CRYSTAL NEUTRON DIFFRACTION
We performed single-crystal neutron diffraction measurements to reduce the ambiguity left by the powder neutron diffraction analysis.This technique exhibits several advantages, such as good sensitivity to thermal parameters and nonoverlapped Bragg peaks, and is therefore often sensitive to the features inaccessible to powder diffraction.Following the same general approach, we first collected the data at 10 K and did structural refinements; 2341 nuclear Bragg peaks were fitted [see Fig. 9(a)].The obtained structural parameters are consistent with the structure at 5.5 K determined using pow-  II).Blue and red symbols represent nuclear and magnetic peaks, respectively.der neutron diffraction; compare Tables III and IX in Appendixes C and D.
At 1.5 K, 229 magnetic Bragg peaks at the wave vector q = (1/2, 0, 0) in the high-temperature parent cell were found.In the following refinements, we used the fixed nuclear structural parameters determined at 10 K, as well as the same set of nuclear peaks for a systematic analysis.Structural refinement was done first using the extended (magnetic) unit cell 2a×a×c.This was done to obtain the accurate scale factor and extinction parameters for the following magnetic refinements.Figure 9(b) shows that these parameters describe the nuclear peaks of the 1.5 K data very well, justifying this approach.We then employed models 1 to 4, described in the previous section, to refine the magnetic structure.Three orientational magnetic domains are possible in a single crystal.The populations of these domains were refined and found to be nearly identical; see Appendix D for the details.
The obtained parameters are listed in Table II for all the models.Reliability parameters, such as GOF, R(All), and wR(All), are very similar in Table II because there is a much larger number of the nuclear Bragg peaks (2341) than magnetic Bragg peaks (229) in magnetic refinements.Thus, R(Mag) and wR(Mag) need to be compared to determine a better magnetic structure in Table II.We note that R(Nuc) and wR(Nuc) values are all the same in Table II, which ensures reliable and TABLE II.Magnetic refinement results from single-crystal neutron diffraction analysis at 1.5 K.The magnetic space group is PC na21.The reliability parameters, R(Nuc) = 4.32 and wR(Nuc) = 5.10, are the same for all magnetic structures presented.

Condition
Model Atom fully-controlled structural calculations during the magnetic refinements in our analysis.While all these fits are comparable, the magnetic reliability factors, R(Mag) and wR(Mag), of model 1 are better than those of models 2−4.The latter models are therefore less favoured by our single-crystal neutron diffraction analysis.We note that the differences in the magnetic reliability factors of models 1 and 2 are marginal.The refinement results for models 1 and 2 are shown in Figs.9(c) and 9(d) for comparison.
We have also tested the previously published magnetic structure [18] using our experimental data.As in the powder diffraction analysis, the published stripy/zigzag initial structures converged to the stripy/stripy order.Specifically, solution 1 from Ref. [18] converged to model 2, and solution 2 converged to model 1.To estimate the goodness of fit differences, we applied the same procedure as in our powder diffraction data analysis.The magnetic connection and the corresponding moment components in the refinements were fixed.The results are shown in Table II.As in the powder diffraction measurements, stripy/stripy models fit the experimental data best.
The refined magnetic structures of models 1 and 2 in Table II show very little difference from the corresponding models obtained in the powder diffraction experiments and shown in Fig. 7. Considering the common features of models 1 and 2, we conclude that the magnetic structure of Ni 2 Mo 3 O 8 exhibits the stripy/stripy magnetic connection in the honeycomb planes, and the magnetic moment of tetrahedral Ni 2+ is larger than that of the octahedral Ni 2+ site.

VII. ORIENTATION-DEPENDENT MAGNETIC SUSCEPTIBILITY
To characterize the magnetic anisotropy of Ni 2 Mo 3 O 8 in the context of its magnetic structure and the observed metamagnetic transitions, we performed orientationdependent magnetic susceptibility measurements on a single crystal.Figure 10(a) shows the magnetic susceptibility for various directions of the applied magnetic field.These directions are described in Fig. 10(b).Figure 10(c) shows the low-temperature part of the data.The largest effect is the clear easy-plane magnetic anisotropy, which is consistent with the refined magnetic structure.The ab-plane anisotropy is much more subtle but still observable.Specifically, the Curie-Weiss fits using the highertemperature data reveal a systematic modulation of the Curie-Weiss temperature Θ CW and effective magnetic moment µ eff .Figure 10(d) shows that Θ CW reaches its maximum value for H || b, while the minimum is ob- served for H ⊥ b. µ eff exhibits the opposite behavior.These systematic tendencies reflect the local magnetic anisotropies and, possibly, those of the exchange interactions.They appear to be consistent with the direction of the magnetic moments at low temperatures; for instance, the in-plane components are perpendicular to the b axis.The full analysis should consider the magnetic domain populations and is therefore the subject of future work.
Another systematic trend for the in-plane magnetic anisotropy is found in magnetic field-dependent magnetization.Figure 10(e) presents the isothermal magnetization M (H) at 2 K for various directions of the external magnetic field.The corresponding differential magnetization dM/dH is shown in Fig. 10(f).Based on the direction of H, the in-plane data can be classified into two categories.One set of M (H) is for H parallel to the hexagonal axes, directions 2 (− a) and 4 (b).The sec-ond set includes the directions 1, 3, and 5, corresponding to [2 1 0], [1 1 0], and [1 2 0], respectively.This difference between the two sets is most clearly observed in the differential magnetization of Fig. 10(f).The two observed peaks correspond to the metamagnetic transitions first discussed in Refs.[17,18].A systematic broadening and shift of the metamagnetic transition fields, H low and H high [vertically dashed lines in Fig. 10(f)] are found.H high , for example, is larger by approximately 1 T for H // a and H // b (directions 2 and 4) than for the other measured directions.It also appears that the higher-field transition is broader for directions 2 and 4.
These observations can be qualitatively understood using the refined magnetic structure.The octahedral moments were proposed to be involved in the low-field transition, while the tetrahedral moment motion was presumed to be dominant in the high-field metamagnetic changes [17].It is well known that the spin-flop transition becomes sharper when the field direction is closer to the easy-axis direction [34].The broader H high for H // a and H // b (directions 2 and 4) is then understood in terms of the flop of the magnetic moments perpendicular to the b axis, as determined in our measurements.This result is valid even if magnetic domains are present because the principal hexagonal axes are perpendicular to the moment direction.The opposite trend observed for H low indicates a different mechanism, such as the establishment of a collinear structure along the c axis by the octahedral moments, as proposed in Ref. [17].While the described features are understood within our magnetic models qualitatively, detailed modeling will need measurements of the magnetic exchange parameters.This will require inelastic neutron experiments.

VIII. IMPLICATIONS ON THE MAGNETOELECTRIC COUPLING AND POLARIZATION
All the refinements of our neutron diffraction data, powder and single crystal, produce the stripy/stripy magnetic connection in the honeycomb planes.That is, all the components of the Ni 2+ spins exhibit the stripy pattern.Previous experiments [18] claimed the stripy/zigzag connection, in which the c-axis spin components order in the zigzag structure in the honeycomb planes.In addition, the interlayer connection of the c-axis spin components in our model is ferromagnetic, while it is antiferromagnetic in the published solutions.The two magnetic structures are compared in Fig. 11, using model 2 from our experiments, and the corresponding solution 1 from Ref. [18].As discussed above, all the starting structures of the stripy/zigzag type converged to the stripy/stripy ones for both single-crystal and powder data.We are therefore confident that the stripy/stripy structure is realized in Ni 2 Mo 3 O 8 .
We now discuss the implication of our results for the magnetoelectric coupling.Ni 2 Mo 3 O 8 exhibits significant changes in the electric polarization with the onset of the magnetic order [Fig.2(f)].The polarization also changes in the applied magnetic field, especially at the metamagnetic transitions [17].This indicates a strong coupling between the magnetic and electric order parameters.
Several mechanisms of such coupling are known.They include the spin-current mechanism that utilizes the Dzyaloshinskii-Moriya (DM) interaction and works in noncollinear systems [35], the exchange striction that works best in the collinear systems [36], and the p-d hybridization between spin and ligand ions mediated by the spin-orbit coupling [37] which works for various spin orders.The type of the magnetic order clearly plays a key role in these mechanisms.The ME properties of Ni 2 Mo 3 O 8 have been recently examined in the framework of these models [17] using a general symmetry-based ME tensor approach [38,39].The stripy/zigzag magnetic order was assumed in those studies.This analysis needs to be revised based on the stripy/stripy magnetic order found in our experiments.
The magnitudes of the terms in the local ME tensor are set by the specific spin arrangement, spin magnitudes, and type of the ME coupling.The global polarization is the summation of the local polarizations generated in the bonds between the Ni ions and in the single-ion [17] p = <i,j> where P αβγ ij are the local ME tensor and (α, β, γ) correspond to the Cartesian orthorhombic coordinate components, which are (x, y, z) in the magnetic phase; they are defined as x = 2a + b, y = b, z = c, respectively.S α i (S β j ) is the α (β) spin component for the i (j) site.<i, j> runs over the Ni-Ni bonds generating the polarization.For i = j, single-spin terms are generated in Eq. ( 1), which therefore combines the two-spin and single-spin ME tensors.
The allowed terms in the P αβγ ij tensor of Ni 2 Mo 3 O 8 [17] are restricted by the point group symmetry (mm2) [39] of the experimentally determined magnetic space group (P C na2 1 ) and is given by , 0, 0) (0, 0, 0) (0, 0, P yy(z) ij where the coordinates inside the parenthesis give the direction of the local polarization and the other two coordinates refer to the spin components.As there are no magnetic moment components along the y axis in the magnetic order of Ni 2 Mo 3 O 8 , Eq. ( 2) is simplified to which consists of two diagonal and two off-diagonal components.Accordingly, p ij in Eq. ( 1) is expressed as Here p ij gives the polarization of the corresponding Ni-Ni bond for the two spins (i ̸ = j), or the single spin (i = j).Since the magnetic space group is the same as the one analyzed in Ref. [17], the conclusions given there also apply to our structure.Both the spin current (from the noncollinear spins along the x axis in the xz plane) and the p-d hybridization (the single-ion terms) mechanisms may contribute to the observed ME effect in Ni 2 Mo 3 O 8 .While the same two mechanisms are allowed for the previously published and revised magnetic structures, the signs of the components in Eq. ( 4) differ for the two structures.In Eq. ( 4), all the signs except the one for the first term are the opposite.This can be seen in Fig. 11, which shows the connections between the ab-plane and c-axis spin components for these structures.The sign of the fourth term involving P zz(z) ij is important for understanding the experimentally measured polarization along the c axis.
For the full analysis of the microscopic interactions underlying the ME effect, one needs to know all the terms in the local ME tensor.So far, experiments have been done only for the c-axis polarization direction.Importantly, the signs of the ME tensor terms should be determined correctly.This involves the determination of the absolute direction of the polar c axis.This has not been done so far, and therefore, the published changes of the c-axis polarization, in fact, give only the absolute values.A detailed study of the ME tensor in Ni 2 Mo 3 O 8 , including the measurement of P x in Eq. ( 4), is highly desired for a better understanding of the microscopic interactions underlying the ME effect.Combined with the magnetic structure determined here, it should lead to a more comprehensive understanding of the microscopic origin of the peculiar ME properties of Ni 2 Mo 3 O 8 .

IX. DISCUSSION
Our measurements show that the magnetic structure of Ni 2 Mo 3 O 8 is noncollinear, of the stripy/stripy type in the honeycomb planes, the ferromagnetic along the c axis, and that the magnetic moment M T of the tetrahedral Ni 2+ site is larger than the moment M O on the octahedral site.The numerical values of the ratio, M T /M O , vary in a certain range, reflecting the systematic errors of the methods used, such as electron spin resonance (ESR) measurements [18] and powder neutron diffraction [18].Specifically, M T /M O values determined in our neutron experiments for model 1 are 1.83 and 2.14 for the powder and single-crystal experiments, respectively.The corresponding values for model 2 are 1.37 and 1.59.On the other hand, the ESR measurements give M T /M O =1.78, and powder neutron diffraction measurements from the literature give 1.21 and 2.24 for solutions 1 and 2, respectively.A comparison of the ESR and the neutron results does not therefore favor either of the two models.We note that the ratios obtained from our neutron diffraction analysis are at least consistent with those from powder neutron diffraction refinements of the reference (see Table XIII for a detailed comparison).
The observed dominant easy-plane magnetism can be explained in the framework of crystal field theory.Reference [40] presented a theoretical model that reproduces the dominant in-plane directions of the magnetic moments.The calculated single-ion anisotropies are different for the tetrahedral and octahedral Ni sites.The observed difference in the corresponding magnetic moment values also obviously results from the different Ni environments.A more detailed analysis requires experimental measurements of Ni crystal field levels by inelastic neutron or Raman spectroscopy.
The noncollinear magnetic order found in Ni 2 Mo 3 O 8 is rather complex, and the spin canting angles out of the ab plane are large.A recent theoretical work [40] introduced the spin Hamiltonian containing, in addition to the Heisenberg terms, bilinear interactions, single-ion anisotropy, and the DM interaction.It was proposed that the DM interaction is responsible for the spin canting.However, given the large observed values of the spin canting angles, anisotropic magnetic interactions and local anisotropies may also need to be considered to explain the observed noncollinear order.Inelastic neutron scattering measurements on single crystals will eventually be required to determine the spin Hamiltonian experimentally.
Finally, Ni 2 Mo 3 O 8 is a member of the compound family in which exotic phenomena as topological magnons, nonreciprocity, and unusual Hall effects are either observed or expected.Theoretical analysis of the revised magnetic structure of Ni 2 Mo 3 O 8 is needed to establish whether they should be sought in this compound.The A 2 Mo 3 O 8 compound family provides many intriguing opportunities for further study.The presence of two distinct magnetic sites is especially useful because it makes selective doping possible, creating various sublattices in the honeycomb structural motif.Nonmagnetic Zn doping is already known to create axion-type coupling, magnetic transitions, and diagonal magnetoelectricity in Fe

X. CONCLUSIONS
The magnetic structure of the pyroelectric honeycomb antiferromagnet Ni 2 Mo 3 O 8 was determined using combined powder and single-crystal neutron diffraction.The structure is noncollinear.The magnetic moment of Ni 2+ on the tetrahedral site is significantly larger than the moment on the octahedral site.The magnetic order in the honeycomb planes is of the stripy type for all the spin components, which is different from the previously proposed stripy/zigzag connection.The ferromagnetic interlayer order of the c-axis spin components in our model is also distinct.In our study, the referential magnetic orders converge into our models for both the powder and the single-crystal data.In single crystals, we found a subtle, but clearly detectable, magnetic anisotropy in the honeycomb planes, using orientation-dependent magnetic susceptibility measurements.It manifests through the systematic modulation of the Curie-Weiss temperature and the effective magnetic moment and through variation of the magnetic fields at which spin-flop transitions occur, consistent with the neutron diffraction analysis.Using our magnetic model, we revised the analysis of the magnetoelectric tensor in Ni 2 Mo 3 O 8 .Our results provide key input for future studies of the Ni-based compounds of the A 2 Mo 3 O 8 family, in which topological excitations, nonreciprocity, unusual magnetoelectricity, and other exotic effects are under investigation.
Polycrystalline Ni 2 Mo 3 O 8 is synthesized using a solidstate reaction.High-purity NiO:ZnO:Mo:MoO 3 powders (1.9:0.1:1:2) were mixed and pelletized.Five percent nonmagnetic ZnO is required to initiate the reaction; pure Ni 2 Mo 3 O 8 polycrystalline samples cannot be grown without it in the conventional solid-state reaction.The pellet is prepared and sealed in a vacuum quartz tube and sintered at 900 • C, 1000 • C, and 1100 • C for 5, 5, and 10 h, respectively, with several intermediate grindings to improve the powder quality.
Single crystals of Ni 2 Mo 3 O 8 are grown using a chemical vapor transport method.Three grams of NiO:Mo:MoO 3 = 2:1:2 mixture powder are placed in a vacuum quartz tube (diameter = 15 mm, length = 200 mm).Then, 0.1 g TeCl 4 is added as a transport agent.The quartz tube is placed in a tube furnace with a hot zone at 1000 • C and a cold zone at 850 • C for three weeks.Single crystals with masses of approximately 20 mg are observed in the cold zone.
Magnetic susceptibility measurements are performed using superconducting quantum interference device magnetometry between 300 or 350 K and the base temperature (2 or 3 K) under 0.1 T applied magnetic field with a Quantum Design MPMS-XL7 and PPMS-14 [Figs.2(a) and 2(b)]; zero-field-cooled and field-cooled measurements are performed when necessary.The temperature dependence of the dielectric constant, electric polarization, specific heat, and magnetic susceptibility measurement of single crystals are measured using a PPMS-9 [Figs.2(c)-2(f)], along with the magnetic susceptibility measurement of polycrystalline powder (Fig. 12).Orientation-dependent magnetic susceptibility measurements (Fig. 10) are performed in the temperature range of 2 -300 K under an applied magnetic field of 0.1 T using a PPMS-14 equipped with a vibrating sample magnetometer.Magnetization data are collected at 2 K up to 14 T by rotating the field direction by 30 • within the ab plane, as well as along the c axis.
A time-of-flight neutron diffraction measurement on polycrystalline samples is conducted on WISH at ISIS, United Kingdom, to determine nuclear and magnetic structures and examine a magnetic transition.Multibank detectors are utilized at different scattering angles; 2.54 g of the powder sample are packed in a cylindrical vanadium can (diameter = 6 mm, height = 26 mm).The entire vanadium cell is placed in a standard cryostat.The data are collected from 1.5 to 10 K with a step of 0.5 K going through an antiferromagnetic transition temperature at around 5 K, followed by a measurement at 100 K. Longer measurements at 1.5 and 10 K are conducted for 1 h to obtain high-quality data for detailed refinements.
Single-crystal neutron diffraction experiments are performed on SXD at ISIS, United Kingdom, where the time-of-flight Laue technique is used to access large threedimensional volumes of reciprocal space in a single measurement.The data are collected by rotating mounted single crystals.We measured the data at two rotation angles at 10 K (the paramagnetic state) and two identical rotation angles at 1.5 K (a magnetically ordered phase) and added three more rotation angles at 1.5 K to increase the statistics.One rotation of the data is collected for 6.5 h.
Nuclear and magnetic refinements are performed using the JANA2006 package [22].As extinction correction is key for reliable refinements using single-crystal neutron diffraction data, we test all isotropic extinction models implemented in the JANA2006 and use the best model for the refinements; that is, the Becker & Coppens model (a mixed Lorentzian model) [22,41] was used for structural refinements (unless other-wise specified in this paper) at the paramagnetic temperature, which was consistently used for magnetic refinements at lower temperatures.We confirmed a single structural domain for measured single crystals and three nearly equally populated magnetic domains in the magnetic phase.Figure 12 shows the temperature-dependent magnetic susceptibility of the polycrystalline sample at H = 0.1 T. We used this powder for the powder neutron diffraction measurement presented in Sec.V. A tentative antiferromagnetic transition from Fig. 12(a) at approximately 5.5 K is confirmed by powder neutron diffraction, as presented in Fig. 20 and Sec.V. From the fits using the Curie-Weiss law, we obtain Θ CW = −29.8(6)K and µ eff = 3.57 µ B , which is somewhat smaller than previously reported values [17,18].The powder susceptibility can be estimated from the single-crystal data (Fig. 2) as χ powder = (2χ ab + χ c )/3, which is located between χ ab and χ c (not shown).We note that a larger susceptibility of the powder sample measured could be attributed to the instrumental difference in the PPMS (powder; Fig. 12) and the MPMS-XL7 [single crystals; Figs.2(a) and 2(b)], different diamagnetic backgrounds, or the sample dependence between the polycrystalline powder and single crystal.

Appendix C: Powder neutron diffraction
This section provides details for powder neutron diffraction analysis.In the reduction of the data, a pair of bank data (in the left and right direction of the scattering) is summed to increase the signal-to-noise ratios of the peaks; they are labeled as bank 1 to bank 5 detectors in this paper (2θ = 27.08 • , 58.33 • , 90 • , 121.66 • , and 152.83 • from bank 1 to bank 5, respectively) with increasing order of the scattering angle.
In the data, we found tiny impurity peaks.Like for three impurity phases, shown in Fig. 4, we did not mark positions of the nuclear Bragg peaks of impurity phases (NiO, MoO 2 , NiMoO 4 ) in this paper for simplicity (except in Fig. 13) because their mass fractions are tiny.However, we included them in all our refinements with the powder neutron diffraction data shown in this paper to determine more precise nuclear and magnetic structures.For completeness, we explicitly compare the structural refinements with and without these three impurity phases in Fig. 13; this comparison clearly reveals a tiny difference in the fits.Table III presents extracted structural parameters from the neutron diffraction analysis at 5.5 K shown in Fig. 13(b).
Among maximal subgroups for the magnetic structural determination (Fig. 14), we describe how we rule out P C ca2 1 in detail.As discussed in the main text, we performed magnetic refinements in two magnetic space groups (P C ca2 1 and P C na2 1 ) because they can support magnetic moments along the a or c axis, unlike the other two.Indeed, in Fig. 6, the magnetic structure in P C na2 1 explains the 1.5 K data very well.However, the magnetic refinements using P C ca2 1 do not work (see Fig. 15).Also, for the more intuitive comparison, we set M c = 0 in the refined magnetic structures, by illustrating a qual-  itative difference in the fitted magnetic structures in P C na2 1 and P C ca2 1 (Fig. 16).The simplified magnetic structure in P C na2 1 provides a stripy magnetic structure both in the ab plane and along c axis; the interlayer coupling is antiferromagnetic in the fitted structure [Fig.16(a)].However, P C ca2 1 cannot support it by symmetry (Table IV); instead, it allows for either a stripy order in the ab plane with ferromagnetic interlayer coupling [Fig.16(b)] or a zigzag order in the ab plane with antiferromagnetic interlayer coupling [Fig.16(c)].This qualitative distinction clearly explains why the magnetic refine-FIG.14.Diagram of the possible subgroups of the paramagnetic parent space group P 63mc with a given q.The diagram is generated using the k -SUBGROUPSMAG software available on Bilbao Crystallographic Server [25].The symmetries of the four maximal subgroups, which are emphasized by red circles, are listed in Table IV.
ment did not work completely with magnetic structures in P C ca2 1 , as explicitly illustrated in Fig. 15.A representative list of refined magnetic structures in P C ca2 1 is presented in Table V to compare the reliability factors with those from P C na2 1 (Table I).We now provide full details about magnetic refinements using models 1 and 2 and solutions 1 and 2. In the main text, we showed that the magnetic structure of Ni 2 Mo 3 O 8 is stripy/stripy in the plane with the fer- romagnetic interlayer coupling, which differs from the previously reported stripy/zigzag structure in the plane with the antiferromagnetic interlayer coupling [18].In Ref. [18], two statistically identical stripy/zigzag magnetic structures were found.To test their solutions with our neutron diffraction data, we first tried to fit our data starting with the stripy/zigzag structures, but the fit always converged to the stripy/stripy model after running a few refinement cycles.This means that the stripy/stripy model is more stable than the stripy/zigzag solution.Thus, we simulated the reported solutions 1 and 2 with constraints in the fit, and slightly adjusted the moment values to obtain the best possible results, labeled as solutions 1# and 2#.The agreement factors for all four cases are listed in Table I for direct comparison (also see Fig. 8).A comparison of the fitting result and agreement factors indicates that the stripy/stripy models yield slightly better fitting results, compared to the stripy/zigzag solutions.Moreover, the same conclusion was reached when we repeated magnetic refinements using the pure magnetic signal at 1.5 K (see Fig. 17)).Furthermore, a consistent result was obtained from the single-crystal neutron diffraction analysis (Table II).
For completeness, we present magnetic refinement results of additional magnetic structural models tested with powder neutron diffraction data collected at 1.5 K, as summarized in Table VI.In Table I, we show that the magnetic space group P C na2 1 fits all observed magnetic peaks well in the stripy/stripy model.We note that the M c value of the Ni2 (Ni1) site in model 1 (2) is very small (i.e., 0.18 µ B and 0.10 µ B , respectively).This may imply an absence of M c in the accurate magnetic structure.Therefore, we tested this hypothesis by imposing M c = 0 in models 1 to 4, labeled models 5 to 8 in Table VI.Then, refinements became marginally worse, as explicitly compared in Fig. 18.Both the (1, 0, 0) and (1, 1, 0) peaks are slightly less fitted in models 5 and 6, and the agreement factors (GOF, wRp, and Rp) increase compared to those in models 1 and 2. This means that such small and finite M c components can be determined in models 1 and 2 within the resolution of our high-quality data.
In addition, as cross-checks, we repeated the magnetic refinements for models 1 and 2 by using pure magnetic signals at 1.5 K, as shown in Fig. 17; the insets confirm that fits work better in models 1 and 2 than solutions 1# and 2#, consistent with all other analysis.Also, we describe how our high-quality powder neutron diffraction data resolve the previous confusion over the indexing of the magnetic Bragg peak.As shown in Fig. 19(b) (see also Fig. 6), we can clearly separate nuclear and magnetic Bragg peaks located nearby in the narrow range of the d spacing.At 1.5 K, the intensity near the (200) nuclear peak (d ∼ 4.96 Å), which was originally indexed as (100) in the paramagnetic phase, increases because of the appearance of new magnetic peaks of (111) and (111) reflections, which our high-resolution data can resolve in magnetic refinements.In a previous report [18], this magnetic peak was indexed as (004), and it was taken as evidence of a significant perpendicular magnetic moment M c , and the irreducible representations Γ1 and Γ3 were disregarded because they cannot have the M c magnetic moment component.However, we cannot discard these irreducible representations based on their arguments because the (004) reflection is not a magnetic peak.Also, even if (004) was a magnetic reflection, we could not confirm the presence of M c because neutrons can see the magnetic moment only perpendicular to the wave vector.Therefore, the presence of "magnetic" (00L) Bragg peak in the neutron diffraction data would mean a possible magnetic moment within the ab plane.Note that we did not observe a magnetic (001) Bragg peak (see Table VII).
Further, we collected neutron diffraction data with increasing temperature to study the evolution of magnetic Bragg peaks towards the paramagnetic phase.Figure 19 shows the temperature-dependent data between 1.5 K and 10 K.It shows the data for the selective range of the d spacing for clarity, where we can see that the magnetic Bragg peaks are gradually suppressed with heating and become absent at around 5.5 K.For a more accurate de- termination of magnetic moments with temperature (using model 1 in Table I), the all-bank and bank 2 data are used separately in magnetic refinements.The transition temperature extracted from the neutron diffraction data via the fit is between 5 K (Fig. 20) and 5.5 K (not shown), and it is consistent with T N found from magnetic, thermodynamic, and electric measurements (Fig. 2).We also attempted to fit the data above 5.5 K, although it does not reveal magnetic Bragg peaks, to estimate the error bars of fitted moments.Then, we set the zero magnetic moment values for the data above 5.5 K in refinements because no magnetic peaks are observed above 5 K.In Table VIII, we present the results of the magnetic refine-ments, which are used for the critical analysis, shown in Fig. 20.
Figure 20 presents the change in the fitted magnetic moments upon a change in temperature, which clearly demonstrates the reliably fitted moments in both types of the data used (see Table VIII for refinement results).Figure 20(a) shows that the fitted magnetic moments of Ni 2+ at the tetrahedral (Ni1) site are almost twice those of the octahedral (Ni2) site.However, the normalized magnetic moments at both sites show a comparable temperature dependence, as illustrated in Fig. 20(b).The solid lines are from the critical analysis using M (T ) = M 0 (1 − T /T c) β [42].The temperature dependence of Enlarged views for selective magnetic Bragg peaks in (a) reveal a strong disagreement in the fit when using PC ca21, compared to I calc from PC na21 (Fig. 6).Model A (stripy/stripy) of Table V was used as an example.
the Ni1 and Ni2 moments exhibits a power-law behavior with critical exponents of β = 0.31(1) and 0.21(2), respectively.It might be understood by a conventional interpretation based on the criticality theory; the β value at the Ni1 site is close to the value of the three-dimensional (3D) Heisenberg (β = 0.36) or 3D XY (β = 0.35) interactions, whereas its value at the Ni2 site is close to that of the two-dimensional XY interaction (β = 0.23) [43,44].However, we should recall that this could be misleading because there could be a non-negligible coupling between two Ni sites.This Appendix provides details about single-crystal neutron diffraction analysis.First, the structural parameters extracted from refinements using the neutron diffraction data at 10 K are listed in Table IX.Also, we present refined magnetic structures in Figs.21(a) and 21(b) which are consistent with those from the powder result.Our analysis found three nearly equally populated magnetic domains, as shown in real space in Fig. 21(c), and the corresponding reciprocal lattices are shown in Fig. 21(d).We observed all three types of magnetic Bragg peaks in the single-crystal neutron diffraction data measured at 1.5 K from three magnetic domains.We tested a single magnetic domain in magnetic refinements, but the test did not work because of the finite intensity of the magnetic Bragg peaks belonging to other domains.
In the magnetic refinements, shown in Figs.9(b)-9(d) in the main text, we use the extended unit cell, 2a×a×c.For completeness, we present a repetitive structural refinement using the 10 K data with the extended unit cell (Fig. 22) to obtain accurate structural parameters such as extinction and scale parameters, which are fixed for more reliable magnetic refinements with the 1.5 K data in our analysis [Figs.9(b)-9(d)].
We present the structural refinements results at 1.5 K using only the allowed Bragg peaks.Our aims are twofold: (i) to confirm that the structural parameters extracted from the 10 K data can be reliably used for magnetic refinements with 1.5 K data and (ii) to see whether our data are sensitive to any possible structural distortion from a hexagonal to an orthorhombic structure via a magnetic transition.We used the only identical nuclear Bragg peaks (2341 nuclear reflections) collected from two identical rotation angles at 1.5 and 10 K while fixing all parameters except the atomic position and thermal parameters in the refinements to determine whether we can obtain any hint of a structural transition.We use nuclear space group P 6 3 mc and magnetic space group P C na2 1 in the extended unit cell (2a×a×c), as shown in Figs.22(c) and 22(d), respectively.The obtained nuclear structures are comparable.Thus, we did not observe any evidence of the structural distortion in the nuclear structure from 10 to 1.5 K within the resolution of our single-crystal neutron diffraction data.
We discuss some excluded forbidden peaks in Figs.9(c)-9(e) in the main text because they are not fitted even with the lowest magnetic space group.Fig- ures 23(a) and 23(b) present the magnetic refinement results after including the additional peaks, where one can see weaker structure factors F calc in both models 1 and 2. If they are magnetic peaks, their existence may imply the symmetry of a true magnetic space group is lower than P C na2 1 .We fitted the data to the lower-symmetric magnetic space groups to test this possibility.Lowering the symmetry of all four orthorhombic maximal subgroups resulted in six possible monoclinic candidates [P A c (5), P a 2 1 (9), P C c (6), P C c (7), P a 2 1 (10), and P C m (8)] and one triclinic candidate [P S 1 (11); Fig. 14].The number in parentheses denotes the subgroup number in the group analysis, as depicted in Fig. 14.Using the two-rotation data, we tested all possible subgroups originating from P C na2 1 , e.g., P A c (5), P C c (6), P a 2 1 (10), and the least symmetric magnetic space group, P S 1, based on model 1.However, none of them explained the additional peaks.Thus, these additional reflections at 1.5 K are not sensitive to the magnetic order.They could be attributed to extrinsic effects, such as low statistics and/or multiple diffraction [26], which cannot be considered in the simple  extinction model in our analysis.We note that the forbidden signals at 10 K were at least already subtracted from the 1.5 K data for identical peaks before any magnetic refinements were done.We obtain the same result that the additional peaks from all five-rotation data are not explained by the lowest-symmetric magnetic space group (Fig. 23).Therefore, we removed these Bragg peaks from  For completeness, we present results of magnetic structural refinements using all five-rotation data collected only at 1.5 K, which include the two identical rotation data at both temperatures; in detail, we collected additional three-rotation-angle data only at 1.5 K to maximize the signal-to-noise ratio of the magnetic Bragg peaks.The same analysis as described in the main text using two-rotation single-crystal data was repeated with the five-rotation single-crystal data (13126 nuclear and 537 magnetic reflections).We fitted the 13126 nuclear reflections collected at 1.5 K with an extended unit cell, 2a×a×c, to obtain accurate structural parameters, including extinction parameters.The obtained structural parameters were then fixed for subsequent magnetic refinements.
The fitted magnetic moments with agreement param-    Those from the electron spin resonance (ESR) [18], powder neutron diffraction (PND) refinements of the reference (solutions 1 and 2) [18], PND refinements of our data, single-crystal neutron diffraction (SND) with the two-rotation data (TR), and with the five-rotation data (FR) are compared.Solutions 1 and 2 are obtained directly from the reference, while the solutions marked with # are the slightly modified solutions 1 and 2 that fit our pND data.

FIG. 2 .
FIG. 2. (a) In-plane and c-axis magnetic susceptibility versus temperature.The open symbols and solid lines are the data from different crystals.The single crystal corresponding to the solid-line data is used for the orientation-dependent susceptibility measurements shown in Fig. 10.The solid-line data for H || ab is taken for the [2 1 0] field direction, corresponding to direction 1 in Fig. 10.(b) Inverse magnetic susceptibility.The dashed lines are fits to the Curie-Weiss law.(c) Magnetic specific heat divided by the temperature.(d) Magnetic entropy extracted from the data shown in (c).(e) Dielectric constant along the c axis as a function of temperature.(f) Temperature dependence of the polarization change ∆P [P (T ) -P (T = 30 K)] and the pyroelectric current density J.

Figure 2
Figure 2 presents the magnetic, thermodynamic, and electric properties of Ni 2 Mo 3 O 8 single crystals.Figure 2(a) shows the magnetic susceptibility χ(T ) for 5 K, signifying the antiferromagnetic long-range ordering.χ(T ) is well reproduced by the Curie-Weiss law with Θ CW = −124.22K and µ eff = 3.93 µ B for H || c, and Θ CW = −49.59K and µ eff = 4.02 µ B for H || ab.The obtained Θ CW values are consistent with the previous single-crystal results [17].The large difference between these values again indicates the substantial magnetic anisotropy of Ni 2 Mo 3 O 8 .

FIG. 7 .
FIG. 7. Refined magnetic structures of (a) model 1 and (b) model 2 at 1.5 K (Table I, −0.1 ≤ z ≤ 0.1).Black (gray) symbols indicate the octahedral (tetrahedral) Ni sites.The parallelogram-shaped black solid lines delineate the magnetic unit cell in the nonstandard setting, 2a×a×c, used in our work.Here a and c represent the crystallographic axes of the high-temperature hexagonal unit cell.Red and blue lines represent the connection of the ab-plane and c-axis spin components, respectively.

FIG. 8 .
FIG. 8. Comparison of selected magnetic refinement results at 1.5 K using (a) model 1, (b) model 2, (c) solution 1#, and (d) solution 2#.Note that models 1 and 2 are stripy/stripy and solutions 1# and 2# are stripy/zigzag.All plots are based on results in Table I for consistent comparison.

FIG. 9 .
FIG. 9. Magnetic refinement results for single-crystal neutron diffraction data.F calc and F obs are the calculated and observed structure factors, respectively.(a) Structural refinement at 10 K in the high-temperature space group.(b) Structural refinement at 1.5 K in the extended unit cell obtained as described in the text.(c) and (d) Magnetic refinements using models 1 and 2, respectively (see TableII).Blue and red symbols represent nuclear and magnetic peaks, respectively.

FIG. 10 .
FIG. 10.Orientation-dependent magnetic susceptibility.(a) Magnetic susceptibility versus temperature for the different directions of applied magnetic field (0.1 T).(b) The field directions with respect to the crystal.The corresponding crystallographic axes are described in the text.(c) The lowtemperature region of the data shown in (a).In (a) and (c), the data variation is smaller than the symbol size for the in-plane fields.(d) Angular dependence of the Curie-Weiss temperature ΘCW (pink circles) and the effective magnetic moment µ eff (black squares) extracted from (a).Dashed vertical lines denote the crystallographic directions a and b.(e) Isothermal magnetization at 2 K for various field directions.(f) The first derivative of the isothermal magnetization at 2 K. Gray and black dashed vertical lines indicate the metamagnetic transition fields, H low and H high , respectively.

FIG. 11 .
FIG. 11.Comparison of (a) the stripy/stripy and (b) the stripy/zigzag magnetic structures of Ni2Mo3O8.The stripy/stripy structure is shown using the refined model 2 for the single crystal obtained in our experiments.The stripy/zigzag structure is represented by the corresponding solution 1 directly from Ref. [18].Red and blue lines represent the connection of the ab-plane and c-axis spin components, respectively.

FIG. 12 .
FIG. 12. (a) Magnetic susceptibility of the polycrystalline sample Ni2Mo3O8 used for neutron diffraction (H = 0.1 T).(b) Inverse magnetic susceptibility.A dashed red line notes the Curie-Weiss fit.
Appendix B: Powder magnetic susceptibility

FIG. 15 .
FIG. 15.(a) Results from magnetic structural refinements in PC ca21 (model A of Table V) using the 1.5 K data.(b) -(e)Enlarged views for selective magnetic Bragg peaks in (a) reveal a strong disagreement in the fit when using PC ca21, compared to I calc from PC na21 (Fig.6).Model A (stripy/stripy) of Table V was used as an example.

FIG. 16 .
FIG. 16.Refined magnetic structures after setting Mc = 0 to compare the intralayer (within the ab plane) and interlayer (along the c axis) coupling in a simpler way.(a) Stripy order in the plane with an antiferromagnetic interlayer coupling (simplified from model 1 in Table I in PC na21).(b) Zigzag order in the plane with a ferromagnetic interlayer coupling (model A of Table V in PC ca21).(c) Zigzag order in the plane with an antiferromagnetic interlayer coupling (model B from Table V in PC ca21).The moment values of (b) and (c) are scaled at both Ni sites for a visualization purposes.Red lines represent the connection of the ab-plane spin component.

FIG. 17 .
FIG. 17.Comparison of results from magnetic structural refinements using the powder data in models 1 and 2 and solutions 1# and 2#.Only magnetic signal at 1.5 K (after the 10 K data are subtracted) is used.(a) Model 1.(b) Model 2. (c) Solution 1#.(d) Solution 2#.Insets enlarging the magnetic (100) peak clearly confirm the better fit with models 1 and 2 compared to solutions 1# and 2#.

FIG. 19 .
FIG. 19.Powder neutron diffraction data of Ni2Mo3O8 collected from 1.5 to 10 K. Positions of Bragg reflections are obtained from fitting the 1.5 K data.

Figs. 9
Figs. 9(b)-9(d) in the main text for simplicity.Magnetic refinement results are presented in Tables X and XI for the two-rotation-angle and five-rotation-angle data, respectively, for comparison.For completeness, we present results of magnetic structural refinements using all five-rotation data collected only at 1.5 K, which include the two identical rotation data at both temperatures; in detail, we collected additional three-rotation-angle data only at 1.5 K to maximize the signal-to-noise ratio of the magnetic Bragg peaks.The same analysis as described in the main text using two-rotation single-crystal data was repeated with the five-rotation single-crystal data (13126 nuclear and 537 magnetic reflections).We fitted the 13126 nuclear reflections collected at 1.5 K with an extended unit cell, 2a×a×c, to obtain accurate structural parameters, including extinction parameters.The obtained structural parameters were then fixed for subsequent magnetic refinements.The fitted magnetic moments with agreement param-

FIG. 20 .
FIG. 20.Evolution of fitted magnetic moments with temperature using the all-bank and bank 2 data separately.(a) Magnetic moments at the Ni1 and Ni2 sites.(b) Normalized magnetic moments.Solid red lines in (a) are from the fits (see texts).Model 1 of Table I is used in the refinement as a representative structure.

FIG. 21 .
FIG. 21.(a) and (b) Refined magnetic structures from Figs. 9(c) and 9(d), respectively.The parallelogram-shaped black solid line denotes the magnetic unit cell.The crystallographic a and b axes based on the extended unit cell are presented.Black (gray) ions are octahedral (tetrahedral) Ni sites.Red and blue lines represent the connection of the abplane and c-axis spin components, respectively.(c) Three equivalent magnetic domains in a collinear stripy order with Mc = 0 for simplicity.(d) Reciprocal space of (c).The blue, red, and green dashed lines show the reciprocal lattice, corresponding to the magnetic unit cells of domains 1, 2, and 3 in real space, respectively. .
2−x Zn x Mo 3 O 8 [6-8].However, studies of doping effects in other compounds of this family are scarce.Ni 2 Mo 3 O 8 is a unique member of this series, showing a rather distinct magnetism.This calls for examining the doping effects in this compound, for instance, in Ni 2−x Mg x Mo 3 O 8 and Ni 2−x Zn x Mo 3 O 8 . ).

TABLE IV .
Symmetry relations of the four maximal magnetic space groups (MSGs) related to Ni2Mo3O8.

TABLE V .
Representative magnetic structural models tested in PC ca21 with obtained reliability parameters.

TABLE VI .
Results from magnetic structural refinements in PC na21 using the all-bank data at 1.5 K. Reliability parameters for each bank data are provided separately for a more detailed comparison.

TABLE VII .
List of magnetic Bragg peaks observed in the powder neutron diffraction data at 1.5 K (bank 2) based on magnetic refinements using the pure magnetic signal.I obs represents the observed intensity.

TABLE XI .
Magnetic refinements results in PC na21 with the 1.5 K single-crystal neutron diffraction data (the five-rotation data, with 537 forbidden peaks).

TABLE XII .
Magnetic refinements results in lower-symmetry magnetic space groups with the 1.5 K single-crystal neutron diffraction data (the five-rotation data, with 537 forbidden peaks).

TABLE XIII .
Comparison of magnetic moments of two Ni sites of Ni2Mo3O8 and their ratios in our study and the literature.