Quantum turnstile operation of single-molecule magnets

The time-dependent transport through single-molecule magnets coupled to magnetic or non-magnetic electrodes is studied in the framework of the generalized master equation method. We investigate the transient regime induced by the periodic switching of the source and drain contacts. If the electrodes have opposite magnetizations the quantum turnstile operation allows the stepwise writing of intermediate excited states. In turn, the transient currents provide a way to read these states. Within our approach we take into account both the uniaxial and transverse anisotropy. The latter may induce additional quantum tunneling processes which affect the efficiency of the proposed read-and-write scheme. An equally weighted mixture of molecular spin states can be prepared if one of the electrodes is ferromagnetic.


Introduction
The single-molecule magnets (SMMs) are foreseen as building blocks of organic spintronic devices [1,2]. Such systems generally behave as magnetic cores with a large localized spin and display slow relaxation of magnetization at low temperature mostly due to presence of the anisotropy-induced magnetic barrier [3]. Similar to quantum dot physics, two-terminal steady-state transport measurements performed on SMMs revealed charging effects such as Coulomb blockade, sequential tunneling or negative differential conductance [4,5]. In the spin sector, Kondo related features were observed experimentally [6,7,8,9] and investigated theoretically [10,11]. The exchange interaction between the local molecular moment and the delocalized spins tunneling through the molecular orbital might be exploited to control the quantum state of the local moment, i.e. to 'write' and 'read' its quantum state [12].
On the experimental side, various techniques [13] are currently used to attach the orbitals (ligands) surrounding the molecular magnetic core to the source and drain probes. Unlike standard transport setups used in quantum dot devices, molecular electronics requires more careful handling of the contact regions. The difficult task of isolating a single molecule between source and drain electrodes is nowadays realized by using more advanced methods as electromigration, mechanically controlled break junctions [14] or spin polarized STM [15]. Recently several groups pushed even further these techniques and reported controlled time-dependent transport measurements for such SMMs when the contacts were switched on and off by varying the substrate-STM tip spacing [16,17] or by bending break junctions [18], and transient currents arising when a molecular tail couples to an STM tip have been recorded [19]. SMMs have also been integrated into carbon nanotube transistors to serve as detectors for the nanomechanical motion due to the strong spin-phonon coupling [20,21]. In two cornerstone experiments Vincent et al. [22] and Thiele et al. [23] detected the nuclear spin of a single Tb 3+ ion embedded in a SMM together with the Rabi oscillations.
These promising experiments motivated us to investigate the transient transport properties of SMMs, with special emphasis on the regime when the couplings to the source/drain electrodes are switched on and off periodically. The transport regime we are interested in is similar to the so called turnstile pumping setup which represents a long-standing [24] asset of pumping or pump-and-probe experiments with quantum dots [25]. Along this periodic pumping, the source and drain tunneling barriers open and close consecutively, such that a single electron is transmitted across the sample; details on the turnstile operation will be presented in Section 3. To our best knowledge the transient and the turnstile regimes have not been theoretically investigated so far in the context of transport across SMMs.
On the theoretical side the magnetic interactions in SMMs are described by effective giant spin Hamiltonians [3], mostly because of the large value of the localized magnetic moment. Using this description, in Ref. [26] the authors investigated the role of relaxation on inelastic charge and spin transport across a SMM weakly coupled to metallic gates. In Ref. [27] transport across a SMM coupled to two ferromagnetic leads was investigated and it was found that the spin current across the SMM can reverse the localized spin if the leads have opposite spin polarizations. Memristive [28] and thermoelectric [29] properties of SMMs were also investigated. All these studies convey a similar message: The current can induce magnetic switching of the localized magnetic moment if the applied bias voltage exceeds the gap between the ground and excited states.
In the absence of transverse anisotropy the effect of quantum tunneling of magnetization (QTM) is negligible and the full magnetic switching requires the transient occupation of all excited (intermediate) states with magnetic quantum numbers in the range [−S, S]. When present, the QTM might leave its fingerprint on the transport properties at resonant values of an applied magnetic field [30].
A complementary approach to transport properties of SMMs relies on density functional theory (DFT) [31,32]. In the DFT approach the molecular structure and the contact regions are carefully taken into account, while the many-body correlations within the SMM are accounted for within various approximations. A detailed ab initio Hubbard many-body model for molecular magnets has been recently implemented [33] and allows the calculation of magnetic interactions.
In the present work we investigate transient transport and turnstile pumping across a SMM embedded between magnetic and non-magnetic electrodes. As we are interested in the time dependent evolution of the currents and the accumulation of the charge and spin on the SMM we rely our investigation on the generalized master equation (GME) technique [34,35]. Let us stress that, to capture the turnstile regime, one has to go beyond the steady-state rate-equation approach. We find that by setting the SMM in the quantum turnstile (QT) configuration with ferromagnetic leads one can address two new issues which are relevant for the use of molecular states as magnetic qubits: (i) the one-by-one all-electrical writing and reading of excited molecular states with spin S − 1, .., −S + 1 (S being the molecular spin of the initial ground state) and (ii) the controlled preparation of statistical mixtures of such intermediate states. So far, the stepwise magnetic switching protocols for excited molecular states that we propose here have not been investigated. In fact, previous studies (see e.g. Refs. [27,26]) were focused only on investigating the full magnetic switching. The second issue was partially addressed by Tejeda et al. [36] some time ago. Their proposal concerns the preparation of equal weight superposition of states (e.g. |ψ = 1/ √ 2(|S + |S − 1 )) using at least two molecular clusters embedded in micro-SQUIDs.
The rest of the paper is organized as follows: In Section 2 we present the theoretical framework by introducing the model Hamiltonian and giving a summary of the GME method. Section 3 presents the main results of our work while in Section 4 we give the conclusions.

Model Hamiltonian
The setup that we consider here consists of a SMM coupled to two external electrodes (see the sketch in Fig. 1). We investigate time dependent transport in the sequential tunneling regime, in which the electrons tunnel one by one from the left (source) electrode to the unoccupied molecular orbitals of SMM and then escape to the right (drain) electrode. In the present work we neglect all other possible transport mechanisms, such as the cotunneling processes which are responsible for the Kondo effect [10].
The model Hamiltonian contains several terms describing the SMM itself (H M ), the left (H L ) and right (H R ) electrodes and the time-dependent tunneling part H T : In general, SMMs are characterized by a large spin S > 1/2. They do also present transverse anisotropy and an easy axis of magnetization [3]. Assuming that the relevant contribution to the transport comes only from the lowest unoccupied molecular orbital (LUMO) the SMM can be modelled by an effective, minimal Hamiltonian [26,27] of the form: In its simplest form the LUMO orbital consists of a single spinful interacting level with energy ε, on-site Coulomb energy U and occupationn =n ↑ +n ↓ , and is coupled to the localized spinŜ through an exchange interaction with a coupling strength J. The fourth and the fifth terms in Eq. (2) describe the easy-axis and transverse anisotropy with the corresponding constants denoted by D and E. For certain molecules [38] E ≪ D but this is not always true, as E can grow up to D/3 in other situations [3].
In the presence of an external magnetic field B pointing in the z-direction, a Zeeman term is supplemented in Eq. (2), with g and µ B being the gyromagnetic factor and the Bohr magneton, respectively. In view of further discussions we single out the transverse anisotropy term and write H M in Eq. (2) as The reason behind this separation is that, contrary to H M , H 0 M has an extra abelian U(1) symmetry generated by the z-component of the total spinŜ t z =Ŝ z +ŝ z . Consequently, the eigenstates of H 0 M can be organized according to the eigenvalues m ofŜ t z [39]. On the other hand the transverse anisotropy term does not commute withŜ t z and needs to be treated separately. We shall discuss in more detail the eigenstates of H 0 M in section 2.2.
The source and drain electrodes are modeled as spin-polarized one-dimensional discrete chains which in the momentum space representation are described by the Hamiltonians: Both leads present an energy dispersion law of the form ε qσ = 2τ cos q + ∆ σ , with τ the effective hopping between the nearest neighbor sites in the leads and ∆ σ a rigidband spin splitting that describes the polarization of the leads. In Eq. (4) a † qασ creates an electron with momentum q and spin σ in the lead α = {L, R}.
The last term in Eq. (1) describes the hybridization of the SMM with the contacts where a † σ creates an electron with spin σ on the LUMO orbital and V α σ is the hopping amplitude of a tunneling process between the LUMO and the majority (σ = +) and minority (σ = −) electron states in the lead α. The coupling of the SMM to the contacts in the case of collinear magnetic configuration, and in the absence of the switching protocol (χ α (t) = 1), is described by Γ α σ = 2π|V α σ | 2 ̺ ασ (0), where ̺ ασ (0) is the spin density of states at the Fermi surface for electrons in lead α.
In view of further investigations we allow tunable spin polarizations in the leads and define P α : . Note that the Hamiltonian H T (t) contains two time-dependent dimensionless functions χ α (t) (α = L, R) which simulate the switching of the contacts between the molecule and the leads. As we are interested in the turnstile pumping it is enough to consider them simply as rectangular periodic pulses (see Fig. 1)

Energy eigenstates
In this section we discuss the energy spectrum and the organization of the eigenstates of the SMM Hamiltonian H M introduced in Eq. (2). We shall start by discussing first the spectrum of H 0 M . When J = 0, the LUMO orbital gets decoupled from the local spin and the Hamiltonian H 0 M has three U(1) symmetries generated by the local charge Q accumulated on the LUMO and by the z-components of the LUMO and local spins,ŝ z andŜ z . Consequently, {Q,ŝ z ,Ŝ z } provide the quantum numbers according to which the multiplets of the Hamiltonian are classified. Notice that this low symmetry classification is valid for finite magnetic fields. When B = 0, the Hamiltonian H 0 M has a much higher symmetry, i.e. U Q (1) × SUŝ(2) × SUŜ(2) in the charge and spin sectors, but this situation shall not be discussed here, as we always assume a finite magnetic field. In the case we consider here, the classification of the states is rather trivial and we can simply denote the eigenstates as follows: |0, 0, S z , |1, ↑, S z , |1, ↓, S z , and |2, 0, S z , with S z = −S, −S + 1, .., S. In the presence of Coulomb interaction the double occupied states |2, 0, S z have an energy of the order ∼ U, which is the largest energy scale in the problem, and in view of the discussion that follows, shall not contribute to transport. Therefore, to simplify the notations, it is enough to relabel the states and keep track of theŝ z andŜ z quantum numbers. In this new notation we have |0, 0, S z → |0, S z and |1, σ, S z → |σ, S z .
A finite exchange coupling J breaks the three U(1) symmetries down to U Q (1) × UŜ t z (1) generated by LUMO charge Q and the z-component of the total spinŜ t z . Still, the Hamiltonian H 0 M can be diagonalized exactly and the states constructed in an analytical fashion in terms of the states introduced previously for J = 0, by using the Clebsch-Gordan construction [26]. Now the new states |Q, m shall be classified by the molecular charge Q, and by the z-component of the total spin, m.
For m ∈ [−S + 1/2, S − 1/2] the single-particle states (Q = 1) are given by: and their associated eigenvalues E ± 1,m read as: where The coefficients C ± mσ in Eq (6) are simply the Clebsch-Gordan coefficients. The states |0, S z are not affected by the exchange coupling and one has |0, m = |0, S z . The corresponding eigenvalue is simply The remaining Q = 1 states are |1, −S − 1/2 and |1, S + 1/2 . For a vanishing magnetic field, B = 0, the states associated to ±m are degenerate and one has So far we have discussed how to construct the states and to compute the energy spectrum for H 0 M . In the rest of this paragraph we shall address the role of the transverse anisotropy term. The transverse anisotropy term ∼ (Ŝ 2 x −Ŝ 2 y ) ∼ (Ŝ 2 + +Ŝ 2 − ) does not commute withŜ t z and induces transitions [3] between the states of H 0 M with the selection rule |m − m ′ | = 2.
As the molecular charge is a good quantum number even in the presence of the transverse anisotropy, the eigenstates of the total molecular Hamiltonian H M can be classified by the molecular charge Q only. We shall label them |ϕ Q,ν , where Q = {0, 1} (states with molecular charge Q = 2 are disregarded) while ν is an internal label that indexes the states within a multiplet. In the presence of the transverse anisotropy E, the 'empty' molecular states (EMS) can be written as: with m running over all allowed values in the range [−S, S]. For half integer S and a small magnetic field, the transverse anisotropy plays a minor role in the mixing of the states |ϕ 0,ν , as the transition amplitudes between the empty molecular states {|0, m } are negligible (see also the discussion following Fig. 2.) In contrast, the transverse anisotropy couples the degenerate, one-particle states (Q = 1) with opposite m's. The strongest mixing is expected for the pairs |1, 1 ± and |1, −1 ± as the off-diagonal matrix element ± 1, 1|H M |1, −1 ± is linear in E. Higher order mixing effects become important as the ratio E/D increases and one can generally write: The eigenvalues E Q,ν of H M and the coefficients in Eqs. (10) and (11) can be found only by numerical diagonalization. More details on the spectral properties and state mixing will be given in Section 3. Finally, we write down the matrix elements of the tunneling Hamiltonian H T with respect to the eigenstates of H 0 M and derive the selection rules for molecular transitions due to electronic back-and-forth processes namely (λ ↑ = 1, λ ↓ = −1): This equation describes the tunneling of one electron with spin σ on the SMM orbital, when the number of electrons in the molecule increases by one while the total magnetic quantum number m changes by ±1/2.

Generalized master equation approach
The GME approach which we use to investigate the time-dependent transport relies on the partitioning approach [37]. More precisely, transient currents develop in the source and drain electrodes as they are contacted to the molecule at some initial instant. The leads are viewed as non-interacting particle reservoirs with chemical potentials µ α (α = L, R), and at equilibrium described by the Hamiltonian (4). This setting is suitable for perturbative calculations with respect to the lead-molecule couplings Γ α ± , and allows us to compute transient currents in the presence of time-dependent modulation of the contacts, as in the turnstile regime. The GME method essentially provides the SMM reduced density operator (RDO) ρ defined as the partial trace over the leads' degrees of freedom ρ(t) = Tr F el {W(t)}. Here W (t) is the density operator of the whole structure which solves the Liouville -von , and the trace is over F el -the Fock space of the non-interacting electronic reservoirs. In the sequential tunneling regime considered here the master equation takes a rather compact form (for a full derivation see Ref. [34]): where we introduced the "free" evolution operator of the disconnected system U t = e −i(H M +H L +H R )t/ and the equilibrium distribution of the leads ρ el [40]. The dissipative operator D t collects all sequential tunneling processes from the switching instant t 0 = 0 to the current time t. We solve numerically Eq. (13) with respect to the fully interacting states of H M and obtain the populations associated to a given state |ϕ Q,ν as Once ρ(t) is known one can calculate average values of molecular observables by performing a trace over the Fock space F M of the molecule. For instance, the total charge accumulated on the orbital involved in transport is given by where the total electronic occupationN =n ↑ +n ↓ and e is the electron charge. The continuity equation then becomes [35] where ν Q is the set of states with charge Q . By inserting the tunneling Hamiltonian H T into the double commutator D t given in Eq. (14) one identifies J L and J R from the RHS of Eq. (17). It is straightforward to show that the 'empty' molecular states |ϕ 0,ν do not contribute to the currents. Similarly, one can calculate the total spin S t z = Tr F M {ρ(t)Ŝ t z } as well as the spin currents. In this work the relaxation of the excited molecular states via phonon emission is not considered. This is a good approximation as long as the timescale on which the quantum turnstile operates is much smaller than the relaxation time which is of order of 10 −6 s (see e.g Refs. [41,42]). In fact previous studies [43] reported that the current-induced magnetic switching is stable against intrinsic spin-relaxation processes.

The transport configuration and tunneling processes
The numerical simulations were performed for molecules with spin S = 5/2 but our conclusions remain valid for larger half-integer values of S. For reasons that will become clear below in this work we restrict ourselves to SMMs with small transverse anisotropy, that is E ≪ D ≪ J (in fact we allow a maximum ratio E/D = 1/25 at fixed easy-axis anisotropy constant D). We shall investigate two spin configurations for the leads. In the so called antiparallel (AP) configuration the left lead carries only spin-down electrons and the right lead is spin-up polarized. For the second configuration the left lead is non-magnetic (i.e. P L = 0) and the right lead remains ferromagnetic. We label this configuration as normal-ferromagnetic (NF).
The chemical potential of the leads µ L,R are set such that only the states |ϕ 0,ν and |ϕ 1,ν contribute to the tunneling processes, while double occupied states, |ϕ 2,ν , have a much higher energies, E 2,ν > µ L , and do not contribute to transport.
Let us first discuss the energy spectrum and the relevant lead-molecule tunneling processes in the absence of transverse anisotropy. In the following discussion we shall use the basis {|Q, m } of H 0 M = H M (E = 0). As we are interested in the pumping mechanism at finite transverse anisotropy E, we shall switch later to the basis |ϕ Q,ν of the full H M . In that situation, with certain modifications, a similar turnstile scenario holds. Figure 2(a) shows the energy levels of H 0 M for a given set of parameters, and in the absence of the external magnetic field. Figure 2(b) schematically shows the charge Q = 1 (integer m's) and Q = 0 (half integer m's) branches of the spectrum; in view of further discussion the double arrow and the dotted lines mark some of the quantum tunneling of magnetization (QTM) processes induced by a nonvanishing E. The states connected by the double arrow are strongly mixed by E, while the ones connected by the dashed lines are only weakly coupled. As charge Q is conserved, it implies that direct transitions between Q = 0 and Q = 1 branches are forbidden by symmetry. This is only possible through processes involving states in the leads that do not conserve the charge. For example in Fig. 2(c), we show such processes (blue and red arrows). The same figure also shows how the SMM evolves from an initial 'empty' molecular state |0, 5/2 to the to the next EMS |0, 3/2 via tunneling processes. We shall call the transitions m → m − 1/2 'forward' processes (they follow the full line arrows in Fig. 2(c)) as they contribute to the magnetic switching m = 5/2 → m = −5/2. On the other hand the transitions m → m + 1/2 compete for the total spin reversal and can be regarded as 'backward' processes (they follow the dashed lines in Fig. 2(c)).
Furthermore, we distinguish two types of 'forward' transitions: (i) 'absorption' of spin-down electrons from the leads, i.e. the charging of the molecular orbital along the transitions |0, m → |1, m − 1/2 ± (full red arrow in Fig. 2(c)) and (ii) tunneling of spin-up electrons from the molecular orbital, i.e. a depletion process associated to the transitions |1, m − 1/2 ± → |0, m − 1 (full blue arrow). Similarly one defines charging and discharging 'backward' processes (associated with the dashed lines in Fig. 2(c)). We find this analysis useful as it provides hints for a write-and-read scheme of states with well defined molecular spin |0, m when operating the SMM in the turnstile regime. Such a protocol will be discussed in the next subsection.

The turnstile protocol
A turnstile pumping cycle entails two steps: (i) the charging of the molecular orbital from the left lead while the drain contact is closed (χ L (t) = 0, χ R (t) = 0) and (ii) discharging/depletion through the drain lead (χ L (t) = 0, χ R (t) = 0). We simulate the turnstile operation by appropriately tailoring the switching functions χ L and χ R in the tunneling Hamiltonian H T .
The main idea behind the proposed operation is the following: use the charging sequences to prepare intermediate one-electron states via 'forward' tunneling from the left lead and write 'empty' molecular states |0, m along 'forward' depletions to the right lead. To be more precise, let us discuss a single turnstile cycle at E = 0 in the NF configuration, starting from the initial state |0, 5/2 . The associated transitions are depicted in Fig. 2(c). By opening the source (left) contact, the states |1, 2 ± become populated by absorbing one spin -↓ electron (forward tunneling). In turn, if a spin -↑ electron is absorbed then the rightmost state |1, 3 becomes activated (backward tunneling). Then, the left contact closes and the drain (right) electrode comes into play. Now, the orbital is depleted through forward tunneling |1, 3 → |0, 5/2 and |1, 2 ± → |0, 3/2 , as spin -↑ electron tunnels out into the right lead. An accurate operation would lead to the preparation of a single EMS or to an equally weighted mixture of EMS, but this scenario is not expected to work if the transverse anisotropy induces strong mixing of states |0, m .
In view of this analysis, let us now discuss how this picture gets modified in the presence of the transverse anisotropy. We start by describing the construction of the empty molecular states {|ϕ 0,ν }. We find that if gµ B B ≪ D the mixing of 'empty' molecular states |0, m is negligible since the QTM between the states |0, ∓1/2 and |0, ±3/2 is very weak (see the dotted lines in Fig. 2(b)). In this case we find a one to one correspondence between ν ↔ m as for any ν in Eq.(10) one can find a single m such that |ϕ 0,ν ≈ |0, m .
This simple correspondence fails as the magnetic field approaches resonant value gµ B B res = −D(S z + S ′ z ), and the Landau-Zenner tunneling processes between E 0,m and E 0,m ′ =m±2 become important and lead to strong mixing of the states. Such a resonant regime will not be discussed in the present work.
We now turn to Q = 1 states. For B = 0, the states |1, m ± and |1, −m ± in Eq. (11) are mixed by the transverse anisotropy term as E ± 1,m and E ± 1,−m are degenerate (see Eq. (8)). The mixing is indicated by the double arrow in Fig. 2(b). However, even a small magnetic field lifts this degeneracy and one finds a rather small mixing of the states |1, m ± and |1, −m ± for E = 0. Once again, for each |ϕ 1,ν there is a single state |1, m ′ s of H 0 M whose weight |c s ν,m ′ | 2 in Eq.(11) is by far the largest one. Under these conditions the one-to-one correspondence between the index ν and the quantum number m is preserved for all states and allows us to index them as: We shall close this section by noticing that although this representation works, ν will always be read as an index, and not as a quantum number. Consequently, in the |ϕ s 1,m basis the populations of the states will be denoted by P |ϕ s 1,m and P |ϕ 0,m .
The writing of a single EMS depends crucially on the spin polarization of the leads. Fig. 3(c) indicates that the 'backward' processes (both charging and relaxation) are forbidden, as the state |ϕ 1,3 is not available along the charging sequence because the left lead provides no spin-up electrons. For the same reason there are no transitions from |ϕ ± 1,2 back to |ϕ 0,5/2 on the discharging sequence. In fact, on the charging cycle the system occupies only two states with imbalanced populations P |ϕ + 1,2 > P |ϕ − 1,2 and which eventually deplete in favor of |ϕ 0,3/2 . It is important to observe that on the depletion cycle the average total spin S t z = 3/2 and coincides with the molecular spin of the |ϕ 0,3/2 EMS (see Fig. 4(a)).
The reverse magnetic switching can be also implemented by simply reversing the bias (µ L ↔ µ R ) while keeping both contacts closed and then repeating the turnstile operation with the new initial state |ϕ 0,3/2 . Then the system returns to |ϕ 0,5/2 . This is a classical NOT operation, as the system evolves from |ϕ 0,5/2 to |ϕ 0,3/2 and then back to |ϕ 0,5/2 without passing through a superposition of these states.
Let us emphasize that the preparation of a pure excited molecular state |ϕ 0,m cannot be achieved in the standard transport regime. In that case the charge flows simultaneously to and from the leads and one cannot completely deplete the molecule and therefore ϕ 0,m |ρ(t)|ϕ 0,m < 1. It should be mentioned that for larger S the time needed to achieve the full magnetic switching m = S → −S also increases as the system must visit all the intermediate states [26]. This fact suggests that the pair of consecutive states (|ϕ 0,S , |ϕ 0,S−1 ) might be more appropriate for faster manipulation of magnetic qubits.
Given these results one can ask about the time evolution of the total spin under repeated pumping cycles and on the possibility to read the states prepared along the turnstile operation by measuring currents. Fig. 4 summarizes our main results on transient currents and spin evolution along few turnstile cycles for the AP and NF configurations. The time-dependent occupation of the molecular orbital (the blue line in Fig. 4(a)) has a typical charging/relaxation pattern, with quick orbital filling and slightly slower depletion. This can be seen by comparing the abrupt increase (in less than 1/2 ns) of the population at the beginning of the charging cycles (e.g. t = 4, 8, 12ns) to the smooth tail of discharging which extends over 1ns (e.g the time range [6,7]ns).
The total spin average S t z presented in Fig. 4(a) displays a step-like structure in both configurations. The steps scan both integer and half-integer values of S t z , the last step for the AP configuration corresponding to S t z = 3 being reached after t ≃ 18 ns (not shown). In the AP configuration the onset of half-integer steps corresponds to the depletion of the molecular orbital (Q = 1 → Q = 0), whereas the transition between half-integer to integer steps is associated to the charging process (Q = 0 → Q = 1).
The NF configuration presents different features: there are fewer but longer steps of S t z , but no clear correspondence can be made between these steps and the behavior of the total charge Q. One notices that in this case the integer steps extend on some charging sequences and that half-integer values are encountered even if the orbital is empty. Fig. 4(b) displays the expected series of spikes for the transient currents J L,R in the AP configuration. The period of the pumping cycles must be chosen appropriately in order to ensure full charging and discharging of the molecular orbital (we find the minimal period to be ∼ 1 ns).
We have found a similar behavior (not shown here) for the transient currents in the NF configuration. The input current J L vanishes when the orbital is fully occupied (Q = 1), whereas on the discharging sequence J R drops to zero as the orbital depletes. The amplitudes of J L and J R are different because the charging process is faster than the depletion (see Fig. 4(a)). By inspecting Figs. 4(a) and (b) one observes that in the AP configuration we have a one to one correspondence between the average S t z and the  peak-to-peak sequence in the transient currents: S t z acquires half-integer values only between a depletion peak and the next charging peak (e.g. for t ∈ [2, 4]ns S t z = 3/2), while between a charging peak and the next depletion peak the average spin is an integer. This means that the AP configuration can be used to record experimentally the initialization of a given 'empty' molecular state |ϕ 0,m .
To this end it is sufficient to know the initial state of the molecule and to carefully 'count' the transient peaks of J L and J R . We need to keep in mind though that Fig. 4(a) shows the average value of the total spin, which does not guarantee that along half-integer steps of S t z the system is in a pure state characterized by the RDO ρ ∼ |0, m 0, m|, especially for larger values of the transverse anisotropy when one expects stronger mixing of states. To further investigate the role of the anisotropy, we calculated the populations P |ϕ 0,m of several 'empty' molecular states |ϕ 0,m for different values of the ratio E/D, at fixed magnetic field. Fig. 5(a) confirms that at E/D = 1/250 the k-th depletion cycle is well described by a single state |ϕ 0,m=5/2−k . This proves the stepwise all-electrical writing of EMS (i.e. point (i) in the Introduction).

Transverse ansiotropy effects
By increasing the transverse anisotropy constant such that E/D = 1/75 we notice in Fig. 5(b) the emergence of a 2nd EMS on the depletion cycles. Nevertheless, the population of the dominant 'empty' molecular state exceeds 0.9 so we can still associate a well defined molecular state to each of the depletion cycle. This no longer holds for E/D = 1/25. Fig. 5(c) reveals that the weight of the state |ϕ 0,−3/2 and |ϕ 0,−5/2 on the 2nd and 3rd depletion cycle increases up to 0.25, reducing the efficiency of the quantum turnstile protocol. Moreover, one can easily see that along the 4th cycle (i.e t ∈ [14,16]ns) P |ϕ 0,−3/2 + P |ϕ 0,1/2 < 1 which suggests that other states have to be populated. We have found that the state |ϕ 1,−3 , corresponding to Q = 1 gets populated after the third cycle due to the forward transition |ϕ 0,−5/2 → |ϕ 1,−3 via spin -↓ tunneling into the SMM. In order to explain the coexistence of two EMSs on the same depletion sequence when the transverse anisotropy increases we have to analyze the QTM between nearly degenerate Q = 1 states. By looking at the off-diagonal matrix element ± 1, 1|(Ŝ 2 + +Ŝ 2 − )|1, −1 ± we infer that by increasing E the hybridization of |1, 1 ± and |1, −1 ± in the fully interacting states increases as well. We find that the weight of the 'minority' state |1, −1 ± in |ϕ ± 1,1 increases from ∼ 10 −2 for E/D = 1/250 to ∼ 10 −1 for E/D = 1/25. The mixing of the states |1, 2 ± and |1, −2 ± arises to the second order in E and is still negligible. As a consequence the accuracy of the first turnstile cycle is preserved even for E/D as large as E/D ≃ 1/25 and that P |ϕ 0,3/2 ≈ 1. This is confirmed by the results presented in Fig. 5(a)-(c). In order to recover 'clean' EMSs on each depletion cycle for larger values of E/D one could slightly increase the magnetic field. The latter lifts even more the degeneracy of the states |1, 1 ± and |1, −1 ± and reduces therefore their mixing in the presence of E/D. are less responsive to charging (this processes correspond to the dashed line in Fig. 6(b)) because spin down tunneling is allowed only through the states |1, 1 ± whose weights are small. By similar arguments, one can see that the discharging process activates two EMSs, namely |ϕ 0,1/2 and |ϕ 0,−3/2 .
Both of these states acquire important weights in the RDO for E/D ≃ 1/25 so the average total spin can no longer be associated to a well defined value of the molecular spin. We therefore conclude that the enhanced QTM between Q = 1 states damages the efficiency of the turnstile protocol even if the EMSs involved in transport are not mixed.
Let us note that the possibility to prepare a single EMS is not obvious as an open system is generally described by a mixed state. Our simulations also show that in the quantum turnstile regime one controls the transitions between any pair of intermediate molecular states (|ϕ 0,m , |ϕ 0,m−1 ) along a pumping cycle, in contrast to the full magnetic switching which involves only the pair (|ϕ 0,S , |ϕ 0,−S ).
Finally, we mention that if the SMM has an integer spin one cannot expect an accurate turnstile operation because the transverse anisotropy induces strong mixing between quasidegenerate EMSs (e.g. between |ϕ 0,1 and |ϕ 0,−1 ).

Conclusions
In the present work we address the transient transport regime and turnstile pumping across a single-molecule magnet coupled to external leads. The time-dependent evolution of the molecular states has been discussed in detail and signatures of the electrically induced magnetic switching on the transient currents were predicted. For ferromagnetic leads with antiparallel spin polarizations the turnstile protocol allows the stepwise writing and reading of excited molecular states.
The evolution of the states along the turnstile operation can be 'read' indirectly from the behavior of the transient currents. More precisely, by recording the charging and discharging currents one can monitor the evolution of the system and identify the regimes where its density matrix is described by a single empty molecular state. This is somehow in contrast to the situation when the leads are simply normal metals and where the control of the excited spin states cannot be achieved as the molecular spin is reversed continuously.
We show that the transverse anisotropy leads to the hybridization of nearly degenerate one-particle states which subsequently relax to empty molecular states with different values of the total spin. However, this dephasing effect can be reduced by applying a moderate perpendicular magnetic field.
Another useful application of the turnstile regime that we address here is the possibility to mix several excited spin states (viewed as magnetic qubits) during the discharging cycles when the source electrode is normal and the drain electrode is ferromagnetic. Note that the short rise time of the switching functions used in our simulations is not essential for the turnstile operation and slower switching functions could be in principle selected to achieve a better resolution of the transient peaks. Our predictive simulations clearly emphasize the potential of the molecular turnstiles as promising candidates for molecular spintronics. As a method approach we have used the generalized master equation formalism adapted to the turnstile configuration.