Intertwined Spin and Orbital Density Waves in MnP Uncovered by Resonant Soft X-ray Scattering

Unconventional superconductors are often characterized by numerous competing and even intertwined orders in their phase diagrams. In particular, the electronic nematic phases, which spontaneously break rotational symmetry and often simultaneously involve spin, charge and/or orbital orders, appear conspicuously in both the cuprate and iron-based superconductors. The fluctuations associated with these phases may provide the exotic pairing glue that underlies their high-temperature superconductivity. Helimagnet MnP, the first Mn-based superconductor under pressure, lacks high rotational symmetry. However our resonant soft X-ray scattering (RSXS) experiment discovers novel helical orbital density wave (ODW) orders in this three-dimensional, low-symmetry system, and reveals intertwined ordering phenomena in unprecedented detail. In particular, a ODW forms with half the period of the spin order and fully develops slightly above the spin ordering temperature, their domains develop simultaneously, yet the spin order domains are larger than those of the ODW, and they cooperatively produce another ODW with 1/3 the period of the spin order. These observations provide a comprehensive picture of the intricate interplay between spin and orbital orders in correlated materials, and they suggest that nematic-like physics ubiquitously exists beyond two-dimensional and high-symmetry systems, and the superconducting mechanism of MnP is likely analogous to those of cuprate and iron-based superconductors.

Unconventional superconductors are often characterized by numerous competing and even intertwined orders in their phase diagrams. In particular, the electronic nematic phases, which spontaneously break rotational symmetry and often simultaneously involve spin, charge and/or orbital orders, appear conspicuously in both the cuprate and iron-based superconductors. The fluctuations associated with these phases may provide the exotic pairing glue that underlies their high-temperature superconductivity. Helimagnet MnP, the first Mn-based superconductor under pressure [1][2][3], lacks high rotational symmetry. However our resonant soft X-ray scattering (RSXS) experiment discovers novel helical orbital density wave (ODW) orders in this three-dimensional, low-symmetry system, and reveals intertwined ordering phenomena in unprecedented detail. In particular, a ODW forms with half the period of the spin order and fully develops slightly above the spin ordering temperature; their domains develop simultaneously, yet the spin order domains are larger than those of the ODW; and they cooperatively produce another ODW with 1/3 the period of the spin order. These observations provide a comprehensive picture of the intricate interplay between spin and orbital orders in correlated materials, and they suggest that nematic-like physics ubiquitously exists beyond two-dimensional and high-symmetry systems, and the superconducting mechanism of MnP is likely analogous to those of cuprate and iron-based superconductors.
In strongly correlated electron systems, the kinetic and interaction energies of electrons compete, which often couples the charge, spin, and orbital degrees of freedoms, resulting in a variety of complex quantum phases. Nematic order is one of the most important forms of intertwined order, which exists in both the two known families of high temperature superconductors -cuprates and iron-based superconductors. In iron pnictides, e.g. LaFeAsO and BaFe 2 As 2 [4][5][6], the nematic orbital order often emerges at a higher temperature than the collinear antiferromagnetic order, and with half the period, indicating that the nematicity may be strongly tied to magnetic fluctuations. In the stripe phase of the cuprates, the charge order forms with half the period of the spin order [7,8]. Moreover, unlike magnetic order, the nematic order is time-reversal invariant, allowing it to directly couple to orders arising from charge and orbital degrees of freedom. Furthermore, it also appears that superconductivity is maximized in the region with the strongest nematic fluctuations. In the iron pnictides, for example, the orbital and/or spin fluctuations associated with the nematic order are suggested to promote superconductivity [9,10]. All these essential phenomena suggest that the electronic nematic phase stands as a bridge between high temperature superconductivity and magnetism.
Under pressure, the helical magnetic compound CrAs and MnP have been recently found to be the first Crand Mn-based superconductor, respectively [1,11]. For example, under ambient pressure, MnP first enters a ferromagnetic state at T C = 290 K, then a metamagnetic transition at T S = 50 K switches it into a double helical magnetic state with moments lying in the ab plane ( Fig.  1(a, b)). Such a complex magnetic behavior suggests that MnP is most likely an unconventional superconductor like the cuprate and iron pnictides. However, not only has the material a helical magnetic state that differs from those of cuprates and iron-based superconductors ( Fig. 1(b)), but its lattice structure is three dimensional and does not possess any high symmetry rotational axis. Moreover, the spin, charge and orbital degrees of free- -pol. The reference spectra were obtained on beam line 8.0.3 at the Advanced Light Source (ALS) [15] and are plotted here offset by 1.55 eV in order to match our MnP spectrum obtained at the SSRL. The Mn L-edge XAS of MnP is characteristic of Mn 2+ ions. (d) Scattering peaks around q1=(0 0 0.116), q2=(0 0 0.234), and q3=(0 0 0.352) at resonant energies (black dots) and non-resonant energy (red dots, E=655 eV) at 20 K. q1 and q2 signals were measured by a Schottky barrier photodiode, and the q3 signals were measured by a quantum-efficient channeltron. The three resonant peaks represent the strongest resonant signals at each q position (see Fig. 3). (e) Schematic diagram of the location of the diffraction peaks in MnP including lattice peaks (red stars), magnetic peaks (blue dots), and the newly discovered resonance peaks in this paper (green and magenta dots).
dom in MnP could be highly interconnected. Therefore, MnP provides a novel playground to probe the relation between unconventional superconductivity, magnetism, and possibly nematicity in the vicinity of a helical spin order on a low-symmetry lattice.
Resonant soft X-ray scattering (RSXS), which can be viewed as a combination of x-ray absorption (XAS) and x-ray emission (XES) spectroscopies with x-ray scattering, provides a direct and powerful probe of the ordering in 3d transition metal compounds [12][13][14]. To elucidate the electronic state of MnP, we first present its XAS spectrum at Mn L-edge in Fig. 1(c), together with that of MnO as a fingerprint of the Mn 2+ valence state [15], Mn 2 O 3 as a typical Mn 3+ spectrum, and MnO 2 representing Mn 4+ . Evidently, MnP has the typical ab-sorption spectrum of a Mn 2+ state with the half filling high spin configuration. To our knowledge, this is the first XAS measurement on MnP. The Mn 2+ state in MnP is contrary to the previously speculated Mn 3+ valence in this material [16]. The spin moment on each Mn site is only ∼ 1.3µ B [17]; this value is significantly reduced from the localized spin-only moment ∼ 5µ B for the 3d 5 configuration in Mn 2+ by Hund's rule. This low spin moment may arise from quantum fluctuations and hybridization [18,19], as is the case for the iron-based superconductors [20,21]. It should be noted that for a half-filled electronic system, the orbital moment usually should be quenched, thus prohibiting spin-orbital coupling. However, orbital anisotropy/ordering have been found in several half-filled systems due to hybridization

T (K)
Wavevector (r.l.u.)  resonance peaks) or Lorentz (for q3 resonance peaks) fittings with linear background. The shadowed area in the upper panel of (c) is to indicate the weak peak intensity at the highest measured temperature. (e) Temperature dependence of the peak area of q1, q2 and q3 near the helical transition, revealing their differing transition temperatures. The data were obtained in the cooling processes. Solid lines are guides to the eye. (f) Temperature dependences of the q1, q2/2, and q3/3. Here, we used a photodiode detector for q1, a channeltron detector for q2, and a Greateyes CCD for q3. (g) and (h) are the temperature dependences of the peak area and 1/FWHM for q1 (black circles) q2 (red squares) resonance peaks in the cooling measurements, respectively. Log scale were used for the vertical axes and the solid lines are guides to the eye.
anisotropy [22] and vacancy modulation [23]. Therefore, orbital angular momentum may be partially restored. Figure 1(b) shows the scattering geometry, where bc was used as the scattering plane. The incident X-ray is either horizontally (π) or vertically (σ) polarized. In this configuration, σ polarization has electric field along the a axis and momentum transfer along the (0 0 L) direction. When helical magnetic order emerges below T S = 50 K, its spin moments (blue arrows) lie in the ab plane with magnetic propagation wavevector (0 0 0.117) [24]. In order to avoid specular reflection, we used a MnP single crystal with a (101) surface and tilted the sample by 48.4 • in our RSXS experiments.
We made an exhaustive search for scattering signals at the Mn L-edge along (0 0 L) at 20 K. By varying L, X-ray polarization, photon energy, and detection mode, we discovered resonance peaks at q 1 =(0 0 0.11634±0.00002), q 2 =(0 0 0.23470±0.00003), and q 3 =(0 0 0.35239±0.00004). q 2 and q 3 wavevectors approximately double and triple that of q 1 , respectively. Figure 1(d) shows the scans for the strongest resonance peaks at q 1 (E = 637 eV, π polarization, black dots), q 2 (E = 648.2 eV, σ polarization, black dots), q 3 (E = 641.2 eV, π polarization, black dots), and the corresponding non-resonant scans at 655 eV (red dots). q 1 is perfectly consistent with the magnetic diffraction peak q m = (0 0 δ) in MnP, thus should correspond to helical magnetism. It should be noted that the magnetic diffraction peaks of MnP only appear at [H K L ± δ] positions in which [H K L] represent the Bragg peaks from the lat-tice, a law universal to helical magnets that has been verified by neutron diffraction in MnP [25], FeP [26], and CrAs [27]. It is thus clear that the new diffraction peaks at q 2 and q 3 discovered by RSXS are not from magnetism. The diffraction positions of lattice, magnetic, and our newly discovered peaks in reciprocal space are illustrated in Fig. 1(e).
Next, we study the temperature dependences of the diffraction peaks. First, we take the (0 0 2) lattice diffraction peak as a reference. It does not show any observable temperature dependence across T S and down to 20 K ( Fig. 2(a)), indicating no structure transition in this temperature range. In contrast, q 1 , q 2 , and q 3 peaks exhibit drastic temperature dependences. Figures 2(b-d) show the three diffraction peaks taken with their corresponding resonant scattering conditions, at q 1 (E = 637 eV, π polarization), q 2 (E = 648.8 eV, π polarization), and q 3 (E = 636.6 eV, π polarization) at different temperatures in a cooling sequence, respectively. As can be seen, not only does the peak intensity rapidly grow across T S , but the propagation wave vector moves to higher q with decreasing temperature. The peak areas of the q 1 , q 2 and q 3 diffraction peaks are plotted as a function of temperature in Fig. 2(e), all showing a jump just around 50 K, consistent with a first order transition (the hysteresis behavior of the q 1 peak can be found in Supplementary Fig. S4(a)). Interestingly, the q 1 peak is fully developed at a slightly lower temperature than the q 2 peak, while the full-development temperature of the q 3 peak lies in between or similar to the q 2 peak. This observation can be justified by more temperature dependence measurements ( Supplementary Fig. S4(b)). In Fig. 2(f), we show the wavevector evolution with temperature, in which q 2 and q 3 are divided by 2 and 3, respectively, in order to scale with q 1 . The propagation wavevectors of the three resonance peaks all show pronounced temperature dependence below T S , which is typical for an incommensurate electronic order. Throughout the measured temperature range, q 2 and q 3 are approximately at the 2q 1 and 3q 1 positions within the experimental accuracy, respectively, indicating that these electronic orders are interconnected with each other. According to the upper panels of Fig. 2(b-d), there are detectable scattering intensities above T S . The q 1 and q 2 peaks persist up to the highest measured temperatures, i.e., 58.5K and 54.5 K, respectively (upper panels of Fig. 2(b) and Fig. 2(c)). However, q 3 peak intensity is not detectable above 52 K (upper panel of Fig. 2(d)) due to its weak intensity and insufficient sensitivity of the detector. To reveal the detailed evolution at high temperatures, Fig. 2(g) plots the peak areas as a function of temperature in the cooling process in log scale, and Fig. 2(h) shows the temperature dependences of the correlation length, ξ i , which is defined as 1/FWHM, reflecting the average domain size along c. Above 50K, there is a similar slow-growing behavior for both the q 1 and q 2 peaks in Figs. 2(g) and 2(h). The peak intensities are two orders of magnitude lower than their full values, thus they are due to fluctuating orders, and the long tails into high temperatures are likely related to local strain distributions. A similar behavior has been observed for the nematic order in iron pnictides under uniaxial strain [5]. With decreasing temperature, ξ 1 and ξ 2 quickly increase almost identically before they saturate at T F 1 and T F 2 , respectively, indicating the peaks are from the same domain and have the same onset temperature, T E , which can be higher than 58.5K, the highest measured temperature. Below T F 2 , ξ 1 continues to increase until T F 1 , so it is larger than ξ 2 at low temperatures. Meanwhile, the peak intensities show sudden jumps to their fully-developed values just before T F 1 and T F 2 as well. Therefore, T F 1 , T F 2 and T F 3 are defined as the temperatures that the q 1 , q 2 and q 3 peaks are fully developed, respectively. Note that, between T F 2 and T F 1 , the q 1 peak intensity is low, indicating that the spin order is still fluctuating while the order corresponding to q 2 is already static. In addition, we found that ξ 3 is almost identical to ξ 2 . These behaviors contain rich information on the evolutions of the orders corresponding the diffraction peak, which will be discussed later.
To comprehensively elucidate the nature of the three resonant peaks at q 1 , q 2 , and q 3 , we have plotted in Fig. 3 the resonant profiles around the three q positions with T = 20 K, i. e., the scattering intensities as a function of reciprocal lattice (0 0 L), X-ray energy, and incident photon polarization. Since there was no polarization analyzer before the detector, the polarizations of the scattered photons were not distinguished. That is, using selfexplanatory subscript convention, the detected scattered intensities for two different incident photon polarizations I π = I ππ + I πσ I σ = I σπ + I σσ .
Clearly, the resonant profiles of q 1 , q 2 , and q 3 are very different from each other in maximal intensity, resonant energy, and polarization dependence. For example, the resonant energies for the maximum peaks at q 1 , q 2 , and q 3 are 637 eV (π polarization), 648.2 eV (σ polarization), and 641.2 eV (π polarization), respectively. The maximal intensity at q 3 is ∼ 23 times weaker than that at q 2 , and the maximal intensity for q 2 is ∼30 times weaker than that of q 1 , based on the data of the resonance peaks taken with the same detector. The q-integrated resonant profiles are compared in Fig.  4, as a function of energy and polarization. The resonant profile of the q 1 peak behaves similarly in both polarizations, consistent with its origin in resonant magnetic scattering. As illustrated by a simple analysis, magnetic scattering will result in ππ, πσ and σπ scatterings with comparable intensity [28,29]. The magnetic scattering matrix element is related to the dipole selection rules, which have the same origin as soft x-ray magnetic dichroism. This explains it having the strongest intensity amongst all three peaks. For the q 2 peak, I σ is about   (c-f) q2 and q3 resonance profiles measured by a channeltron detector, due to their relatively weak intensity. For q2, the Iσ maximum is ∼ 3 times of the Iπ maximum. For q3, the Iσ maximum is about 80% of the Iπ maximum. All data shown here are raw scattering intensity without background subtraction or absorption correction. The color bars indicates scattering intensity in arbitrary unit.
triple I π . Moreover, the energy positions of the q 2 profile differ between the σ and π polarizations, which is unlikely from magnetic scattering. The resonance profile of the q 3 peak shows moderate polarization dependence. The saturation of the q 1 peak intensity below T F 1 represents the full development of the incommensurate helical spin density wave. Remarkably, the weak q 2 peak is fully developed slightly above T F 1 , when the q 1 peak intensity is still two orders of magnitude smaller than its full value. It implies that this is not a simple second-order harmonic of the spin density wave but an indication of a hidden order which is induced by magnetic fluctuations. Applying Landau theory to the phase transition and using general symmetry analysis, we extract the nature of these orders by assuming that the Landau free-energy functional depends only on the fundamental Fourier components of the two different electronic orders. The q 1 peak corresponds to the double helical spin density wave, S q1,a ∝ (i, 1, 0) where a = 1, 2 denote the two helices. As q 2 ≈ 2q 1 , the hidden order, Π q2 , must be coupled to a quadratic term of S q1 . If we take Π q2 to be a scalar, a natural choice of Π q2 is the charge density wave, ρ q2 . The lowest order coupling in Landau theory can be written as H c = ab λ c,ab ρ q2 ( S * q1,a · S * q1,b ) + h.c.. However, such an unequal charge distribution at different Mn sites could not be found in our LDA calculations shown in the supplementary materials. Instead, an orbital distribution with half the periodicity of the magnetic order can be explicitly obtained in the calculation.
Here, we suggest that the hidden order must be a orbital density wave (ODW), which can be induced by the orbital redistribution in developing the double helical spin density wave. In a simple helically ordered state, it is well known that an orbital redistribution can be induced by spin-orbital coupling in the presence of the crystal field, and such an orbital redistribution only depends on | S q1 | [30], giving the X-ray magnetic linear dichroism effect (XMLD). In this case, spin-parallel and -antiparallel Mn atoms would have identical orbital (wave function) distribution, which explains why q 2 corresponds to an order with a period half that of the magnetic peak. This situation resembles the stripes in cuprates, but in MnP both orders are incommensurate with respect to the lattice. We observe resonance at the Mn 2p to 3d transition, and orbital order of 3d electrons gives an electronic orbital order, Π αα . Thus a ODW is simultaneously developed in the helical spin state. Now, in the case of double-helical spin order, we argue that a tiny orbital modulation can be induced by spin fluctuations even before the static helical spin order is developed. In general, the coupling between the ODW and the double-helical spin order in Landau theory can be written as With this coupling, there is a new phase in which < S α * q1,a >= 0 but < S α * q1,1 S α * q1,2 > = 0, which describes the locking of the magnetic fluctuations between the two helices. In general, this phase could exist slightly above T s , thus explaining the intriguing full development of the q 2 peak at a slightly higher temperature than the q 1 peak. The linear coupling between the ODW and S α * q1,1 S α * q1,a must result in < Π αα,q2 > = 0 if < S α * q1,1 S α * q1,2 > = 0. This argument is generally known as "order by disorder" and has been used to explain the nematicity in FeAs-based superconductors in which a similar phenomenon has been observed in the parent compounds, e.g. BaFe 2 As 2 [5,6], where the nematic orbital order emerges at a slightly higher temperature than the collinear spin order due to spin fluctuations.
Since the spins lie in the ab plane in the double helical phase, the orbital moment should also be in the ab plane from the coupling H Q . This explains the observed polarization dependence of the q 2 peak and also suggests that the charge redistribution mainly occurs in the d x 2 −y 2 and (a) The orbital density wave with 1/3 the period is represented by an exaggerated mixture of dxz and dyz orbitals that rotates three times as fast as the spin. (d) Cartoon showing that short-ranged and fluctuating SDW (represented by the blue area) and ODW domains (represented by the texture) share the same region in MnP for T E ≥ T≥ T F 2 . (e) When T F 2 ≥ T≥ T F 1 , the ODW is fully developed and stops to grow with decreased temperature, while SDW domains continue to grow. There is regions with fluctuating spin order but without orbital order. (f) When T ≤ T F 1 , the spin order is fully developed.
d xy orbitals (the orbital basis in the octahedral coordination with these two orbitals in the ab-plane is used here for the ease of explanation, while we note Mn ions are in a tilted and distorted octahedrons made of P ions).
The observation of the weak q 3 peak also lends strong support to the presence of hidden ODW order since it represents the harmonics generated by the coupling between the ODW and the helical spin order.
The intensity of the q 3 peak is less than 1/600 of that of the SDW peak, yet it saturates at a slightly higher temperature. This remarkable behavior suggests that the q 3 peak is also non-magnetic in nature. Moreover, since it is slightly stronger for π-polarized than σ-polarized incident photons (Fig. 4(c)), it should involve d xz , d yz , and possibly d 3z 2 −r 2 orbitals. The simplest possibility for the q 3 peak would the rearrangement or modulation of these orbitals induced by the coupling between the ODW (q 2 ) and the helical SDW (q 1 ) as a third-order effect. This represents another unique ODW that is observed here for the first time.
In Figs. 5(a-c), we summarize the observed orders. The helical spin density wave is shown in Fig. 5(a), represented by spins rotated by equally-spaced in-plane angles. In Fig. 5(b), the in-plane ODW follows the SDW with half the period, represented by a charge distribution whose axis follows the orientation of the spin. In Fig. 5(c), the out-of-plane ODW rotates with one third of the SDW period. These helical orbital density waves in MnP are discovered for the first time by RSXS, and they are intertwined with the helical spin density wave. It is particularly noteworthy that the temperature dependencies in Fig. 2 further illustrate the intricate relation between the spin and orbital orders. The short-ranged and fluctuating orbital and spin orders share the same domain as they both emerge and grow upon lowering the temperature (Fig.5(d)). This starts from quite high temperatures, likely related to strain distributions in the system. Upon further cooling, as shown in Fig. 5(e), the ODW order freezes and its domains are fixed after passing its fully-developed temperature T F 2 . Meanwhile, the spin order is still fluctuating and its domains keep expanding, since the intensity of the spin order peak is still a few percent of its full value at low temperature. As a result, there are regions with fluctuating spin order but no orbital order in this temperature range. Because the orbital order is fairly weak, as shown by the weak diffraction peak intensity, the orbital distribution can be influenced or pined by local strain or defects. Therefore, the orbitals in these regions are disordered, and its configuration could not follow the helical spin rotation. When temperature is further lowered below T F 1 , the spin order is fully developed and its domains stop to grow, as illustrated in Fig. 5(f). We note that because the out-ofplane ODW in Fig. 5(c) is the descendant of the SDW and the in-plane ODW, its domain is identical to that of the in-plane ODW.
A similar observation was made in Pr 0.6 Ca 0.4 MnO 3 , whose spin order correlation length is longer than its charge/orbital order correlation length, which was attributed to likely decoupling between spin and charge/orbital orders [31]. The diffraction peak intensity of its orbital order is also much weaker than the spin order in Pr 0.6 Ca 0.4 MnO 3 , although its spin ordering temperature is lower than that of the charge/orbital order. So the picture revealed in Figs. 5(d-f) through our detailed temperature dependence are likely ubiquitous for systems with both spin and orbital ordering, which explains the difference in the correlation lengths in Pr 0.6 Ca 0.4 MnO 3 .
Our findings provide an unprecedented picture on the intricate interplay between spin and orbital orders in MnP, and show that intertwined ordering and nematiclike ordering is ubiquitous to the phase diagrams of unconventional superconductors, even in the case of low symmetry and incommensurate ordering. The extraordinary helical ODWs found here may provide a foundation for understanding the complex behaviors of spin/orbital order in helimagnets and correlated systems in general, and for understanding the unconventional superconductivity in MnP and other related materials.
We After sealing in a quartz ampoule under vacuum, the following synthesis procedure was conducted: first raise to 650C in 10h, stay for 8h, and then raise to 1150C, stay for 6h, and then slowly ramp to 600C with a cooling speed of 3C /h. The liquid Sn flux was then immediately centrifuged. The filtered MnP crystals are needle-shaped and the typical crystal length is 2~7mm. The single crystals were characterized by Physical Properties Measurement System (PPMS) and X-ray diffraction (XRD). As shown in Fig. 1(a) in the main text, the magnetic susceptibility of MnP shows a sharp transition into the helical magnetic phase at 50 K and a ferromagnetic transition at 290 K, its resistance is also consistent with previous reports 1 . XRD pattern from the main plane of the MnP crystal used in experiments is shown in Fig. S1.

Resonant soft-X ray elastic scattering (RSXS) experiments.
RSXS and XAS experiments on a MnP single crystal at the Mn L-edge were performed on an in vacuo four-circle diffractometer at the BL13-3 station at the Stanford Synchrotron Radiation Lightsource (SSRL) and the REIXS beamline at the Canadian Light Source (CLS). In our experiment, the long direction of the needle-like MnP sample in RSXS experiment is along the b axis and the flat shiny surface is a (1 0 1) plane. We tilted the sample using a wedge with an angle of 48° in order to make the momentum transfer Q direction along the c axis and use (0 K L) as the horizontal scattering plane. In this way, the specular reflection background can be effectively suppressed. The experimental geometry is shown in Fig. S2. The beamlines are equipped with Elliptical Polarized Undulators (EPU) which can control the σ(π) incident photon polarization. The sample was aligned by the (0 0 2) lattice peak with X-ray energy at 2.18 keV and the lattice constants used here are a = 5.26 Ǻ, b = 3.17 Ǻ, and c = 5.90 Ǻ. The base experiment temperature is 20 K and the ultrahigh vacuum pressure was below 10 -9 torr. The momentum resolution of the RSXS instrument is roughly Q ~ 0.0005 Ǻ -1 , far smaller than the FWHM of the resonance peaks in MnP. The detectors used in the RSXS experiment are a Schottky barrier photodiode, a quantum-efficient channeltron, and an in-vacuum CCD (charge coupled device, Greateyes GmbH) camera. The sample's X-ray absorption spectra were measured by the total electron yield (TEY).

The charge/orbital density of MnP
Here we study the relationship between the magnetic period and the charge/orbital density period by density functional theory (DFT). Because the incommensurate ordering requires a large unit cell that is beyond our computing capability, a realistic calculation is currently unfeasible. Nevertheless, insight can be gained by using artificial commensurate lattice. Here, we design two antiferromagnetic (AF) states, one has the 1 × 1 × 2 supercell as shown in Fig. S3 (a) and the other has the 1 × 1 × 4 supercell as shown in Fig. S3 (b). We calculate the charge density difference in this two AF state. The charge density difference is defined as the difference between the selfconsistent charge density and the non-self-consistent charge density.
Our DFT calculations employ the projector augmented wave (PAW) method encoded in Vienna ab initio simulation package(VASP) [2][3][4] , and generalized-gradient approximation (GGA) 5 for the exchange correlation functional is used. Throughout this work, the cutoff energy of 600 eV are taken for expanding the wave functions into plane-wave basis. The number of these k points are 9 × 15 × 5 for 1 × 1 × 2 supercell case and 9 × 15 × 3 for 1 × 1 × 4 supercell case in the calculations of the magnetic structures.
In the designed model, MnP is in antiferromagnetic states with spin alignment along the c axis. The local magnetic moments are indicated by arrowhead in Fig. S3. So the magnetic unit cell are 1 × 1 × 2 supercell in Fig. S3(a) and 1 × 1 × 4 supercell in Fig. S3(b). In both cases, there is a rotated charge density distribution, indicative of orbital order. Its unit cell, as indicated with red frame, is half of the magnetic unit cell.
In summary, the period of orbital order is half of the period of the magnetic order in MnP. Fig. S3: Charge density difference maps of MnP in two antiferromagnetic state visualized by VESTA 6 . The 1 × 1 × 2 supercell case (a) and 1 × 1 × 4 supercell case (b). The isosurface electron density difference distribution of MnP is plotted at 0.02 e Bohr -3 . As shown by the rotated charge density distribution, there is an orbital order, and the red frame indicates its unit cell. Figure S4(a) shows the q1 resonance peak (E = 637 eV, π polarization) area as a function of temperature at the cooling (black dots) and warm up (red dots) measurements. There is a clear hysteresis loop between the cooling and warm up curves around TS, consistent with a first order transition.

More temperature dependence data of diffraction peaks
More data on the temperature dependence of the three resonance peaks, i.e., peak area as a function of temperature are shown in Fig. S4(b). By over-plotting several temperature curves at L2 or L3 edges, it apparently that the fully developed temperature of q2 (red dots) is slightly higher than the that of q1 (black dots).