Evolution of the magnetic and polaronic order of $\rm{Pr_{1/2}Ca_{1/2}MnO_3}$ following an ultrashort light pulse

The dynamics of electrons, spins and phonons induced by optical femtosecond pulses has been simulated for the polaronic crystal $\rm{Pr_{1/2}Ca_{1/2}MnO_3}$. The model used for the simulation has been derived from first-principles calculations. The simulations reproduce the experimentally observed melting of charge/orbital order with increasing fluence. The loss of charge order in the high-fluence regime induces a transition to a ferromagnetic metal. At low fluence, the dynamics is deterministic and coherent phonons are created by the repopulation of electronic orbitals, which are strongly coupled to the phonon degrees of freedom. In contrast to the low-fluence regime, the magnetic transitions occurring at higher fluence can be attributed to a quasi-thermal transition of a cold-plasma-like state with hot electrons and cold phonons and spins. The findings can be rationalized in a more complete picture of the electronic structure that goes beyond the simple ionic picture of charge order.


I. INTRODUCTION
Unlike traditional spectroscopy, ultrafast pump-probe spectroscopy is a powerful tool to dynamically track the interactions in a strongly correlated system. Specifically, manganites are a suitable model system, because different types of correlations between electrons, spins and phonons have similar strength. A number of different electronic ground states can be realized just by changing temperature or doping. In these materials, high-resolution ultrafast pump-probe spectroscopy experiments have unraveled interesting physical effects, such as photo-induced phase transitions [1][2][3] and transient "hidden" phases [4].
In manganite perovskites, a growing number of ultrafast experiments provide access to the dynamics on different time scales from femto-to nanoseconds. The dynamics depends on the phase of the selected manganite as well as on the photon energy and intensity of the pump pulse. In GdSrMnO 3 close to half doping, a photo-induced transition from a charge-order phase to a ferromagnetic metallic phase within 200 fs has been observed [2]. An ultrafast metal-insulator transition has been induced in Pr 0.7 Ca 0.3 MnO 3 by selectively exciting phonon modes at 625 cm −1 [5]. On shorter timescales, coherent oscillations in the sample resistance in the insulator-to-metal dynamics with 30 THz are interpreted as orbital waves. [6] Several ultrafast optical-pump terahertz-probe studies in manganites focused on probing the nature of the quasi particles and their dynamics within a given phase [5][6][7].
involves thermalization of electronic subsystem and its energy redistribution with the lattice subsystem. The slower component on the 20-200 ps timescale, with a T cdependent lifetime, is attributed to the spin-lattice relaxation [3,10,11].
Recently, long-lived polaron-type optical excitations with lifetime of 1-2 ns are observed in the charge-ordered phase of Pr 0.7 Ca 0.3 MnO 3 [7]. Another long-lived intermediate state, with a mixture of ferromagnetic metallic and charge-ordered nanoscale domains, is observed in La 0.325 Pr 0.3 Ca 0.375 MnO 3 [12].
Depending on the phase and excitation intensity, coherent acoustic phonons as well as oscillating strain waves are observed on time scales of several 10 ps [13,14]; for manganites see [15].
The goal of this work is to augment previous experimental studies with a detailed description of the microscopic processes occurring during the first few picoseconds. For this purpose, we perform simulations, that are verified by comparing our findings with experimental observations. The work presented here is parameter free in the sense that the model parameters have been extracted from first-principles calculations [16].
We simulate the photo-excitation of Pr1 /2 Ca1 /2 MnO 3 by a femto-second light pulse and the subsequent relaxation of the magnetic and polaronic microstructure for the first few picoseconds. Ehrenfest dynamics [17] is adopted to propagate wave functions, spins and atoms. Peierls substitution has been employed to incorporate the external light pulse. A systematic study is performed to investigate the relaxation process by varying the intensity of the light pulse.
The paper is organized as follows: The methods of the paper are covered in section II: The tight-binding model used for the electron, spin and phonon degrees of freedom is described first, followed by the dynamical equations of motion, the treatment of the optical excitation and quantities used in the analysis. In section III, we revisit the arXiv:1909.03412v4 [cond-mat.str-el] 8 Mar 2020 ground state of Pr1 /2 Ca1 /2 MnO 3 , which experiences the optical excitation. Thereafter, we describe in section IV the results of our simulations and discuss the underlying mechanisms. Finally we summarize the findings in section V.

A. Tight-Binding model
To investigate the electronic, atomic and magnetic microstructure of manganites, we employ a tight-binding model [18][19][20]. The selection of energy terms and the parameter values have been extracted from the firstprinciples calculations [16].
The model describes the correlations between electron, spin and phonon degrees of freedom. The explicit electronic degrees of freedom describe the Mn-d-electrons with e g -character. The spin degrees of freedom describe the three majority-spin d-orbitals of each Mn ion with t 2g character. The phonon degrees of freedom are two Jahn-Teller active vibration modes and one breathing mode of each MnO 6 octahedron. In addition, we allow for a global expansion of the lattice.
The potential-energy functional of the system is where the E e , E S and E ph are the energies of the isolated electronic, spin and phonon subsystems, respectively. E e−ph is the electron-phonon coupling of Mne g electrons with the Jahn-Teller active modes as well as with the breathing mode of the MnO 6 octahedra. E e−S is Hund's coupling between the e g electrons and the spins of the Mn-t 2g electrons. Our model avoids the common infinite-Hund'scoupling limit, [21] and uses a more realistic description with explicit minority-spin electrons and a fully noncollinear treatment of the electron spin. Furthermore, we include the strong cooperativity of Jahn-Teller distortions and octahedral breathing modes by expressing them in terms of the explicit oxygen positions, which are shared each by two adjacent MnO 6 octahedra.
The Mn e g electrons are described by a Slater determinant formed by the one-particle wave functions |ψ n . The latter are expressed as superposition of local orbitals |χ σ,α,R with complex-valued coefficients ψ σ,α,R,n |ψ n = σ,α,R |χ σ,α,R ψ σ,α,R,n . ( The local orbital |χ σ,α,R is a spin-orbital with e g character at the Mn site R. It is a spin eigenstate with spin σ ∈ {↑, ↓} and spatial orbital character α ∈ {a, b} (where a=d x 2 −y 2 and b=d 3z 2 −r 2 ) [16]. The wavefunctions are Pauli spinor wavefunctions that account for the non-collinear nature of the magnetization. The three majority-spin t 2g electrons of a Mn site with site index R are described by a spin vector S R [16]. The two Jahn-Teller active phonon modes Q 2,R and Q 3,R for each MnO 6 octahedron as well as the breathing mode Q 1,R are expressed by the displacements of oxygen ions along the Mn-O-Mn bridge using Eq. 22 of [16].

B. Dynamics
The dynamics of the optical excitation and its relaxation processes are described by Ehrenfest dynamics [17,22]: That is, the electrons evolve under the timedependent Schrödinger equation, while the atoms are treated as classical particles and obey Newton's equations of motion.
i ∂ t ψ σ,α,R,n = ∂E pot ∂ψ * σ,α,R,n R j are the structural degrees of freedom of the oxygen ions and M O is their mass. m S is the absolute value of the magnetic moment of the t 2g -spins and the quasimagnetic field B R is Further details are given in appendices A and B. We also allowed for a strain dynamics, which does not have notable consequences on the results presented here. For the sake of completeness, it is described in appendix C.

C. Light Pulse
The light pulse is described by a spatially homogeneous, but time-dependent electromagnetic field (see appendix D) where A 0 is the amplitude of the vector potential and ω is the photon energy. c is the speed of light. The polarization of the electric field and of the vector potential is the unit vector e A . The temporal profile of the laser pulse is described by a Gaussian The intensity, which is proportional to |g(t)| 2 , has the full-width-at-half-maximum (FWHM) of 2c w √ ln 2. The light pulse is implemented with the Peierlssubstitution method [23,24]. Details are given in Appendix D.

D. Parameters of the simulation
Tab. I summarizes the relevant parameters used in this paper. We use a Cartesian coordinate system with the coordinate axes e x , e y and e z pointing along the Mn-Mn nearest-neighbor distances. The vectors e j are defined with length 1.
The Mn-Mn nearest-neighbor vectors are g x d M n−M n e x , g y d M n−M n e y , and g z d M n−M n e z with d M n−M n from Tab. I and scale factors g x , g y and g z . With N x , N y and N z , we denote the number of Mn sites in each of the three spatial directions in the supercell used in the calculation.
To describe perovskites, one usually refers to the lattice vectors of Pbnm unit cell, an non-standard setting of space group 62. The lattice vectors a, b and c of the Pbnm unit cell are The pulse length has been chosen consistent with the laser pulses used in ultrafast pump-probe experiments [7].

E. Diffraction patterns
In order to link our results with diffraction experiments, we inspect the intensities of diffraction for charges, orbitals and spins.
The intensity of diffraction of an observableX, such as number of electrons, orbital occupations, or spins, with density x( r) is [25] I X ( q) := I X,0 ( q) where the integration is performed over the illuminated sample volume V . I X,0 ( q) is the intensity of diffraction of a point scatterer x( r) = Xδ( r).
For a periodic lattice of atom-centered distributions (11) with the correlation function [26] The distributions x R ( r) are placed at the lattice sites R R + t, where R R is the position of an atom inside the first unit cell and t is the lattice translation vector of a specific unit cell. N is the number of atoms in the unit cell, Ω T is the unit-cell volume, Ω G = (2π) 3 /Ω T is the unit-cell volume of the reciprocal lattice and G are the corresponding general reciprocal-lattice vectors. N ∈V is the number of sites in the illuminated region. X R ( q) := d 3 r x R ( r)e −i q r is the form factor of x R . Note that the correlation function is meaningful only at the reciprocal lattice vectors G. Specifically, we address the following diffraction patterns: • The charge-order correlation function [26] probes the deviation of the electron density from its mean value, i.e X R = n R − n , where n R = α,σ ρ σ,α,R,σ,α,R is the number of e g -electrons on Mn-site R and n = 1−x with the doping x = 1 2 is the average number of electrons on Mn sites. The one-particle-reduced density matrix ρ σ,α,R,σ ,α ,R =: n f n ψ σ,α,R,n ψ * σ ,α ,R ,n (14) is given by the wave-function coefficients ψ σ,α,R,n and the occupations f n .
• The orbital-order correlation function probes difference between the orbital occupations X j,R = n x,R − n y,R , where are calculated for the set of orthonormal orbitals |θ j with j ∈ {x, y}, which are nearly axial in the x, respectively the y direction. These orbitals are defined in terms of the original basis states |d x 2 −y 2 and |d 3z 2 −r 2 as and We skip spin and site indices of the Wannier-like orbitals for scalar products, where they are identical.
• The spin-correlation function [26] probes the total spin X R = S R + s R of the Mn-sites, where S R is the spin of the t 2g electrons and s R is the spin of the e g electrons at Mn-site R.
In order to account for the blurring of the diffraction peaks due to fluctuations, we included for the timedependent spin correlation functions, in Figs 13, 19 and 23, the contribution from a (3 × 3 × 3) set of reciprocal lattice vectors of the super cell centered at the specified reciprocal lattice vector.

III. EQUILIBRIUM
Before investigating the optically-induced dynamics, let us remind of the salient features of the spin, charge, orbital, and lattice order of Pr1 /2 Ca1 /2 MnO 3 in equilibrium. The low-temperature phase serves as initial state for the excitation and it determines the time-evolution of the system. Pr x Ca 1−x MnO 3 has a perovskite lattice formed by a network of corner-sharing MnO 6 octahedra. Large cations such as Ca 2+ and Pr 3+ fill the voids in between the octahedra. The octahedral network distorts to optimize the Coulomb energy between the ions, which results in a characteristic pattern of the octahedral tilts. This tilt pattern fits into the orthorhombic Pbnm unit cell, which holds four octahedra.
The Mn-ions occur in the formal Mn 4+ and Mn 3+ oxidation states, with spin-aligned d-electrons on each Mn site. The additional electron of Mn 3+ produces a Jahn-Teller distortion, which is highly cooperative.
At half doping, Pr x Ca 1−x MnO 3 has a lowtemperature phase, which is described as charge and orbital ordered. [27][28][29] The low-temperature phase of Pr1 /2 Ca1 /2 MnO 3 is shown schematically in Fig. 1. It has a CE-type antiferromagnetic order exhibiting ferromagnetic zig-zag chains in the ab-plane, which proceed along the b-direction. These zig-zag chains are antiferromagnetically coupled among each other. Along the zig-zag chain, we can distinguish alternating central and corner sites. The central sites are described formally as Mn 3+ ions and exhibit Jahn-Teller distortion, while the corner sites are formally Mn 4+ ions with a negligible Jahn-Teller distortion. Left: CE-type magnetic order and orbital order in the ab-plane of the ground state of Pr1 /2 Ca1 /2 MnO3. Black and white symbols indicate up-and down-spins. The orbitalpolarized central orbitals are indicated by a d 3z 2 −r 2 orbital symbol in the corresponding direction. The corner sites have no orbital polarization and are indicated by circles. The gray square indicates the magnetic unit cell of the material in the ab-plane. Also shown are the lattice vectors a and b of the orthorhombic Pbnm unit cell. Right: Wannier-like orbitals along a trimer of the zig-zag chains of the CE-type magnetic structure. The sign of the orbital-lobes are indicated by black and white. The Wannier-like orbitals are orthogonal within and between trimers. The arrows connect orbitals with dipole-allowed transitions. Transitions between all other orbital pairs within and between trimers are dipole forbidden. The orbital |w1 in the majority-spin direction is filled. The optical excitation lifts electrons from |w1 to |w2 in the majority-spin direction.

Diffraction patterns
A transition at 250 K is attributed to the emergence of charge and/or orbital order [29] from a disordered polaron distribution at higher temperatures. The charge and orbital order has been explored experimentally by X-ray diffraction of the Mn-K-edge [30].
The diffraction patterns for the low-temperature phase are listed under the header CE-type in Tab. II. The diffraction peaks are quoted as relative coordinates (h, k, l) of the reciprocal lattice vectors in the setting of the orthorhombic (Pbnm) crystal structure obtained at room temperature.
Zimmermann et al. [30] exploited that the diffraction spots of charge and orbital order can be distinguished when scanning (h, k, l) = (0, k, 0) along the b-direction, the direction of the zig-zag chain. The dominant peaks of the Mn-K-edge with even integer k are due to the pseudo-cubic atomic lattice of Mn-sites. The charge order introduces additional diffraction peaks at odd-integer k, and the orbital order produces diffraction spots at halfinteger (but not integer) values of k, as shown in Fig. 2 and listed in Tab. II. At 175 K, Pr1 /2 Ca1 /2 MnO 3 undergoes a Néel transition. Neutron-diffraction studies [29] identify the magnetic lines characteristic for the CE-type spin order for Pr1 /2 Ca1 /2 MnO 3 .
Our simulations reproduce the diffraction patterns due to charge, orbital, and spin order for the CE-type ground state as shown in Fig. 2.
The diffraction spots of other magnetic orders listed in Tab. II will be used to characterize the evolution of the magnetic order following the light pulse: B-type refers to a pure ferromagnet, A-type refers to ferromagnetic planes that are stacked antiferromagnetically in cdirection, C-type refers to ferromagnetic Mn-lines running along the c-direction, which are antiferromagnetically aligned with respect to neighboring strands. In a G-type antiferromagnet, the Mn-sites are antiferromagnetic with respect to all their neighbors. The magnetic orders can also be described by their wave vector of the magnetization. Expressed by the pseudo-cubic lattice formed by the Mn-sites, the B-type order refers to k = π d M n−M n (0, 0, 0), A-type refers to k = π d M n−M n (0, 0, 1) C-type refers to k = π d M n−M n (1, 1, 0) and G-type refers to k = π d M n−M n (1, 1, 1).

Charge order, orbital order and Jahn-Teller distortions
The pattern of Jahn-Teller distortions in Fig. 1 has been attributed to checkerboard-like charge order in the ab-plane with Mn 3+ ions at the central sites of each seg- ment and Mn 4+ ions at the corner sites of the zig-zag chains. [28] More recently, this picture of charge order has been challenged: Rather than deducing the charge state from the pattern of Jahn-Teller distortions, experimental techniques such as core level spectroscopy, (XANES, ELNES) and neutron diffraction measurements of the magnetic moments provide a more direct access to the charge on the ions. These experiments rule out a fully ionic picture and indicate the absence of a considerable charge disproportionation [29,[31][32][33][34] This seeming contradiction between atomic structure and charge distribution can be reconciled by considering orbital polarization. [16] A Mn-ion has a complete orbital polarization, when only one of the two spatial e g orbitals is occupied, while the other is empty. Hereby, the shape and spin orientation of the occupied orbital is not relevant. When both spatial e g orbitals are equally occupied, the atom lacks orbital polarization.
To quantify the orbital polarization at site R, we determine the difference |f orb 1,R −f orb 2,R | of the orbital occupations f orb α,R , which are the eigenvalues of the spin-averaged local density matrix ρ orb R with matrix elements The orbital polarization P O is defined as where a and b denote the two Mn e g orbitals. Orbital polarization can be recognized indirectly via the resulting Jahn-Teller distortions. Mn-ions without orbital polarization do not exhibit a Jahn-Teller distortion, irrespective of the number of electrons in the e g shell. Hence, there is a direct link between Jahn-Teller distortions and orbital polarization. The connection to the charge order is, however, indirect. It is present only in the case of complete orbital polarization. This assumption is violated in Pr1 /2 Ca1 /2 MnO 3 [16]. The Jahn-Teller distortions of the corner sites are small, not because of their charge, but because of their lack of orbital polarization. The orbital polarization of the corner sites is small because Wannier-like orbitals from two segments of the zig-zag chain contribute equally to the two e g orbitals. [16] As shown earlier [16], the electronic structure of the low-temperature phase of Pr1 /2 Ca1 /2 MnO 3 can be rationalized using a specific set of Wannier-like states formed from the Mn e g orbitals. These states, shown in Fig. 1, are localized on specific segments of the zig-zag chains, which we denote, in the following, as trimers. The Wannier-like states are constructed as orthonormal eigenstates of a pseudo symmetry of a trimer, namely three orthogonal mirror planes through the central Mn-ion of a trimer. The functional form of the Wannier-like orbitals has been given in an earlier publication [16]. The requirements given above determine the Wannier-like orbitals up to a single parameter, which governs the charge disproportionation between central and corner sites. With a suitable choice of this parameter, the first Wannier-like state |w 1 describes the occupied states almost perfectly. This can be seen in Fig. 3, which shows that the occupied portion of the density of states can be attributed almost exclusively to |w 1 . To be specific, the one-particlereduced density matrix of the e g states is well described by ρ σ,α,R,σ ,β,R = m χ σ,α,R |w σm,1,m δ σ,σm δ σ ,σm × w σm,1,m |χ σ ,β,R where |w σ,α,m is a Wannier-like orbital with spin σ, spatial type j with j = {1, 2, 3, 4} according to Fig. 1 and the index m specifying a particular trimer. σ m denotes the majority-spin direction of the trimer with index m. In our model calculations, the charge on the central site is 3.75 e and that of the corner site is 3.25 e, which corresponds to a charge disproportion of q/e = 3.5 ± δ with δ = 0.25. This value lies within the range of values obtained from various experimental probes as discussed earlier [16].
The shape of the orbital |w 1 is responsible for the orbital order with strong orbital polarization on the central size and negligible orbital polarization on the corner sites. Thus, the electronic structure is consistent with both, the observed Jahn-Teller pattern and the more direct measurements of the charge state [29,[31][32][33][34] When the charge order is described in terms of integral charge states, they should be understood as oxidation states, which, per definition, attribute electrons as a whole to the more electronegative partner. [35] This is, however, a definition rather than a detailed description of an electron distribution. The notion of integral charge states Mn 3+ and Mn 4+ ions shall be understood in this context. The real charge distributions in manganites are more subtle. The photo-excitation in manganites occur through both d-to-d and p-to-d transitions in the spectral energy range ∼0.5-2.3 eV [16,[36][37][38][39]. While the d-to-d transitions occur between Mn-3d states, the p-to-d transitions occur between O-2p and Mn-d states [16,40,41]. The transitions observed experimentally in the ∼0.5-0.75 eV energy range are mainly dipole-allowed transitions from the occupied to the unoccupied e g states. While d-d transition on a single Mn-site are dipole forbidden, there are dipole-allowed transitions, which involve charge-transfer oscillations between different Mn-sites. [16] In this work, we focus entirely on transitions within the Mn e g shell. The charge-transfer transitions from O-p to Mn-d states dominate only at comparatively higher energies [16,[36][37][38][39].
A quantity used to describe the excitation is the photon-absorption density D p , which is the total number of photons absorbed per Mn-site and pulse. We calculate it as from the energy ∆E tot f −i added by the light-pulse to the system with N M n Mn ions (see also Fig. 31). The division by the photon energy ω and N M n provides the number D p of absorbed photons per Mn-ion and pulse.
Another quantity, often used to characterize experiments, is the pump-fluence F p . It is the energy transmitted to the sample per pulse and unit area. The pump fluence determines together with the pulse duration the intensity of the light-field.
The pump fluence F p is where c is speed of light and ε 0 is vacuum permeability.
With A 0 we denote the amplitude of the vector potential (see Eq. 7). The photon energy is ω. Due to the normalization of the pulse-shape function g(t), Eq. 8, the relation Eq. 24 is independent of the pulse duration. The spectral distribution of the photon-absorption density is shown in Fig. 4.
For the simulations discussed below, we have chosen the photon energy equal to the absorption maximum of ω =1.17 eV. The position of the absorption maximum appears to be rather independent of pulse length and light intensity. Fig . 5 shows the photon-absorption density D p as a function of the amplitude A 0 of the light field, as defined in Eq. 7. At the largest fluences, every third electron is excited, which explains the large changes of the magnetic and polaronic microstructure observed in those simulations.
The photon-absorption density in Fig. 5 grows approximately linearly with the amplitude of the light field. This behavior differs from the low fluence regime, where the photon absorption grows quadratically with the amplitude. The approximate linear behavior can be attributed to damping and decoherence. For a two-state system, the optical Bloch equations with damping have a steady state solution with an excited-state population of P e (A 0 ) = A 2 0 /(a+A 2 0 ), where a is constant determined, among others, by detuning and friction parameters. [42] P e has a point of inflection at a population of 25 %, which explains the approximate linear behavior on the amplitude for the fluences studied here. Additional features seen in Fig. 5 can be attributed in parts to strongly damped Rabi oscillations.

B. Regimes with distinct relaxation behavior
The relaxation following the optical absorption depends strongly on the pump-fluence F p . Based on the diffraction patterns, we identify four regimes with distinct relaxation behavior. The diffraction intensities of the characteristic diffraction spots are shown in Fig. 6. The boundaries of the different regimes depend little on the pulse duration.
The nature of these regimes are discussed in detail below. They can be characterized as follows: 1. In regime I, the spin, charge, and orbital-order is preserved. Coherent phonons with long lifetime are observed.
2. In regime II, the spin dynamics sets in, but the spin pattern relaxes back to the original state on a picosecond time scale. The charge order remains unaffected. Coherent phonons are present as in regime I 3. In regime III and beyond, the charge order is disrupted and the system is driven into a photoinduced ferromagnetic state. 4. In regime IV, the system enters a photo-induced anti-ferromagnetic state.
For the demonstration of the characteristic behavior in the four fluence regimes, we selected the fluence values in Tab. III.
The time-dependent distributions of excited electrons and holes are shown in Fig. 7. While the band structure for regimes I and II are qualitatively similar, in regime III the band gap at the Fermi level collapses as a ferromagnetic metallic state is formed. Also the band gap between minority and majority spins collapses. The band gap between majority-and minority-spin orbitals reoccurs in regime IV, where the system evolves into a new antiferromagnetic state. The antiferromagnetism is accompanied by a smaller band width of majority-and minorityspin bands so that the gap between them opens. Like in the ferromagnetic regime III, the system is metallic in regime IV.

C. Regime I
For the weak pump fluences of regime I, the magnetic, charge and orbital orders remain intact. The excitation can be described as formation of electrons and holes in an essentially rigid band structure. The electron-hole pairs are strongly coupled to breathing modes and Jahn-Teller active phonons at the Γ-point. As a consequence, two coherent phonons with a long lifetime are excited.
The excitation can be rationalized using the Wannierlike states introduced previously. They are shown in Fig. 1 for one segment of the zig-zag chain. Unless mentioned otherwise, the electric field of the light wave points along the b-direction of the Pbnm unit cell, that is along to the zig-zag chains of the ground-state magnetic structure.

During the light pulse
In the initial phase of the excitation, i.e. during the 100 fs light pulse, charge and orbital order drop to lower values. Furthermore, long-lived oscillations, discussed below, are initiated. On top of these effects, high-frequency oscillations of the electronic system are induced that, however, decay after few tenths of picoseconds. These oscillations are apparent in the charge-order and orbital-order diffraction peaks in Fig. 8 and Fig. 9. The only dipole-allowed transitions are between the bonding orbital |w 1 and the non-bonding orbital |w 2 as well as between latter, |w 2 , and the antibonding orbital |w 3 with the same spin direction shown in Fig. 1. There are no dipole-allowed transitions to |w 4 . Furthermore, there are no dipole-allowed transitions between Wannierlike orbitals from different segments of the zig-zag chains.
There is only one type of dipole-allowed transitions from the filled states. It lifts an electron from the majority-spin bonding orbital |w 1 to the non-bonding orbital |w 2 of the same segment and with the same spin.
The nature of the high-frequency charge oscillation on a segment of the zig-zag chain is rationalized via the Bloch waves of |w 1 and |w 2 character shown in Fig. 10, which are connected by the optical excitation. The excitation depends on the polarization of the light. Two representative Bloch waves of the initial state with |w 1 character and the corresponding final state with |w 2 character are shown schematically in Fig. 10 for each polarization in the ab-plane. The product of initial and final state wave functions is proportional to the first-order change of the charge density, which is in turn responsible for the dipole oscillation that couples to the light field.
For an electric field along the a-axis, i.e. perpendicular to the zig-zag chains, large dipole oscillations between the two corner sites are visible in Fig. 9. The charges of the two corner sites (red/green in Fig. 9) oscillate outof-phase and with the frequency of the light field. They describe the oscillating charge transfer between the corner sites. This is consistent with the Bloch waves (1) and (2) in Fig. 10 for e A || a. The product of initial, |w 1 derived, states and final, |w 2 derived, Bloch waves lead to charge contributions with alternating sign on the corner sites of a zig-zag chain. The central site (blue in Fig. 9) has a smaller oscillation with twice the frequency of the light-field. This oscillation is due to the rescaling of the |w 1 contribution required to maintain a normalized overall wave function, while |w 2 is mixed in.
When the electric field is polarized along the b-axis, i.e. parallel to the zig-zag chain, we observe in Fig. 9 only oscillations with small amplitude and with the doubled frequency of the light wave. The two corner sites oscillate in-phase with the doubled frequency. The charges on the central sites oscillate out of phase with the corner sites. This is consistent with the Bloch waves (3) and (4) in Fig. 10: The product of |w 1 -and |w 2 -derived waves has contributions only on the corner sites. However, there, the product of two orthogonal orbitals is formed, which does not contribute to the net charge. The net charge dipole coupling to the light field is not apparent in the bulk. It would show up at the surface of the material. Thus, only the second-order charge oscillations describing the charge transfer from the central site to the corner sites is visible in the bulk. The charge oscillations between the corner sites are initially absent and only kick in at later times.

Orbital order
Let us now turn from the high-frequency electronic excitations during the light pulse, to the changes that persist beyond the light pulse.
For the analysis, we expand the time-dependent oneparticle wave functions |ψ n (t) in Wannier-like orbitals The (j ∈ 1, 2, 3, 4) selects one of the four Wannier-states from a given trimer according to Fig. 1.
The instantaneous occupancy F tot j (t) of the j-th Wannier-like state |w j is and its spin polarization is with the number of Mn-sites N M n and Q j,m,σ,σ = n w j,σ,m |ψ n f n ψ n |w j,σ ,m . (28) As shown in Fig. 11, the dipole-allowed optical transitions from |w 1 to |w 2 dominate at the low pump fluences of regime I. Thus, the occupancy of |w 2 grows at the expense of |w 1 during the light pulse, while the occupancies of |w 3 and |w 4 remain small. After the lightpulse, the occupation of |w 2 remain almost constant, which is one sign of the preservation of the ordered state of the material. Often, e.g. [43], the excitation is attributed to an onsite d-to-d transition at the central Mn-site of a trimer segment, which formally is the Mn 3+ -ion. The picture, which emerges from our calculations, is more subtle. Rather than a dipole-forbidden excitation on central Mnion from |w 1 to |w 4 , the excitation is a |w 1 -to-|w 2 The description in terms of Wannier-like orbitals in Fig. 1 explains the strong optical absorption and it has consequences on the coherent modes described below.

Charge order
The |w 1 -to-|w 2 re-population rearranges electrons from the central to the corner sites, which reduces the charge disproportionation between central and corner sites. The reduced charge disproportionation reflects on the charge-order correlation C Q (0, 1, 0) shown in Fig. 8. The light pulse induces a sharp drop of the charge-order diffraction intensity C Q (1, 0, 0) from the initial value. Then, the intensity oscillates around this reduced intensity, with little sign of recovery during our simulation.
In our simulation, the orbital-order peak C O (0, 1 2 , 1) has a characteristic frequency of 10 THz, while the charge-order diffraction peak C Q (1, 0, 0) exhibits one frequency at 10 THz and a second one at 16 THz.

Coherent vibrations
The removal of electrons from the central site during the |w 1 -to-|w 2 transition reduces its Jahn-Teller distortion. The sudden change excites a Jahn-Teller mode with ν = 10 THz, which affects predominantly the central site. The displacements as function of time are shown in Fig. 12.
The charge transfer from central to the corner sites reduces the charge order and thus induces a planar breathing mode on the corner sites with a frequency of ν = 16 THz. Note, that an electron addition to a Mn site populates the Mn-O antibonds, which, in turn, expands the nearest neighbor distances. The expansion of the corner sites also affects the Jahn-Teller vibration on the central site.
On the time scale of a few picoseconds, the vibrations do not dissipate significantly. Furthermore, the phonon modes are fully coherent on the picosecond time scale of our calculation.
We attribute the lack of dissipation to the absence of heat conduction. In an experiment, a limited spot is illuminated and the energy can escape from the illuminated region by heat conduction. In our simulation, this process is prohibited, because the infinite material is illuminated homogeneously.
Another reason for an unexpectedly slow dissipation is the specific non-equilibrium state at hand. The coherent phonons are in contact with a phonon bath that is extremely cold. Thus, the collision probability of the coherent phonon with another one is extremely small.
The two highest frequencies measured in Nd 1/2 Ca 1/2 MnO 3 [44] at 10.2 THz and 14.1 THz agree very well with those in our simulations. This suggests a different assignment of the coherent vibrational modes: the mode previously assigned to a octahedral rotation mode at 10 THz, is in our simulation a Jahn-Teller oscillation on the central site of a segment. The mode assigned as Jahn-Teller mode at 16 THz, is in our simulation a symmetric breathing mode at the corner sites. Overlapping with the 10 THz vibration, we also find the antisymmetric expansion of the corner sites along the trimer axis. This latter vibration, however, is not coupled to the optical excitation. It should be noted that displacements of the oxygen ions in our model, denoted as the Jahn-Teller and breathing modes, also have a small implicit component from octahedral tilting.
Vibrations observed in the low-frequency range 2.4-7 THz [7,[43][44][45][46], have been attributed to A-type ion motion and rotational modes of MnO 6 octahedra. Our model does not describe these low-frequency modes because it does not contain explicit A-type ions. Pure octahedral rotations are not included, because our model does not describe oxygen vibrations perpendicular to the oxygen bridge.
Our description of the coherent phonons differs from that given earlier [44] as being due to an instantaneous melting of the charge and orbital order. The picture emerging from our simulations is that of a mechanistic rather than a thermal process.

Magnetic order
In the fluence regime I, the magnetic order is preserved: The spin angles deviate by about 10 degree from the ideal ground-state arrangement.

D. Regime II: Transient changes of the magnetic order
In the fluence regime II, the Jahn-Teller-active phonons respond similar to regime I. After an initial period of about 200 fs, however, also the spin order is perturbed. The spin system relaxes back into the ground state within few picoseconds, while the reduction of charge and orbital order persists much longer.

Magnetic order
As shown in Fig. 13, the spin correlations of the initial CE-type order drop to lower values at about 0.2 ps, while a prominent but short-lived signal of an A-type magnetic order shows up. This signal lives for about 0.2 ps before it dies out again. The spin-diffraction pattern during this period, which reflects the superposition of CE-type and A-type spin patterns, is shown in the left graph of Fig. 14.
As shown in Fig. 15, the angle of antiferromagnetic neighbors changes by up to 90 • . Besides some fluctuations, the ferromagnetic order within the zig-zag chain is preserved in the fluence regime II. What is affected most, is the spin correlation between neighboring chains in the ab plane. The onset time decreases with increasing fluence.
We attribute the response of the spin system to the optically-induced intersite spin transfer (OISTR) [47] caused by the coupling of the majority-spin |w 2 states with minority-spin states |w 1 and |w 4 of a neighboring zig-zag chain: the time-dependent populations of the rel-  Fig. 16. The |w 2 state is populated by the photo excitation. It is located at the corner sites and has lobes pointing towards the central atom of a neighboring zig-zag-chain. Thus, there is a spatial overlap of the |w 2 orbitals with the minority-spin |w 1 and |w 4 orbitals of a neighboring chain. The excitation into the |w 2 orbital is therefore accompanied by a spin transfer between neighboring zig-zag chains in the ab-plane. The delocalization of electrons among the antiferromagnetically coupled zig-zag chains changes the magnetization of the e g electrons, which in turn acts onto the classical spins describing the t 2g elec- In the period from 0.2 to 0.4 ps in Fig. 16, we observe a transfer of weight from |w 1 orbitals of one chain to the majority spin |w 4 of a neighboring chain. This is a consequence of the transient magnetic transition, which takes place during this period. The hybridization between these orbitals becomes possible because the spin orientation of neighboring chains deviate from 180 • .
In order to obtain a better understanding of these spin fluctuations, we investigated the spin structure as function of the ratio of excited electrons. For this purpose, we reduce the occupation for the 1 2 N M n occupied states from 1 to 1−δ and we increase the occupation of the first 1 2 N M n unoccupied states from 0 to δ. Then, we investigate the ground state as function of δ, where all degrees of freedom except strain are relaxed.
The original CE-type spin order is preserved up to δ = 11 %. For larger δ, a non-collinear but co-planar spin-wave structure emerges, where every two antiferromagnetic zig-zag chains in the ab-plane pair up and form a finite spin-angle with the spin axis of the next pair. The ferromagnetic spin order within the zig-zag chains and the strict antiferromagnetic coupling in c-direction remain intact. The angle of the spin axes grows with increasing δ until the angle approaches 90 • . At this point, δ ≈ 16 %, the system collapses into a ferromagnetic metallic state.
Hence, we interprete the spin fluctuation observed in the simulations as the onset of a Néel transition, which, in this case, is driven by photo-excited rather than by thermal |w 1 → |w 2 electron-hole pairs. In this context, the OISTR concept may be generalized to a temperatureinduced intersite spin transfer (TISTR).

Charge order
Charge and orbital order, Fig. 17, as well as coherent phonons, Fig. 18, show the same behavior as in regime I, albeit with larger amplitude. The oscillations of the phonon displacements and the correlation functions for charge and orbital order are long-lived. Unlike regime I, a slow decay of the oscillations is noticeable in regime II. The coherent phonons, charge order, and orbital order are only little affected by the transient change of the magnetic structure. The reason is that the charge order and orbital order remain intact during the transient change of the magnetic correlations.
Coherent phonons, charge order, and orbital order, are strongly coupled via electron-phonon coupling. They are due to the same physical mechanism, namely the |w 1 to |w 2 excitation. Thus, they are expected to exhibit similar decay properties.  In regime III, the system undergoes a photoinduced phase transition, which converts the CEantiferromagnetic order into a ferromagnetic metallic state without charge and orbital order. Initially, the antiferromagnetic correlation between the zig-zag chains is perturbed rather similar to regime II, leading to an A-type magnetization as seen in Fig. 19.
Unlike regime II, however, the A-type diffraction pattern persists for several picoseconds. During this time, the diffraction pattern of the ferromagnet builds up until it replaces the A-type diffraction pattern altogether.
The ferromagnetic state obtained is not fully established in the simulation: The non-relativistic Schrödinger equation employed in the simulations conserves the total spin. As a result, the system evolves into a state that is better characterized as a spin wave or a lattice of ferromagnetic domains.
Nevertheless, the diffraction pattern obtained is very similar to the ferromagnetic structure. For each diffraction spot of the ferromagnetic structure, we do not obtain a single spot, but a set of two "twin-peaks". The two peaks are located at the supercell reciprocal-space vectors adjacent to those of the ideal ferromagnet as seen in Fig. 20. The displacement of the twin peaks from the diffraction spot of a true ferromagnet is governed by the size of our supercell, which limits the wave length of the spin wave, respectively the domain size.
The ferromagnetic spin correlation function in Fig. 19 exhibits a finite signal by considering the contribution from the immediate neighborhood of the specified reciprocal lattice vector. The signal at the center of the spot is zero.
We envisage that a larger supercell leads to larger domains and thus to twins-peaks that are even closer together, making them indistinguishable by experiment. As shown in Fig. 21, the charge-order correlation C Q and the orbital-order correlation C O are completely wiped out after about 0.2 ps. The loss of orbital order makes the system metallic as seen in Fig. 7. The loss in charge and orbital order is also reflected in attenuation of the phonon displacements shown in Fig. 22.
We attribute the ferromagnetic order to a mechanism in the spirit of the double-exchange picture. [48-  50] The origin of the band gap in the ground state of Pr1 /2 Ca1 /2 MnO 3 is the formation of Zener polarons. These Zener polarons are also the origin of charge and orbital order. Due to the re-population of electrons across the band gap, the stabilization due to Zener polarons is lost and another, competing mechanism can take over. A ferromagnetic alignment of the spins lowers the kinetic energy of the electrons, which can now spread over a large area: An electron with a given spin is effectively excluded from orbitals of Mn-ions with opposite spin. A Mn-ion with opposite spin thus leads to an energy cost. Thus, it is favorable, when all spins align ferromagnetically. In other words, a configuration of ferromagnetic spins produces a larger effective band width of the majority-spin configuration. The larger band width stabilizes electrons which populate the lower half of the majority-spin band formed.

Experiment
A recent ultrafast pump-probe experiment [46] carried out on Pr x Ca 1−x MnO 3 at 100 K with different pump fluences showed that the characteristic charge-and orbitalorder reflection peaks of the CE-type ground state disappear for larger fluences F p >2.5 mJ/cm 2 .
The experimental study with the same material class by Li et al. [51] revealed a photo-induced ferromagnetic state within about 120 fs above the threshold fluence F p = 2.4 mJ/cm 2 . A rise in magnetization has also been measured by Zhou et al. [52].
The measured threshold fluence of F p = 2.5 mJ/cm 2 translates via Eq. 24 into an amplitude A 0 = 1.26 /(ea 0 ). In our simulations, the loss of charge and orbital order sets in with A 0 = 0.58 /(ea 0 ) (F p = 0.52 mJ/cm 2 ) as seen in Fig. 6 from the appearance of the A-type magnetic order, which finally converts into the B-type (ferromagnetic) order. Our simulations and experiments produce the transition in the same fluence range. The remaining difference of a factor five in the threshold fluence may be attributed, for example, to the different photon energies.
It is worth mentioning that the ferromagnetic states observed in our study for region III is expected to persist on longer timescale hinting towards its possible long lifetime. Similar long-lived states are recently observed in Pr x Ca 1−x MnO 3 series within the charge-ordered region of the phase diagram [7].

F. Regime IV: Non-collinear antiferromagnet
In regime IV with the highest fluence, the system evolves first into a G-type antiferromagnet as shown in Fig. 23, rather than forming an A-type antiferromagnet as in regime III. After about 1.5 ps, the diffraction spots of the G-type structure fall off again in favor of a more complicated structure with non-collinear magnetic order. We note that this regime may be difficult to access experimentally due to the limited stability of the material. shown in red, respectively green, are characteristic for the CE-type ground state. The ferromagnetic (B-type) peak CS(0, 0, 0) is shown in orange. The peak with hkl = (1, 0, 1) is characteristic for G-type magnetic structure. Over time, a new, non-collinear magnetic structure (see Fig. 25) emerges evidenced by the occurrence of the peak hkl = (0.5, 0.5, 1).

(color online)
The diffraction patterns calculated for characteristic times along the trajectory are shown in Fig. 24. The spin structure of the final state, which emerges at approximately 1.5 ps after the light-pulse has been extracted on the basis of the real space spin-correlation function. It is shown in its idealized form in Fig. 25. All spins are perpendicular to its neighbors in the ab-plane and antiferromagnetic in the c-direction. This model produces the spin diffraction pattern of the final spin configuration of regime IV. To our knowledge this configuration has not been investigated before in the context of manganites.
The charge and orbital order is destroyed almost immediately, that is during the light pulse as shown in Fig. 26. This destruction of the orbital order excites phonons, that, however, dissipate on a picosecond time scale as seen in Fig. 27. The amplitude of the phonon vibrations is considerably larger than that in regime III.
We attribute the transition with increasing fluence from a ferromagnet in regime III to the non-collinear antiferromagnetic structure in regime IV to a mechanism analogous to that described earlier. [53] The doubleexchange mechanism favors ferromagnetism through the increase in band width only, when a majority of electrons populate the lower half of the majority-spin e g states. Thus, increasing the fluence beyond a certain point switches off the double-exchange mechanism again, so that a antiferromagnetic structure can develop. The re-duction of the band width of majority-spin and minorityspin electrons opens a band gap between them, which is seen in Fig. 7. Compared to the minimal model [53], our system evolves into a more complicated non-collinear antiferromagnetic structure.

G. Thermalization
One of the questions of interest is how thermal equilibrium is established from the excited state.
Therefore, we investigated the evolution of the temperatures of the subsystems.
Below, we consider the quasi temperatures described in appendix D 1. As shown in Fig. 28, the light pulse immediately raises the temperature of the electronic system to high temperatures, i.e. several thousands of Kelvin, while the temperatures of the phonon and spin systems remain low in comparison. Thus, a state far from equilibrium is formed. The state is analogous to that of a non-thermal (cold) plasma, where the electrons reach 10 4 K, while the ions remain near room temperature.
While the temperature of the phonon system remains cold, the coherent phonons of regime I and II are strongly coupled to the electronic subsystem and reach comparable temperatures. When we attribute the complete thermal energy of the ions considered to the two phonon modes, the resulting temperature of the two modes is comparable to the electronic temperature.

Non-equilibrium distribution of the electrons
In particular, the channel of the electrons is of interest because the energy from this channel is most easily put into practical use, such as in a solar cell.
Our simulations may shed light onto the workings of the Boltzmann equation. For this purpose, we inspect the emergence of a distribution, i.e. the occupations, as function of energy and we compare the distribution obtained in our calculation with the Fermi distribution. The approach to a Fermi distribution is one of the common assumptions made for the Boltzmann equation.
In order to explore the approach to the Fermi distribution, we choose a representation of atanh 1 − 2f j versus energy BO j . In this representation, a Fermi distribution maps onto a straight line with slope 1/k B T ψ and zero µ ψ .
The Born-Oppenheimer energies BO j , their occupationsf j and the Born-Oppenheimer wave functions |φ BO j are defined in appendix D 1 b. The occupations for different time slices and for the four regimes are shown in Fig. 29.
Initially, that is 0.2 ps after the maximum of the light pulse, the occupations do not lie on a continuous function of the one-particle energies but scatter wildly. This is expected because the occupations of the electron states are dominated by the ground-state occupations, the photon energy and the dipole matrix elements.
Already after an initial period of about 0.5 ps after the pulse maximum, the occupations form a continuous function of the energy. This indicates that the thermalization between electrons of the same energy is very efficient.
This result is specific to the choice of one-particle orbitals and energies, namely the Born-Oppenheimer states.
However, on the time scale of our simulation, the system does not approach a Fermi distribution. Rather, occupations of electrons further away from the chemical potential deviate further from integral occupations than a Fermi distribution: The electrons further away from the Fermi level seem to be "hotter" than those close to the Fermi level. Surprisingly, the tails of the distribution become even flatter with time, that is they seem to deviate even further from a Fermi distribution.
We attribute this behavior to the strong coupling between different subsystems: Due to the dynamics of the spin and phonon systems, the electrons experience a time-dependent Hamiltonian, that constantly drives the electrons out of their equilibrium distribution. For the one-particle basis |φ BO j (t) used here, the approach to a Fermi distribution is not a requirement.

Cold-plasma model
On the picosecond time scale after an excitation, we find a large disparity between the high temperature of the electrons on the one hand and the low temperatures of spins and phonons on the other hand. This suggests that the optically accessed states are the result of thermodynamic equilibrium of the electron system alone. A quasi-equilibrium state such as this has been assumed earlier [14].
In order to test this conjecture, we investigated the phase diagram by increasing the temperature of the electrons, while spins and phonons are kept at zero temperature. That is, spins and phonons are optimized for each electron temperature. The lattice constants are kept equal to the values before excitation, because they are usually too slow on the short time scale under consideration.
As shown in Fig. 30, we find four different temperature ranges, which, however, do not directly correspond to the ranges of different excitation behavior.
• For T < T 1 with T 1 ≈ 2500 K, we obtain the charge-ordered phase with CE-type magnetic order. Charge-order correlations and the corresponding Jahn-Teller distortions on the central site vanish at T 1 with an approximate square-root behavior ∼ (T 1 − T ) α with 0 < α < 1.
• T 1 < T < T 2 with T 2 = 4200 K: At T = T 1 the charge and orbital order is completely lost, and the system transforms abruptly from the CE-type antiferromagnetic structure into a ferromagnet. The system is a pure ferromagnet only at T = T 1 . For T > T 1 the spin angle φ c between adjacent abplanes increases with an approximate square-rootlike behavior towards increasing temperatures, i.e. φ c ∼ (T − T 1 ) β with 0 < β < 1. The spin orientation alternates between two values from plane to plane.
• T 2 < T < T 3 with T 3 = 5000 K: At T = T 2 , the spins in the ab-planes become non-collinear. The angle φ ab between adjacent spins in the ab-plane increases approximately linearly from 0 • to 180 • as the temperature is raised from T = T 2 to T = T 3 .
• T > T 3 : At T = T 3 , both spin angles, φ ab and φ c , are 180 • , which corresponds to the G-type magnetic order. This is the favorable high-temperature phase for the temperature range explored.
It is important to note that the phase diagram described here has little to do with the equilibrium phase diagram of the material. The phases described above are extreme non-equilibrium states, because spins and phonons are at T = 0.
We can identify the excitation regimes I and II with the temperature range T < T 1 . The ferromagnetic state obtained in regime III can be attributed to the range T 1 < T < T 2 . A non-zero spin angle φ c between the ferromagnetic planes has not been apparent in our timedependent simulation. We expect this to be a fluctuating quantity that is averaged out.
In regime IV, we find configurations which are non-collinear in the ab-plane. Interestingly, the non-collinear state with φ ab = 90 • shown in Fig. 25 is a typical state obtained for a range of fluences, while in Fig. 30 it is just one point in a region with continuously changing angles φ ab . At even higher fluences, also the G-type structure is encountered.

V. SUMMARY
The optical excitation of half-doped Pr 0.5 Ca 0.5 MnO 3 has been simulated to study the physical interplay between electronic, spin and lattice degrees of freedom in response to a femtosecond light pulse. The simulations use Ehrenfest dynamics, in which electrons and spins follow the time-dependent Schrödinger equation while the nuclei proceed on a classical trajectory.
Femtosecond excitations with various intensities and pulse lengths are studied. The pulse acts on the chargeordered, low-temperature phase with CE-type antiferromagnetism.
Four different intensity regimes with qualitative different behavior could be identified.
• In regime I, the electron-band structure remains essentially rigid. The electron-hole distribution excitation transfers weight from the central Mn ions of the zig-zag chain to the corner sites. The dipole oscillations shuffle charge between two adjacent corner sites. Two coherent phonons with long life time are excited as result of the electron-phonon coupling.
• In regime II, the spins react and rearrange into short-lived A-type antiferromagnetic struc-ture. The ground-state CE-type antiferromagnetic structure is recovered within a picosecond. The coherent phonons, present also in regime I, survive this transition.
• In regime III, charge and orbital orders are destroyed within few hundred femtoseconds and a ferromagnet is formed. In contrast to regime I and II, the coherent phonons are damped out rapidly. Due to spin conservation, the ferromagnet is not directly accessible. Rather, an A-type antiferromagnet is formed, which evolves over several picoseconds into a ferromagnet having domains compatible with the size of our simulation cell.
• In regime IV, charge and orbital orders are immediately destroyed as in regime III, but now a G-type antiferromagnet is formed rather than an A-type antiferromagnet. Over time, the system evolves into a new non-collinear spin structure with neighboring spins having 90 • angles in the ab-plane.
The transient magnetic state observed in regime II may shed light onto the thermal Néel transition of Pr1 /2 Ca1 /2 MnO 3 at 175 K. In regime II, the system maintains the orbital and charge order, but it modifies the spin correlations of neighboring zig-zag chains of the CE-type spin structure in the ab-plane. Analogously, the Néel transition may be due to a melting of the antiferromagnetic correlations between the zig-zag chains, while maintaining the ferromagnetic order within the chains. When the ferromagnetic order within the chains melts at higher temperature, the integrity of the chains with their orbital and charge order is destroyed as in regime III and IV.
The long lifetime of the magnetic orders in regime III and IV may qualify for the concept of "hidden phases". Hidden phases [4] are states with unique order which can not be accessed thermodynamically. It must be noted however, that the time scales covered in our simulations are short compared to those studied experimentally.
In order to make contact with thermodynamics, we estimated the temperatures of the individual subsystems, namely electrons, spins and phonons. The temperature of the electronic subsystem raises quickly to several thousand Kelvin, while phonon and spin degrees of freedom remain relatively "cold". An exception are the coherent phonon modes, which initially reach the temperature of the electrons before dissipating their energy into other degrees of freedom.
Following this concept of hot electrons and cold phonons and spins, we have been able to identify the phases accessed by optical excitation with those obtained by raising only the electron temperature.

ACKNOWLEDGMENTS
Financial support from the Deutsche Forschungsgemeinschaft (SFB 1073) through Projects B02, B03 and C03 is gratefully acknowledged. We are grateful to Michael Ten Brink, Salvatore Manmana and Stefan Kehrein for fruitful discussions.

Appendix A: Numerical integration of time-dependent Schrödinger equation
To solve the time-dependent Schrödinger equation for wave functions and spinors, we use the second-order differencing scheme proposed by A. Askar and Cakmak. [54,55].
Given the wave function |ψ(0) at time t = 0 and the time-dependent HamiltonianĤ(t), the wave function |ψ(t) can be obtained as |ψ(t) =Û (t, 0)|ψ(0) using the propagator T D is Dyson's time-ordering operator [56], which rearranges all operators in a product into ascending time order from right to left. With the time step ∆, subsequent wave functions of a time sequence are related by The error is reduced by splitting off the dynamical phase using the corresponding energy expectation value E n (t) := ψ n (t)|Ĥ(t)|ψ n (t) e i∆ En(0) |ψ n (∆) − e − i∆ En(0) |ψ n (−∆) This leads to the following iterative scheme These equations of motion are time-inversion symmetric per construction. However, the equations of motion produce besides the correct solution also a spurious partial solution which changes sign in each iteration. This implies that, over time, the wave function will pick up a contribution from the spurious solution. In order to purify the solution, we interrupt the simulation at regular time intervals and perform a correction step. In the correction step, we filter out the spurious partial solution.
and perform a Gram-Schmidt orthonormalization on the one-particle wave functions for the two subsequent time steps used in the next iteration. We use a time step of ∆ ≈ 10 −3 fs. Correction steps are performed every 20 time steps.
The energy conservation is shown in Fig. 31 for different light amplitudes A 0 . The dynamics of the spins S R describing the t 2g electrons requires special attention. While the spin dynamics is intrinsically of quantum nature, we want to keep all three t 2g electrons of a given Mn ions strictly collinear. For this purpose, we map the spin vector S i onto a normalized, complex-valued, two-component spinor a ↑,R a ↓,R such that The magnetic moment m S of the t 2g shell is anti-parallel to its spin direction, namely m S = −m S The spinors (a ↑,R , a ↓,R ) evolve under the time-dependent Schrödinger equation The summation index R ∈ N N R runs over nearestneighbor sites of site R. The first term in Eq. B4 describes the antiferromagnetic coupling with neighboring spins. The second term in Eq. B4 describes Hund's coupling between t 2g and e g electron on the same site. J AF is the antiferromagnetic spin-coupling parameter and J Hund is the Hund's coupling parameter. The dynamical equation Eq. B3 is equivalent to where × denotes the vector product.

Appendix C: Strain dynamics
In the present study, the scale factors g x , g y and g z are dynamical variables, which describe long-wave length acoustic modes that are responsible for the strain effects in manganites [45,46]. We enforce g x = g y .
In order to describe the sound wave observed in experiment, we introduce a classical kinetic energy In thin-film experiments, the wave vector of a sound wave perpendicular to the film is quantized, which results in a standing wave with a characteristic frequency. The sound wave modulates the optical density of the material, which can be detected by the optical absorption measurements. [13][14][15] We adjusted the fictitious mass M g to simulate this effect. In our model, a sound wave is excited in a material without e g electrons with q = 0 and with polarization along the c-direction. M g is chosen so that our model material oscillates with the same frequency as the film in experiment [7], namely ≈ 25 GHz.
Appendix D: Peierls substitution In this section, we give a brief derivation of the Peierlssubstitution method [23,24].
The electric field E = −∂ t A of the light pulse is expressed by a vector potential e A is the polarization direction of the vector potential and g(t) is an envelope function, which is normalized so that dt |g(t)| 2 = 1 .
The electrons experience a Hamiltonian of the form where q is the electron charge, m e its mass, and V is the lattice potential. The Hamilton matrix elements are evaluated in a basisset of local orbitals centered at positions R α , that have the vector potential explicitly built in. From a regular basisset of local orbitals |χ α , field-dependent basis functions are constructed [23,24]. The integral of the vector potential is path dependent: we choose a straight line from the central atom R α to the position r.
Substituting the above ansatz Eq. D4, we obtain From Eq. D3 and Eq. D5, we obtain is the Peierls phase. Furthermore, we define the small quantity F α,β ( r, t), which appears in the above Eq. D6, as F α,β ( r, t) is a magnetic flux through triangle with corners at R α , R β and r.
Furthermore, the vector potential is approximated by a constant. This is equivalent to the long-wavelength limit. It also excludes dipole-forbidden, but quadrupoleallowed transitions. The latter are not considered relevant in comparison with the strong charge-transfer transitions in the present work.
With these approximations, by exploiting the orthonormality of our basisset, and after ignoring off-site terms of the dipole matrix elements, we obtain The first term on the right-hand side describes chargetransfer transitions, while the second term describes dipole-allowed onsite transitions. The latter vanish in our model and are included here only for the sake of completeness. The Peierls phase only affects off-site matrix elements. In our case, these are the hopping matrix elements. Thus, the only change required to incorporate the excitation is to multiply the hopping matrix elements with the timedependent Peierls phase.

a. Phonon temperature
The temperature T ph of the phonon degrees of freedom has been evaluated from the kinetic energy of Jahn-Teller-active and breathing phonon modes of the oxygen ions. We use the relation where index i runs over all N O oxygen ions in the unit cell and M O is the mass of the oxygen ion. There is only one degree of freedom per oxygen atom in our simulation, because only three phonon modes per formula unit are considered. These are the modes with strong electron-phonon coupling, which receive the energy directly from the excited electrons and holes. Only later, these "hot" modes dissipate their energy into the other phonon modes and the spin system.

b. Electron temperature
The temperature of the electronic degrees of freedom are obtained from the occupations of the Born-Oppenheimer wave functions. For that purpose, we extract the one-particle wave functions |ψ j (t) and the instantaneous one-particle Hamiltonianĥ BO (t) acting on the electrons. The Hamiltonian is the Born-Oppenheimer Hamiltonian for the instantaneous spin distribution and atomic positions.
Let |φ BO j (t) be the eigenstates and BO j (t) the eigenvalues of the one-particle Hamiltonianĥ BO (t). The occupationsf j of the Born-Oppenheimer states |φ BO j (t) are obtained from their projections onto the occupied wave functions |ψ j (t) as Electron temperature T ψ and electron chemical potential µ are determined such that energy and particle number coincide with those of a thermal distribution, i.e.
with B R defined in Eq. B4 and the absolute value m S of the t 2g magnetic moment. The projections of the instantaneous Pauli spinors a R (t) onto the eigenvectors a BO j,R (t) of h BO,S R (t) yield occupationsf j,R (t) = σ∈{↑,↓} a * σ,R (t)a BO σ,j,R (t) 2 (D18) for ground state with j = 0 and excited state with j = 1.
The comparison with the internal energy for noninteracting spins in an external magnetic field provides a relation which is resolved for the instantaneous temperature T S (t) of the spin system. A more detailed derivation of the expressions summarized here is provided in appendix E.

Appendix E: Temperature of the spin subsystem
The temperature of the spin system is extracted analogously to that of the electrons. We consider a system of uncoupled spins in a magnetic-field distribution B R defined by the local Born-Oppenheimer Hamiltonian for the spin system according to Eq. D17. The free energy of this system is where the energy of a spin distribution σ is σ R ∈ {0, 1} characterizes the ground and excited state of the local spin S R of the t 2g electrons at site R. The absolute value of the magnetic moment related to the spin S R of the three t 2g electrons at a given site is m S = |γ 3 2 |. This yields for the free energy The instantaneous temperature of the spin system is extracted by comparing the energy obtained from the instantaneous spin distribution with the energy E T = β∂ β F T of the thermal ensemble, where β = 1/(k B T ).
With the eigenstates a BO j,R and the eigenvalues˜ j,R of the Born-Oppenheimer Hamiltonian h BO,S R , the instantaneous energy Eq. E4 is The frequencies of the coherent modes, present in regimes I and II, have been extracted by non-linear curve fitting of a superposition of a constant and two cosine functions with amplitude, frequency and phase shift as variable parameters. The quality of the fit is shown in Fig. 32. The fit gives a frequency of 9.7 THz and 15.7 THz.
The vibration of 15.7 THz is dominant in the ∆d 4+ y and can be attributed to the planar breathing mode on the corner sites of the CE-type magnetic structure, which couples to the charge transfer from the central to the corner sites.
The lower frequency with ν = 9.7 THz is a Jahn-Teller mode on the central site of a trimer. In order to confirm that the oscillations of the chargeand orbital-order correlation functions are a direct consequence of the coherent phonons, the correlation functions have been fitted with same two frequencies. The perfect agreement shown in Fig. 32 supports our conjecture.