Investigation into the role of the orbital moment in a series of isostructural weak ferromagnets

The orbital contribution to the magnetic moment of the transition metal ion in the isostructural weak ferromagnets ACO$_3$ (A=Mn,Co,Ni) and FeBO$_3$ was investigated by a combination of first-principles calculations, non-resonant x-ray magnetic scattering and x-ray magnetic circular dichroism. A non-trivial evolution of the orbital moment as a function of the $3d$ orbitals filling is revealed, with an anomalously large value found in the Co member of the family. Here, the coupling between magnetic and lattice degrees of freedom produced by the spin-orbit interaction results in a large single-ion anisotropy and a peculiar magnetic-moment-induced electron cloud distortion, evidenced by the appearance of a subtle scattering amplitude at space group-forbidden reflections and significant magnetostrictive effects. Our results, which complement a previous investigation on the sign of the Dzyaloshinskii$-$Moriya interaction across the series, highlight the importance of spin-orbit coupling in the physics of weak ferromagnets and prove the ability of modern first-principles calculations to predict the properties of materials where the Dzyaloshinskii$-$Moriya interaction is a fundamental ingredient of the magnetic Hamiltonian.

The orbital contribution to the magnetic moment of the transition metal ion in the isostructural weak ferromagnets ACO3 (A=Mn,Co,Ni) and FeBO3 was investigated by a combination of first-principles calculations, non-resonant x-ray magnetic scattering and x-ray magnetic circular dichroism. A non-trivial evolution of the orbital moment as a function of the 3d orbitals filling is revealed, with an anomalously large value found in the Co member of the family. Here, the coupling between magnetic and lattice degrees of freedom produced by the spin-orbit interaction results in a large single-ion anisotropy and a peculiar magnetic-moment-induced electron cloud distortion, evidenced by the appearance of a subtle scattering amplitude at space group-forbidden reflections and significant magnetostrictive effects. Our results, which complement a previous investigation on the sign of the Dzyaloshinskii−Moriya interaction across the series, highlight the importance of spinorbit coupling in the physics of weak ferromagnets and prove the ability of modern first-principles calculations to predict the properties of materials where the Dzyaloshinskii−Moriya interaction is a fundamental ingredient of the magnetic Hamiltonian.

I. INTRODUCTION
Recent studies of the weak ferromagnetic carbonates ACO 3 (A=Mn,Co,Ni) 1 and FeBO 3 2 represent the first systematic experimental and theoretical investigation of the changes in the sign and magnitude of the Dzyaloshinskii−Moriya interaction (DMI) across a series of insulating 3d transition metal (TM) compounds. The combination of novel resonant x-ray diffraction technique and modern first-principles calculations revealed a dramatic evolution of the sign of the DMI as the 3d orbitals of the TM are gradually filled with electrons. The ability to accurately model the DMI is essential for the fundamental understanding of a plethora of exotic noncollinear magnetic ground states, such as spin spirals 3 and Skyrmions 4-6 , and their exploitation as candidate materials for spintronics applications.
The DMI has its microscopic origin in spin-orbit coupling (SOC) 7,8 . In the common paradigm of the physics of TM oxides, SOC is regarded as negligible for 3d electrons, where its role is merely as a small perturbation to the ground-state Hamiltonian 9 . This contrasts with the case of heavier (4d and 5d) TM compounds, where SOC competes with the crystal field and other relevant en-ergy scales on an equal footing and gives rise to more exotic ground states 10 . Nonetheless, even for 3d TM compounds, SOC is expected to have a significant impact on the magnetic properties of the system whenever a finite orbital moment is present 9 . A substantial unquenched orbital contribution to the magnetic moment has been indeed reported for several 3d oxides [11][12][13][14][15][16][17][18] . In this case, the coupling between spin and orbital moment caused by the spin-orbit interaction can generally produce a strong magnetoelastic coupling and lead to the appearance of a large single-ion anisotropy and magnetostrictive effects 9 . The magnetic properties will then considerably differ from the case a spin-only system with quenched orbital degrees of freedom.
A careful determination of the strength of the orbital moment and its impact on the magnetic ground state is of particular interest for the weak ferromagnets ACO 3 (A=Mn,Co,Ni) and FeBO 3 , where SOC, and the resulting DMI, underpins one of the most peculiar aspects of the physics of this system, i.e. the existence of a weak net magnetization. Weak ferromagnets also represent an ideal model system for (i) their manageable magnetic unit cell, in contrast to the more complex spin-spiral and Skyrmion states of interest in light of spintronics arXiv:1805.12423v1 [cond-mat.str-el] 31 May 2018 applications and (ii) the availability of state of the art calculations 1,2 which can be conveniently used to predict the relative orbital and spin contribution to the TM magnetic moment.
In this paper we present a detailed investigation into the role of the orbital moment in the isostructural weak ferromagnets ACO 3 (A=Mn,Co,Ni) and FeBO 3 by means of a combination of theoretical calculations, Nonresonant X-ray Magnetic Scattering (NXMS) and X-ray Magnetic Circular Dichroism (XMCD). While MnCO 3 and FeBO 3 behave as almost pure spin systems, a sizeable orbital contribution to the magnetic moment was found in CoCO 3 and NiCO 3 . In particular, an anomalously large orbital moment is present in the Co compound which results in a remarkable coupling between lattice and magnetic degrees of freedom. The latter is unveiled by a sizeable magnetocrystalline anisotropy of the magnetic interactions and, more spectacularly, by the emergence of an unusual space-group-forbidden scattering process.
The paper is organized as follows. A brief description of the samples is given in § II, while the theoretical calculations and the NXMS and XMCD experimental setup are outlined in § III. § IV presents the results of the DFT calculations and the NXMS measurements on the orbital contribution to the magnetic moment across the series. § V includes specific findings on CoCO 3 , in particular: § V A outlines the XMCD measurements used to support the NXMS results on the size of the orbital moment, § V B discusses the role of the magnetocrystalline anisotropy in the NXMS data and § V C deals with the space-group forbidden scattering and its microscopical interpretation based on the multiplet calculations. Finally, the concluding remarks are presented in § VI.

II. SAMPLES
The weak ferromagnetic carbonates ACO 3 (A=Mn,Co,Ni) and FeBO 3 are isostructural compounds, with the trigonal R3c crystal symmetry [19][20][21][22] (Fig. 1). The latter consists of alternating TM and oxygen-carbon/boron layers, such that each TM ion is at the center of a distorted TMO 6 octahedra. The TM magnetic moments of ACO 3 (FeBO 3 ) display an analogous antiferromagnetic (AFM) order at low (room) temperature: the moments lie in the crystal ab basal plane, and are coupled ferromagnetically in each TM layer and antiferromagnetically between adjacent layers. The moments in different layers, however, are not exactly antiparallel one to another: the finite DMI causes the moments to be slightly canted and results in a small net magnetization in the basal plane of the crystal 1,2 (Fig. 1). The single crystals used in the present investigation are the same as Ref. 1 and 2, which the reader is referred to for further details on the crystal and magnetic structures and the sample growth.

A. First-principles calculations
The orbital and spin moments of the selected compounds were calculated using the Vienna ab initio simulation package (VASP) 23,24 within the local density approximation taking into account the on-site Coulomb interaction U and SOC (LDA + U + SO) 25 . The calculations are the same as outlined in our recent resonant scattering investigation 1 , where further details on the calculation methods can be found. The initial magnetisation directions were set to lie along the x direction, with x perpendicular to a 2-fold axis and contained in a c glide plane of the R3c structure. This results in having a canted AFM state, which is the lowest-energy state for all compounds. The results will be discussed and compared to the experiment in § IV. Values of spin and orbital moments reported in the present work are projections of the magnetisation density onto a sphere around the corresponding TM ion. Due to covalent bonding of the TM 3d orbitals with the oxygens 2p states, part of the magnetisation density appears on the ligand sites. The latter also contributes to the net magnetic moment. as derived from DFT calculations and measured by means of NXMS. The xyz reference frame is defined such that x is perpendicular to a 2-fold axis and contained in a c glide plane of the R3c structure and z is parallel to the crystallographic c axis.

B. Non-resonant X-ray Magnetic Scattering
The NXMS 26,27 measurements were performed in vertical scattering geometry at beamline I16 of the Diamond Light Source, Didcot UK 28 . The crystals were mounted on the sample rotational stages of the 6-circle kappa diffractometer with the c axis of the R3c trigonal structure aligned vertically. A standard closed-cycle cryostat was used to cool down the ACO 3 samples below the Néel transition temperature of the canted AFM structure, while the data on FeBO 3 were collected at room temperature. A magnetic field µ 0 H ≈ 35 mT, sufficient to drive the canted AFM structure into a single-domain phase 29,30 , was applied to the ab plane of the crystal using the rotating permanent magnet setup already successfully employed in our recent resonant x-ray scattering measurements 1,2 .
The diffracted signal arising from several space-group forbidden reflections of the type (00L), L = 6n + 3 and (HHL), L = 2n + 1 was measured using linearly polarized radiation (30 × 200 µm 2 spot size), with the electric field vector of the incident x-rays lying in the horizontal plane (referred to as σ polarization following the conventions by Blume and Gibbs 27 ). For each reflection, the scattered signal for both the rotated (σ − π ) and unrotated (σ −σ ) polarization channels was measured as a function of the magnetic field direction (where π and σ denote the polarization of a scattered beam whose electric field vector is perpendicular or parallel to the scattering plane, respectively). The latter is described by the angle η: this is defined such that the field lies in the vertical scattering plane (pointing towards the detector) for η = 0 • and is perpendicular to the latter for η = 90, 270 • .
Polarization analysis of the scattered beam was achieved by means of the (004) reflection from a pyrolytic graphite (PG) single crystal [with the exception of the data presented in § V B, for which the (006) was used]. The total scattered intensity without polarization anal-ysis was also measured for MnCO 3 , FeBO 3 and CoCO 3 in order to correct for the different reflection efficiencies of the PG crystal in the σ − π and σ − σ polarization channels (see Appendix B). The energy of the incident beam was kept fixed to E = 5.223 keV (E = 7.684 keV for the data of § V B), chosen for being away from any sample absorption edges and for minimizing the crosstalk between the two orthogonal light polarizations. For most of the measured reflections, equivalent data sets were collected at several different sample azimuths ψ 30 , whose values were selected to minimize the contribution of multiple scattering to the diffracted intensity. The corresponding results were then averaged. All the ψ values reported in this paper are defined with respect to the (100) azimuthal reference [ψ = 0 • when the (100) reciprocal direction lies in the scattering plane pointing towards the detector].
C. X-ray Magnetic Circular Dichroism XMCD was measured on a single crystal of CoCO 3 at the high-field magnet end station (BLADE) of beamline I10 of the Diamond Light Source. A thin film of Pt (≈ 2 nm) was deposited via sputter coating on the crystal's facet orthogonal to the c axis at the Research Complex at Harwell (Didcot, UK) prior to the XMCD measurements. The purpose of the Pt coating was to create an electrical contact with the illuminated area of the sample. The latter allowed the drain current of photo electrons to be extracted thus making total electron yield (TEY) detection possible despite the strong insulating character of CoCO 3 . The crystal was clamped on a electricallygrounded copper holder and inserted in the UHV sample environment of the I10 superconducting magnet with the coated surface facing the incident beam.
The measurements were performed at a shallow (20 • ) incident x-ray angle, so that the external magnetic field, (Color online) Representative magnetic reflections dependence on the magnetic field direction for two different polarization states of the diffracted x-ray beam in CoCO3. The data points represent the diffracted intensity integrated over a rocking scan while the solid curves correspond to the best fit to Eq. (4). The data were collected at T = 5 − 6 K and ψ = 83 • , 0 • , 108 • , 30 • for the (003), (107), (117) and (009) reflection, respectively. A small constant background originating from residual multiple scattering has been removed from all the data sets. For each reflection, the intensity is normalized to the maximum value across both polarization channels. A schematic drawing of the vertical scattering geometry, along with the definition of the u1u2u3 reference frame used to express the cross sections of Eq. (1) and the magnetic field angle η, is reported in the top panel.
directed along the incident beam wave vector, was almost perpendicular to the c axis. A relatively small field value (µ 0 H = 0.4 T) was used: this was chosen to be sufficiently large to suppress the magnetic domains structure and thus generate a net magnetization along the field while being, at the same time, small enough not to significantly perturb the in-plane canted AFM order of the Co 2+ moments. X-ray Absorption Spectroscopy (XAS) measurements were collected across the Co L 3 (778.2 eV) and L 2 (793.1 eV) edges for opposite helicities of the incident circularly-polarized soft x-ray beam (20 × 100 µm 2 spot size) and opposite directions of the external field. Several XAS spectra were collected and averaged for each permutation of light polarization and field direction and the resulting spectra combined to obtain the XMCD signal. The orbital and spin angular momenta derived from the DFT calculations are summarized in Table I. The calculations predict a negligible orbital contribution to the total angular momentum of the Mn 2+ and Fe 3+ ions. On the other hand, a significant orbital angular momentum is found in CoCO 3 and NiCO 3 . The orbital angular momentum is particularly large for the Co 2+ ion, where it reaches almost 60% of the spin value. This peculiar trend does not find a trivial explanation in a simple isolated-ion picture. Although Hund's coupling applied to Mn 2+ and Fe 3+ (3d 5 , l = 0, s = 5/2) predicts a zero orbital moment, the orbital contribution should be larger in Ni 2+ (3d 8 , l = 3, s = 1) than in Co 2+ (3d 7 , l = 3, s = 3/2). Moreover, despite the nominal 3+ oxidation state of the magnetic ion in FeBO 3 , the calculations predict a covalent, rather than ionic, character for the Fe-O bond: this results in an electronic configuration close to 3d 6 (l = 2, s = 2), which is then expected to host a finite orbital moment.
In order to verify the theoretical predictions, we combined the well-established polarization dependence of NXMS 27 with our novel rotating magnet technique 1,2 . The diffracted intensity arising from the long-range AFM order of the A(C,B)O 3 compounds was probed at several space-group forbidden reflections below the Néel transition temperature. A complete summary of the reflections measured for each compound can be found in the Supplemental Material 30 along with the relevant experimental parameters. The scattered signal is purely magnetic in origin, as proved by the fact that the signal vanishes upon warming above T N following the expected critical behaviour as a function of temperature 30 . The only exception is represented by the (105) and (207) reflections in CoCO 3 , which will be discussed in §V C. For each reflection, the signal was measured in both the σ − σ and σ − π channels as a function of a 360 • rotation of the external field in the basal plane of the crystal. The canted AFM structure rotates in response to the application of the field 1,2 : the corresponding magnetic field dependence of the scattered intensity can then be exploited to extract the relative orbital and spin contribution to the magnetic moment, as described in the following.
Given an incident σ-polarized x-ray beam, the NXMS amplitudes (neglecting a constant imaginary pre-factor) for σ and π -polarized scattering read as follows 27 : (1) where θ is the Bragg angle of the measured (HKL) reflection and L i and S i are the components of the orbital (L) and spin (S) structure factors in the u 1 u 2 u 3 reference frame defined in Ref. 27, respectively. As shown in the schematic of Fig. 2, the latter is defined such that u 3 is antiparallel to the scattering vector Q = k out − k in , u 1 lies in the scattering plane and points towards the detector and u 2 is orthogonal to the scattering plane. The magnetic structure factors represent the Fourier transforms of the orbital and spin magnetization densities and thus directly depend on the direction of the magnetic moments. In the specific case of the magnetic reflections under study, they are given by (see Appendix A for a detailed derivation): ] is the i-th component of the magnetic moment (expressed as unit vector) of the A (B) sublattice along the u i direction of the u 1 u 2 u 3 frame. C L and C S are constants terms, whose ratio depends on the relative magnitude of the orbital (l) and spin (s) angular momenta and the orbital [f l (Q)] and spin [f s (Q)] form factors: where Q is the modulus of the momentum transfer asso- . The data were collected at T = 3 K using a magnetic field µ0H = 0.4 T applied in the ab plane of the crystal. The dashed grey line represents the integrated XMCD signal used for the application of the sum rules: the l/s value refers to corresponding orbital-to-spin angular momenta ratio. The XAS data are normalized such that the post-edge spectral weight is equal to unity.
In the case of negligible magnetocrystalline anisotropy (see § V B for the case when this assumption no longer holds), the sumμ A +μ B of the moments of the two sublattices aligns along the direction of the rotating external field H(η): the differenceμ A −μ B of Eq. (2), perpendicular to the field, is forced to follow and causes the scattering amplitudes to vary accordingly. After inserting Eq. (3) and (2) into Eq. (1), the corresponding diffracted intensities are described by the following relations: where the dependence of the magnetic moments differences (μ ) on the field angle η has been emphasized. The momentum-dependent orbital-to-spin ratio l(Q)/s(Q) can be extracted through a fit to Eq. (4) of the measured dependence of the diffracted intensities on the magnetic field direction in the σ − π and σ − σ polarization channels. Data for two different light polarizations are needed due to the arbitrary scale factor relating the modulus square of the scattering amplitudes to the measured intensities values.
Representative data measured in CoCO 3 for four different magnetic reflections are displayed in Fig. 2 along with the best fits to Eq. (4). The data were collected by measuring the integrated intensity of the diffraction peak over a rocking scan of the sample at each value of the magnetic field angle η. The magnetic intensity displays very well defined 180 • -periodic sinusoidal oscillations, which are out of phase in the two polarization channels. For any given reflection, the l(Q)/s(Q) ratio is encoded in the relative amplitude of the σ − σ and σ − π intensity modulations. It should be noted that the measured σ − σ / σ − π amplitude ratio is also subject to a trivial azimuth-dependent geometrical factor related to the components of the magnetic moments in the u 1 u 2 u 3 reference frame [i.e. the quantities (µ (4)]: this explains why the two symmetricallyequivalent (107) and (117) reflections show significantly different amplitudes in Fig. 2 despite being characterized by the same value for the orbital-to-spin ratio. The data are of extremely high quality since the magnetic field measurements of Fig. 2 are performed without moving the sample. Therefore, variations in the self-absorption and grain hopping caused by the sphere of confusion of the diffractometer, which affect more conventional azimuthal scans, are not present in this case.  Table I along with the corresponding values from DFT calculations.
There is generally a very good agreement between the measurements and the calculations. Importantly, the trend across the series of compounds is confirmed: while MnCO 3 and FeBO 3 behave as almost pure spin systems, a significant unquenched orbital moment is found for both CoCO 3 and NiCO 3 . In particular, the predicted large value of the orbital moment in CoCO 3 is confirmed. This is somewhat consistent with previous studies on crystals 14,15,18 and thin films 31 of CoO, where a large orbital moment was also found. As we will show in § V, the presence of a large orbital moment is confirmed by XMCD measurements and results in the emergence of several interesting phenomena in the physics of CoCO 3 . When comparing the measurements with the calculations, it is important to bear in mind that due to the covalent bonding of the TM 3d orbitals with the oxygen 2p states, part of the magnetisation density appears on the ligand sites. While the NXMS measurements are sensitive to both, the DFT calculations neglect the oxygen contribution. This could explain, for instance, the partial discrepancy (still within the experimental uncertainty) between the measured and calculated value for CoCO 3 . This seems to be confirmed by the XMCD measurements outlined in § V A, which probes selectively the magnetization of the TM ion.

V. THE PECULIAR CASE OF CoCO3
A. XMCD investigation of the orbital moment in CoCO3 One of the main results of the previous section is the presence of an unusually large unquenched orbital moment in CoCO 3 . In order to confirm this finding, we performed XMCD measurements at the Co L edges as described in § III C. The results are shown in Fig. 4, where the absorption spectra obtained by combining the field and polarization reversal measurements are plotted along with the corresponding dichroism. The presence of a significant unquenched orbital moment is immediately evident from the much larger XMCD signal at the Co L 3 edge compared to L 2 . The application of the sum rules for the spin (µ s ) 32 and orbital (µ l ) 33 magnetic moment to the integrated XMCD signal shown in Fig. 4 leads to a value of the orbital-to-spin ratio l/s = 2µ l /µ s = 0.5(1) which is also confirmed by Co L-edge multiplet simulations (l/s ≈ 0.6) as described in § V C 2. This value confirms the one derived from our NXMS measurements within the experimental uncertainty, thus further consolidating our findings. One could argue that the nominal value of the l/s ratio found by means of XMCD is closer to the calculated one (Table I), which neglects the oxygen contribution to the total magnetization density. This is perfectly consistent with the resonant nature of the absorption process which, contrary to NXMS, selectively probes the magnetization density localised on the Co 2+ ions.

B. Single-ion anisotropy
One of the consequences of the large orbital moment in CoCO 3 is the presence of a large (≈ 7 meV/Co 2+ ) out-ofplane magnetocrystalline anisotropy, which confines the magnetic moments in the ab plane of the crystal. The magnetic anisotropy within the ab plane is expected to be significantly smaller. However, as we will show hereafter, its effect on the magnetic field dependence of the scattered intensity is clearly visible. For a crystal of space group R3c the single-ion anisotropy in the ab plane is described by the following energy cost per unit volume 35 : where K 0 (T ) and K(T ) are temperature-dependent constants (in energy per unit volume) which define the strength of the anisotropy and α is the angle describing the magnetization direction with respect to the crystal axes: this is defined such that the net magnetization resulting from the moments canting (which we shall refer to simply with the term "magnetization" from now on) is orthogonal to the [100] crystallographic direction for α = 0 • (see inset in Fig. 5). This 6-fold energy term is minimized for α = 30 • + n60 • (n integer index) and thus defines three main easy magnetization axes along the [100], [110] and [010] crystallographic directions.
In the presence of an external magnetic field H(β) the total energy per unit volume can be written as [dropping the constant K 0 in Eq. (5)]: where M (T ) is the temperature-dependent magnitude of the magnetization. Analogous to the magnetization angle α, β defines the direction of the external magnetic field with respect to the orthogonal to the [100] direction. This is related to the angle η used to express the magnetic field dependence of the scattered intensity through the relation β = ψ − η + 60 • , where ψ is the sample azimuth and the 60 • offset is simply due to the initial magnet position with respect to the crystal axes. In Eq. (5), E f ield represents the Zeeman interaction of the net magnetic moment with the external field. For the case of negligible anisotropy considered in § IV (E anis. ≈ 0), the Zeeman term forces the magnetization to align parallel to the applied field (α = β). In the general case of non-negligible anisotropy, however, M and H will lie along different directions. For any given direction β of the external field, the equilibrium direction α of the magnetization is obtained by minimizing the Hamiltonian of Eq. (6): regardless of the direction of the field and jumps discontinuously from one easy axes to another as the field rotates. A non-trivial 6-fold periodic function is obtained in the intermediate regime.
The solutions α(β) can be used to calculate the magnetic structure factors (2) and simulate the magnetic scattering intensities of Eq. (4) for different h values. The simulations are reported in Fig. 6 for the (009) σ − π intensity and three representative h values. For negligible anisotropy (large h) a smooth sinusoidal oscillation is obtained, analogous to the data shown in Fig. 2. As the anisotropy increases (h decreases) the intensity modulation takes on a peculiar "shark-fin" shape and eventually becomes discontinuous. Exactly the same trend is seen in the measured data as a function of temperature shown in Fig. 7. Increasing the temperature towards the Néel transition has, in this case, the effect of weakening the magnetocrystalline anisotropy: upon warming, the shark-fin shape progressively disappears and symmetric sinusoidal oscillations are recovered close to T N .
The scattered intensity calculated from the solutions of Eq. (7) can be used to fit the experimental data of Fig. 7 leaving the anisotropy constant K(T ) as a free fitting parameter and using the magnetization values of Ref. 36 in the parameter h. The calculations reproduce extremely well the measurements, as shown in the two representative fits of Fig. 8 for data collected well below and close to the magnetic transition. The values of K(T ) obtained from the fits are displayed in Fig. 9. As expected, the basal plane anisotropy constant decreases with increasing temperature: a quadratic dependence of the type ∝ (T N − T ) 2 is observed over most of the temperature range explored, similar to what was reported by Kaczer 34 . We find K = 11(2) neV/Co 2+ at T = 4 K, in agreement with a previous estimate 34 . Although almost 6 orders of magnitude smaller than the out-of-plane value, the effect of a finite basal plane anisotropy is clearly visible in the data and proves the extremely high sensitivity of our novel rotating magnetic technique to small interaction terms of the magnetic Hamiltonian.
The anisotropy-induced distortion at low temperature is not as evident in the data of Fig. 2 collected around T = 5 K. The two sets of data were measured during different experiments on two different crystals and several factors might explain the observed discrepancy. As well as differences in the crystal quality (crystal defects could, for instance, impact the field dependence of the scattered intensity), a significant beam heating has been observed in several occasions and might have also played a role despite the precautions taken to minimize it. The latter is very sensitive to the exact experimental conditions (which were different for the two data sets), such as sample mounting and incident flux: a sample temperature just a few degrees higher than the nominal value could explain the apparent lack of anisotropy in Fig. 2. Moreover, unlike the data presented in Fig. 2, the measurements of Fig. 7 were not collected integrating the intensity over a rocking scan at each value of the field angle. This allowed a much larger number of data points to be collected in a significantly shorter time (minutes compared to hours): both the coarser sampling and the averaging of any long term drifts could be at the origin of the apparent lack of any significant shark-fin distortion. It should also be noted that the magnetic field value depends on the exact position of the permanent magnet used for the measurements. This was fixed on the diffractometer rotational stages manually: a slightly different position between the measurements of Fig. 2 and Fig. 7 is likely to play a role in the observed discrepancy.
C. Forbidden charge scattering

Experimental data and empirical model
One of the most striking manifestations of the large orbital moment in CoCO 3 consists in the appearance of a peculiar interference pattern in the magnetic field dependence of space-group forbidden reflections. As argued hereafter, this is the result of the presence of a subtle, extremely weak, contribution to the scattered intensity induced by the ordered magnetic moment below the Néel  207) intensity magnetic field direction dependence as a function of the sample azimuth ψ for two different polarization states of the diffracted x-ray beam. For each azimuth value and polarization channel, the measured intensity (represented through the color scale) as a function of the field direction is normalized to its maximum value. A detail of the ψ = 65 • data sets (dashed horizontal lines) is shown in Fig. 11. (c-d) Global fit of the (207) magnetic field dependence shown in (a-b). The σ − σ and σ − π data at each azimuth value have been fitted to Eq. (8) and the second of Eq. (4), respectively, as shown in Fig. 11 for ψ = 65 • . The ratio between the magnetic and forbidden charge amplitudes in σ − σ was kept constant across all ψ values, while the multiple scattering amplitude was left free to vary: the resulting values are plotted in Fig. 14. transition. The same effect has not been observed in the other compounds of the series 30 and is thus to be considered a distinctive aspect of the physics of CoCO 3 . Of all the space-group forbidden reflections that we have measured, this interference is clearly visible only for the (207) and (105). Here, as shown in Fig. 10(a) for the (207) reflection, the dependence of the scattered intensity in the σ − σ channel on the magnetic field direction displays abrupt variations with the sample azimuth. A similar effect is also seen by varying the energy of the incident x-rays 30 . This contrasts the scattered intensity in σ − π [ Fig. 10(b)], which, regardless of the ψ and energy value, exhibits the normal 2-fold oscillation seen for the other reflections (Fig 2).
A detail of the (207) magnetic field dependence at ψ = 65 • [dashed lines in Fig. 10(a-b)] is shown in Fig. 11. Magnetic scattering alone cannot account for the nearly 4-fold pattern observed in the measured intensity in σ − σ : an extra contribution, displayed as a green dashed line in Fig. 11(a), must be introduced. The latter is, in turn, the result of two interfering scattering amplitudes of charge origin: (i) a sinusoidally-oscillating forbidden scattering term, C F S σσ (η), plotted in Fig. 12 as a function of the field direction and (ii) a field-independent multiplescattering amplitude, C M S σσ . The total scattered intensity in the σ − σ channel is therefore given by: where I M agnetic σσ = |M σσ (η)| 2 is given by the first line of Eq. (4). The charge origin of C F S σσ (η) is suggested by the absence of the interference effect in the rotated polarization channel and further confirmed by its temperature dependence (Fig. 13): the latter exhibits a critical exponent twice as large as the magnetic one, as expected for magnetic-induced charge scattering 37 .
The forbidden amplitude stems from a peculiar distortion of the Co 2+ electron cloud induced by the magnetic moment in the magnetically ordered phase. While  Fig. 10(a-b)]. The symbols represent the measured diffracted intensity, while the solid line refer to the global fit as explained in the text. In (a) the dashed green line and the dash-dot light blue line correspond to the charge and magnetic contribution to the global fit, respectively. The data were collected at T = 4 K and are normalized to the peak intensity of the σ − σ channel. A constant background originating from multiple scattering has been removed from the σ − π data set. a microscopic description requires detailed calculations of the Co 2+ ground-state wave function (see § V C 2), most aspects of the resulting scattering process are captured by the simple "toy model" sketched in the drawing of Fig. 12. This assumes a small elongation of the Co 2+ electron cloud along the magnetic moment (µ) direction, which is modelled by artificially adding a pair of negative charges to either side of each Co 2+ ion along µ. The electron cloud distortion reduces the symmetry of the crystal in the magnetically ordered phase such that, for an arbitrary field direction, only the inversion centre is left [space group P1 (No.2)]. The two extra charges are set to rigidly follow the rotation of µ as this is dragged around by the external field, thus originating the fielddependent term C F S σσ (η) shown in Fig. 12. The latter interferes with the multiple scattering amplitude and gives rise to the observed magnetic field dependence.
The multiple scattering amplitude C M S σσ plays a key role in our observations. In particular, it explains the azimuthal dependence of the scattered intensity of Fig. 10(a). Contrary to standard Bragg diffraction (also referred to as two-wave diffraction), where the diffracted radiation originates from a single scattering event of the primary beam, in multiple-wave diffraction the secondary beam originated from the scattering of the incident xrays can act as a primary beam for a second scattering process, thus giving rise to a tertiary reflection 38 . This results in additional diffraction peaks, which can appear at nominally forbidden (HKL) values. Although much weaker than Bragg reflections, multiple scattering peaks can have a comparable intensity to the magnetic ones.
The condition for generating multiple-wave diffraction is much more stringent than in the two-wave case: as a result, the multiple scattering amplitude displays a strong dependence on the sample azimuth (see Fig. 14) as opposed to Bragg scattering, which does not depend on ψ. Moreover, while the latter does not change the polarization of the primary beam, multiple scattering can in general give rise to both σ and π -polarized radiation.
The ψ values at which multiple scattering occurs can be calculated 30,39 , and thus avoided, following a simple kinematical approach. Nonetheless, broad tails are also present away from the nominal scattering condition and generally result in a residual contribution in both polarization channels. Because of the inversion centre of the R3c space group, all structure factors, including the multiple scattering one, are real if one considers only Thomson scattering far from any absorption edge. Therefore, multiple scattering interferes with the "forbidden" amplitude (which turns out to be of Thomson nature and is thus real), but not with non-resonant magnetic scattering, which is out of phase by 90 • , and gives rise to the dramatic evolution with the sample azimuth reported in Fig. 10. In the absence of the forbidden amplitude, as is the case for the σ − π intensity of all reflections and the σ − σ one when the additional term C F S σσ (η) is absent (Fig. 2), multiple scattering simply results in a constant background superimposed to the intensity of magnetic origin. A significant multiple scattering background is responsible for the high-intensity streaks visible in the σ − π color map of Fig. 10(b).
The interference between the amplitude induced by the electron cloud distortion and multiple scattering repro-  Fig. 10(a-b), while the solid line represents calculations performed using a mixed kinematical/dynamical approach where the standard kinematical structure factors of the secondary and tertiary reflection are weighted by terms calculated in a dynamical framework.
duces extremely well our data: this is clearly shown by the color map of Fig. 10c, which displays the fit of the measured intensity of Fig. 10(a) to Eq. (8). The fit of the σ − π intensity [ Fig. 10(d)] was performed using the second line of Eq. (4), analogous to § IV. An arbitrary positive scale factor, constant throughout all ψ values, was used for the forbidden charge amplitude of Fig. 12 in the fit of the σ − σ intensity. This is because the empirical model is not capable of reproducing a physicallymeaningful value of the scattering amplitude. On the other hand, the phase of the oscillations (including its sign) is correctly predicted. This is elegantly proved by the values of the multiple scattering amplitude C M S σσ extracted from the fits of Fig. 10(c-d), which are reported in Fig. 14  x is parallel to the a crystallographic axis and z is parallel to the crystallographic c axis. dard kinematical structure factors of the secondary and tertiary reflection were weighted by terms calculated in a dynamical framework. The ψ dependence of the measured amplitude, in particular its sign, is remarkably consistent with the calculations, thus confirming the correctness of the forbidden amplitude phase. The empirical model also grasps other significant features of the forbidden amplitude that are confirmed by our form factor calculations ( § V C 2). In particular, it correctly predicts a vanishing amplitude when (i) the canting of the Co 2+ magnetic moments is set to 0 (perfect AFM alignment) or (ii) a specular (00L) reflection is considered. Within this simple model, a non-vanishing forbidden amplitude is expected to be present also for the two equivalent reflections (107), (117). However, this term appears to be much smaller than the magnetic contribution 30 and is not clearly visible in the measured data. The latter are well described by the magnetic scattering cross sections alone.

Microscopic model and role of SOC
In order to achieve a microscopic understanding of the electron cloud distortion induced by the ordered moment, we derived the 3d electron ground-state wave function for different directions of the external field by means of mul-tiplet calculations. The latter were carried out using a Hartree-Fock method in the mean-field approximation, using the code RCN 40 for the radial part of the Co 2+ wave function and Quanty 41 for the angular part. The ground state was computed separately for two clusters of Co 2+ ions: the latter correspond to the two unique orientations of CoO 6 octahedra of the R3c crystal struc- 3d−3d ), (ii) crystal field (10D q , D σ ), (iii) SOC, (iv) magnetic exchange (H ex ) and (v) Zeeman term of interaction with the external field. H ex is a mean-field term which mimics the effect of the field produced by the ordered moments. The values of the main parameters used for the calculations are summarized in Table II: these were obtained refining the initial atomic values to reproduce the XMCD spectra of Fig. 4 (see Supplemental Material 30 ). The corresponding ground-state Hamiltonian expectation values are also summarized in Table II: despite the absolute  values of the spin and orbital moments are somewhat  different from the ones reported in Table I 30 , a large value of the orbital contribution (l/s ≈ 0.6) is confirmed.
Moreover, as well as reproducing the absorption spectra, the electronic structure of the Co 2+ ion thus calculated is consistent with the one presented in Ref. 42 and reproduces the experimental XMCD spectra well 30 . Further information on the crystal field parameters and additional details on the calculations are reported in the Supplemental Material 30 .
The results of the calculations confirm that, as predicted by our empirical model, the charge density of the valence 3d electrons depends on the magnetic moment orientation. This is shown in Fig. 15, where a real-space representation of the charge density is reported for different directions of the external field. The Co 2+ groundstate wave functions for different field directions can then be used to calculate the corresponding atomic form factor: the resulting (207) and (105) scattering amplitudes show a sinusoidal magnetic field dependence analogous to the one of Fig. 12. Consistent with the empirical model, the amplitude vanishes for specular (00L) reflections and when the magnetic moment canting angle is set to 0. Most importantly, the multiplet calculations show that the amplitude also vanishes when the SOC is artificially switched off 30 . This attributes the magnetic-momentinduced distortion of the Co 2+ electron cloud to the coupling between lattice and magnetic degrees of freedom driven by SOC and further highlights the fundamental role played by the large unquenched orbital moment in the physics of CoCO 3 . . The deformation of the in-plane lattice parameters resulting from the fit is shown in Fig. 17.

D. Magneto-striction
Another distinctive evidence of the elongation of the Co 2+ electron cloud along the magnetic moment direc- tion is the resulting expansion of the unit cell in-plane lattice parameters. This is revealed by the angular shift of the Bragg peak of symmetry-allowed reflections as a function of the field direction in the magnetically-ordered phase (Fig. 16). Given an expansion ∆l > 0 of the unit cell along µ, the Bragg angle θ magnetic field dependence is correctly described by the following field-dependent lattice parameters distortion: where a 0 (b 0 ) is the value of the lattice parameter a (b) when the magnetic field is orthogonal to the a (b) axis [α = 0 • (60 • )]. The magnitude of the unit cell distortion ∆l can be obtained by fitting the measured shift of the Bragg peak to the one calculated through the lattice parameters of Eq. (9). This is shown by the red solid line of Fig. 16. The resulting unit cell deformation along the a and b axes is plotted as function of the field direction in Fig. 17. The deformation amounts to ≈ 70 ppm at T = 5 K (which correspond to a change of ≈ 35 fm in the lattice parameters) and decreases upon warming towards T N following the same critical behaviour of magnetic scattering 30 . This further confirms the magnetostrictive origin of the Bragg peak oscillations of Fig. 16 and constitutes further evidence of the coupling between crystallographic and magnetic properties induced by the large unquenched orbital moment.

VI. CONCLUDING REMARKS
In conclusion, our combined DFT, NXMS and XMCD investigation of a series of isostructural weak ferromagnets led to the following findings: • A non-trivial evolution of the orbital contribution to the magnetic moment with the filling of the TM ion 3d orbitals is present across the series. In particular, the value of the orbital moment was found to be particularly large for CoCO 3 as confirmed by both NXMS and XMCD.
• In CoCO 3 , SOC couples the large orbital moment and the spin of the Co 2+ ion and results in a strong single-ion uniaxial anisotropy and a much smaller, although still clearly visible in our NXMS data, basal plane anisotropy.
• SOC is also responsible for a distortion of the Co 2+ 3d electron cloud in the magnetically order phase: the latter is evidenced by a sizeable magneto-striction and, more spectacularly, by the appearance of a forbidden scattering amplitude at space-group forbidden reflections.
Our results combined together highlight the importance of SOC in the physics of weak ferromagnets and show how, even in the case of 3d transition metal oxides, SOC can have a significant impact on the magnetic properties of the system whenever the orbital degrees of freedom are not quenched. Finally, our investigation also proves the ability of modern first-principles calculations to predict the properties of materials which exhibit magnetoelectric coupling, Skyrmion lattices and other non-collinear magnetic ordering.
this case, the magnetic form factors can be expressed as follows [43][44][45] : Here, j 0 and j 2 are radial integrals of the type j k nl (Q) = ∞ 0 R 2 nl (r)j k (Qr)r 2 dr where R nl (r) is the radial part of the magnetic ion wave function and j k (Qr) is the spherical Bessel function of order k. The radial integrals (espressed as a function of the normalized momentum tranfer s = Q/4π = sin θ/λ, being θ the Bragg angle of the magnetic reflection and λ the wavelength of the incident x-ray beam) can be approximated by the following 43 : j 0 (s) = Ae −as 2 + Be −bs 2 + Ce −cs 2 + D j 2 (s) = s 2 (A e −as 2 + B e −bs 2 + C e −cs 2 + D ) (A5) where the values of the coefficients are tabulated for different oxidation states of each element in Ref. 43.
Appendix B: Details on the data treatment As mentioned in § III B, as well as the intensity in the σ − σ and σ − π channels, the total NXMS intensity was also measured in order to correct for the different reflection efficiencies of the PG analyser crystal in the two polarization channels. This mainly originates from the different beam divergence in the vertical and horizontal plane: as a result, the intensity detected in σ − σ and σ − π will generally be different even in the case of an equal distribution of σ and π polarization. Using the total intensity I tot , the measured values of the intensity in σ − σ and σ − π can be corrected by introducing a compensation factor f , which is defined by the following relation: Here, f = C σπ C σσ is the ratio of the two arbitrary scale factors which link the intensity measured in the two polarization channels to the total one recorded through the area detector. In the ideal case in which the reflection efficiencies for the two polarization channels are equivalent, only one scale factor would be necessary, which corresponds to having f = 1. In practice, the compensation factor f can be extracted using Eq. (B1) to fit the total intensity as a function of a 360 • rotation of the magnetic field measured with the 2D detector from the σ − σ and σ − π magnetic field dependences. The measured data can then be corrected by multiplying the intensity in σ − π by f : finally, Eq. (4) can be used to fit the corrected values.