Effects of spin–orbit coupling and many-body correlations in STM transport through copper phthalocyanine

  1. ,
  2. and
Institut für Theoretische Physik, Universität Regensburg, D-93040, Germany
  1. Corresponding author email
Guest Editor: J. M. van Ruitenbeek
Beilstein J. Nanotechnol. 2015, 6, 2452–2462. https://doi.org/10.3762/bjnano.6.254
Received 14 Aug 2015, Accepted 27 Nov 2015, Published 22 Dec 2015
Full Research Paper
cc by logo


The interplay of exchange correlations and spin–orbit interaction (SOI) on the many-body spectrum of a copper phtalocyanine (CuPc) molecule and their signatures in transport are investigated. We first derive a minimal model Hamiltonian in a basis of frontier orbitals that is able to reproduce experimentally observed singlet–triplet splittings. In a second step SOI effects are included perturbatively. Major consequences of the SOI are the splitting of former degenerate levels and a magnetic anisotropy, which can be captured by an effective low-energy spin Hamiltonian. We show that scanning tunneling microscopy-based magnetoconductance measurements can yield clear signatures of both these SOI-induced effects.


Spin–orbit interaction (SOI) can play a major role in molecular spintronics. For example, in combination with the configuration of the non-magnetic component (organic ligand), it is known to be essential in establishing magnetic anisotropy in high-spin molecular magnets [1], and it is quite generally expected in metalorganic compounds. Effective spin-Hamiltonians are commonly used to describe this anisotropy, and usually capture well the low energy properties of these systems [1-3]. Such effective Hamiltonians have been derived microscopically for widely studied molecular magnets such as Fe8, Fe4 and Mn12 [4]. Recently, magnetic anisotropy effects could be directly probed by magnetotransport spectroscopy for Fe4 in quantum-dot setups [5,6]. An interesting question is hence if other classes of metallorganic compounds, such as the widely studied metal phthalocyanines [7,8], exhibit sizeable magnetic anisotropy induced by the interplay of electronic correlations and SOI. Indeed, in an X-ray magnetic circular dichroism (XMCD) analysis copper phthalocyanine (CuPc) was found to exhibit enormous anisotropies in both spin and orbital dipole moments [9]. Furthermore, recent experimental findings for cobalt pththalocyanine in a scanning tunneling microscopy (STM) setup [10] suggest that many-body correlations play an important role in the interpretation of the transport measurements. In a recent work [11], we have explicitly investigated long-range and short-range electron–electron correlation effects in CuPc and found a singlet–triplet splitting of the former anionic groundstate of about 18 meV, and thus a triplet as anionic ground state.

In this work we add the SOI to our analysis. We find that it further removes the triplet degeneracy by inducing splittings of few tenths of millielectronvolts. Moreover, in combination with exchange correlations, it produces a magnetic anisotropy which can in turn be captured by an effective spin Hamiltonian.

In general, the accurate calculation of the many-body properties of metallorganic molecules, such as molecular magnets or our CuPc, is a highly nontrivial task. In fact, the number of their atomic constituents is large enough that exact diagonalization is not possible and standard density-functional schemes have difficulties in capturing short ranged electron–electron correlations [4]. In order to reduce the size of the many-body Fock space, we use a basis of frontier molecular orbitals as the starting point to include electronic correlations [11,12] and construct a generalized Hubbard Hamiltonian. Furthermore, the symmetry of the molecule greatly helps to reduce the number of matrix elements one has to calculate in this basis.

To probe both SOI-induced splittings and magnetic anisotropy, we further investigated the current characteristics of a CuPc molecule in an STM configuration similar to the experiments in [13,14]: The molecule is put on a thin insulating layer grown on top of a conducting substrate. The layer functions as a tunneling barrier and decouples the molecule from the substrate. Hence the CuPc molecule acts as a molecular quantum dot weakly coupled by tunneling barriers to metallic leads (here the STM tip and the substrate). This quantum dot configuration should be favourable to experimentally probe SOI splittings and magnetic anisotropies when an external magnetic field is applied to the system, in analogy to the experiments in [6]. Indeed, we demonstrate that experimentally resolvable SOI splitting should be observed at magnetic fields of a few teslas.

The paper is organized as follows: We first derive a microscopic Hamiltonian for CuPc in the frontier orbital basis which includes exchange correlations and the SOI. This Hamiltonian is diagonalized exactly and used in further spectral analysis and transport calculations. Its spectrum is also used to benchmark the prediction of an effective spin Hamiltonian that captures well the low-energy properties of CuPc both in its neutral and anionic configurations. Finally, transport calculations with and without magnetic fields are presented and SOI-induced signatures are analyzed.

Results and Discussion

Microscopic model Hamiltonian for CuPc

The focus of this section is the establishment of a minimal model Hamiltonian for an isolated CuPc molecule capable to account for both electron–electron interaction and spin–orbit coupling effects. As discussed below, parameters are fixed such that experimental observations for the singlet–triplet splitting [8] as well as positions of anionic and cationic resonances [14] are satisfactorily reproduced. In its most general form and for a generic molecule such Hamiltonian reads


where the single-particle Hamiltonian of the molecule is given by [Graphic 1], [Graphic 2] describes electronic interactions and [Graphic 3] accounts for the spin–orbit interaction (SOI).

Single-particle Hamiltonian for CuPc

The one-body Hamiltonian [Graphic 4], written in the atomic basis [Graphic 5], reads


where α is a multi-index combining atomic species and orbital quantum number at position rα, see Figure 1a. For the ligand we consider the set of all 2s (1s for hydrogen), 2px and 2py orbitals as the σ-system, and consequently the set of 2pz orbitals as the π-system. On the metal, the 3dxy, [Graphic 6], [Graphic 7] and 4s orbitals contribute to the σ-system, while the 3dzx and 3dyz belong to the π-system. This basis yields a total of 195 valence electrons for neutral CuPc. Atomic on-site energies εα and geometrical parameters were taken from [7,15]. The hopping matrix elements bαβ in Equation 2 are obtained by using the Slater–Koster [16] and Harrison [17] LCAO schemes, similar to [18]. Numerical diagonalization of [Graphic 8] finally yields single particle energies εi, see Figure 1b, and molecular orbitals [Graphic 9], cf. Supporting Information File 1.


Figure 1: (a) Geometry and atomic composition of CuPc. (b) Single particle energies of relevant molecular orbitals. Black (grey) circles depict the π (σ) character of the corresponding orbital. The color (diameter) of the inner circles characterizes the type (weight) of the metal orbital contribution on the corresponding molecular orbital. (c) Depiction of the four frontier orbitals retained in this work: SOMO (S), HOMO (H) and LUMOzx/yz (Lzx/yz).

Stemming from Hartree–Fock calculations for isolated atoms [15], the atomic on-site energies εα do not take into account the ionic background of the molecule and crystal field contributions. Therefore, molecular orbital energies εi have to be renormalized with parameters δi to counteract this shortage, yielding (cf. Supporting Information File 2)


In this work we use a constant shift δi = δ = 1.83 eV.

Due to the odd number of valence electrons, in its neutral configuration CuPc has a singly occupied molecular orbital (SOMO). When the molecule is in its anionic groundstate, this orbital does not become doubly occupied [7]. Hence, the orbitals most relevant for transport (frontier orbitals) are the SOMO (S), the HOMO (H) and the two degenerate LUMOs (Lzx/yz), which transform according to the b1g, a1u and eg irreducible representations of the point group of CuPc (D4h), respectively. They are depicted in Figure 1c. The LUMO orbitals in their real-valued representations, [Graphic 10] and [Graphic 11], have equal contributions cL ≈ 0.097 on both 3dzx and 3dyz orbitals of the metal. Due to their degeneracy, they can be transformed into their complex, rotational invariant representations:


where [Graphic 12] is the n = 3 metal orbital with angular momentum [Graphic 13] and magnetic quantum number m = ±1. To distinguish contributions from the pure phthalocyanine (Pc) ligand and the copper (Cu) center, we introduced [Graphic 14] and [Graphic 15], respectively. Likewise, with cS ≈ 0.90, we can write for the SOMO:


where [Graphic 16] is the n = 3 metal orbital with angular momentum [Graphic 17] and projection m = ± 2 onto the z-axis. Finally, the HOMO has no metal contributions and thus we have trivially [Graphic 18]. The representations introduced in Equation 4 have the advantage that the four frontier orbitals can then be characterized by the phases φi acquired under rotations of π/2 around the main molecular symmetry axis. For the SOMO φS = π, for the HOMO φH = 0 and for the two LUMOs φL± = ±π/2.

Many-body Hamiltonian in the frontier orbitals basis

In order to set up a minimal many-body Hamiltonian, we restrict the full Fock space to many-body states spanned by the SOMO (S), the HOMO (H) and the two LUMO (L±) orbitals and write Equation 1 in this basis. Hence, for neutral CuPc the number of electrons populating the frontier orbitals is N0 = 3.

We exploit the distinct phases acquired by the frontier orbitals under 90° rotations to determine selection rules for the matrix elements Vijkl in [Graphic 19],


namely [Graphic 20] if [Graphic 21], [Graphic 22], cf. Supporting Information File 2. Equation 6 in this basis then reads


where the indices i, j, k now run over the set of frontier orbitals, and the notation [ijk] means that the sum runs only over unlike indices, i.e., i, j and k are different from each other in the corresponding sum. The abbreviations we introduced in Equation 7 are the orbital Coulomb interaction Ui = Viiii, the inter-orbital Coulomb interaction Uij = Viijj, the exchange integral [Graphic 23], the ordinary pair hopping term [Graphic 24] and the split pair hopping term [Graphic 25]. Contributions with four different indices are found to be very small (of the order of microelectronvolts) and thus omitted in this work. The matrix elements Vijkl are calculated numerically using Monte Carlo integration [19] and renormalized with a dielectric constant εr = 2.2 in order to account for screening by frozen orbitals [12]. A table with the numerically evaluated interaction constants is found in Supporting Information File 2.

Spin–orbit interaction (SOI) in the frontier orbitals basis

A perturbative contribution to the bare one-body Hamiltonian [Graphic 26] relevant in molecular systems is provided by the SOI. In the following we derive an effective spin–orbit coupling operator acting on the subset of frontier orbitals. The atomic SOI operator reads


where α and [Graphic 27] run over all atoms and shells, respectively. By evaluating Equation 8 only on the central copper atom, i.e., [Graphic 28] and α = Cu, [Graphic 29] in second quantization is given by


where [Graphic 30] creates an electron with spin σ on the copper atom in the orbital specified by [Graphic 31]. For an electron in the 3d shell of Cu we use ξCu ≈ 100 meV [20]. Projecting Equation 9 onto the minimal set of frontier orbitals then yields:


where λ1 = 1/2ξCu|cL|2 = 0.47 meV and λ2 = ξCu(cScL)/[Graphic 32] = 6.16 meV are now effective spin–orbit coupling constants. A similar analysis of SOI in CuPc, laying more focus on the central Cu atom, can be found in [21].

Finally, many-body eigenenergies ENk and eigenstates [Graphic 33], labeled after particle number N and state index k, are obtained by exact numerical diagonalization of [Graphic 34] in the frontier orbitals basis. Despite numerically tractable, the problem described by [Graphic 35] is still highly intricate, as the Fock space has dimension 44 = 256. In reality, though, only few low-lying many-body states are relevant at low energies. This enables further simplification and even an analytical treatment, as discussed in the next subsection.

Low-energy spectrum of CuPc and effective spin Hamiltonian

In the following we will analyze the neutral and anionic low-energy part of the many-body spectrum of CuPc and establish an effective Hamiltonian which enables us to analyze the low-energy behaviour in a more lucid way. To this extent, we start by observing that [Graphic 36] (in the considered particle number subblocks) contains different energy scales, in particular, U > J > λ, which suggests a hierarchy of steps. We use U, J and λ to denote the set of all Hubbard-like parameters (Ui, Uij), all exchange parameters ([Graphic 37]) and all SOI parameters (λi), respectively. As a first step we set both the exchange (J) and SOI (λ) contributions to [Graphic 38] to zero and determine the neutral and anionic groundstates. In a second and third step exchange and SOI are added, respectively.

Neutral low-energy spectrum

In the neutral low-energy part of the spectrum, we retain the two spin-degenerate groundstates of [Graphic 39],


with corresponding energy [Graphic 40]. We defined [Graphic 41]. The groundstates in Equation 11 are neither affected by [Graphic 42] nor by the exchange terms in Equation 7. Trivially, the effective Hamiltonian in the basis of [Graphic 43] reads:


In principle Equation 7 also contains terms that act on the neutral groundstate, such as for example pair hopping terms proportional to [Graphic 44], and cause admixtures with other many-body states. However, according to our full numerical calculations, these admixtures are rather small and do not affect transitions between neutral and anionic states.

Anionic low-energy spectrum

Continuing with the anionic low-energy part of the spectrum of [Graphic 45], we find an eightfold degenerate groundstate:


with corresponding energy [Graphic 46]. The eightfold degeneracy comes from the two unpaired spins in either SOMO or LUMO and the orbital degeneracy of the LUMO orbitals. In order to make the anionic eigenstates also eigenstates of the spin operators [Graphic 47] and [Graphic 48], they can be rewritten as


The orbital degeneracy of the LUMOs, expressed by the index τ, is responsible for the two sets of singlets (total spin S = 0) and triplets (total spin S = 1). Considering exchange interaction in a second step, we find that only the [Graphic 49] term in Equation 7,


directly determines the low-energy structure of the anionic low-energy part because of the singly occupied SOMO and LUMOs: The degeneracy between singlets and triplets is lifted, see Figure 2, and we obtain


for the singlets and triplets, respectively.


Figure 2: Lowest lying anionic states of CuPc, together with their grade of degeneracy d. Without exchange and SOI, the anionic groundstate is eightfold degenerate. When exchange interaction between SOMO and LUMOs is introduced, the degeneracy is lifted, yielding two triplets and two singlets because of the orbital degeneracy of the LUMO. SOI further splits the triplet states, generating a twofold degenerate anionic groundstate consisting of the states [Graphic 50] and [Graphic 51].

Finally, to analyze in a third step how [Graphic 52] affects the low-energy part of the anionic part of the spectrum, in particular which degeneracies are lifted, we treat it as a perturbation and apply second-order perturbation theory to obtain the energy shifts. To this end, some additional states have to be considered. They are listed in Supporting Information File 3.

The states [Graphic 53] and [Graphic 54] experience a downshift due to [Graphic 55] and become the groundstates. Measuring energies with respect to [Graphic 56], we get


see Figure 2. Note that in our numerical calculations [Graphic 57] and [Graphic 58] are mixed and the degeneracy of the resulting states is lifted by a small shift in the range of some μeV. A more detailed discussion concerning the mixing of [Graphic 59] and [Graphic 60] can be found in Supporting Information File 3. The next states are [Graphic 61] and [Graphic 62] with


Due to their quadratic dependence on λ1 and λ2, these states change very little with [Graphic 63]. The degeneracy of the states [Graphic 64] and [Graphic 65] is lifted by the mixing of these states through [Graphic 66]. We find


where for [Graphic 67] we omitted smaller additional contributions from other states. The energies change according to


For further details we refer to Supporting Information File 3. Finally, the singlets S+ and S, similar to [Graphic 68] and [Graphic 69], change very little (with respect to [Graphic 70]):


By introducing [Graphic 71], an approximate Hamiltonian up to first order in [Graphic 72] can be given for the N0 + 1 particle subblock:


Equation 24 is one major result of this work. It shows that, similar to the well-studied molecular magnets [3-6], the interplay of spin–orbit coupling and exchange interactions yield magnetic anisotropies that can be captured by effective spin Hamiltonians. Noticeably, because Equation 24 was derived from the microscopic molecular Hamiltonian [Graphic 73], it was possible to check that deviations are in range of microelectronvolts and only of quantitative nature by comparison of the spectrum to the numerically evaluated one. Another source of magnetic anisotropy is the Dzyaloshinskii–Moriya interaction [22,23]. Although the latter is also linear with respect to the SOI, it does not appear in our model. The fundamental reason for neglecting it is the large ratio between the hopping integrals (of the order of electronvolts) and the SOI (ξCu ≈ 100 meV), which also justifies our perturbative analysis in terms of molecular orbitals. However, for molecular quantum dots with comparable SOI and hopping integrals the Dzyaloshinskii–Moriya interaction is sizeable and produces interesting effects on magnetization [24] and transport characteristics [25].

Interaction with magnetic fields

An experimentally accessible way to probe magnetic anisotropies is to apply external magnetic fields. In order to account for interactions of orbitals with magnetic fields, the atomic hopping matrix elements bαβ in Equation 2 have to be corrected with Peierls phase factors,


where, using the gauge [Graphic 74], the phase is given by


Here (xα, yα) are the in-plane atomic coordinates. Owing to the planar geometry of CuPc, [Graphic 75] depends only on the z-component Bz of the magnetic field B. In Figure 3 we show the dependence of the energies of the frontier molecular orbitals on the strength of the magnetic field in z-direction, Bz. For the two LUMOs we observe a linear dependence on the magnetic field, yielding an effective orbital moment of μorb = 33.7 μeVT−1, while the LUMO−(+) goes down (up) in energy with Bz (Figure 3a). The energies of the HOMO and the SOMO, however, scale quadratically with the magnetic field at a much lower scale (Figure 3b). This behaviour is expected, since the a1u and b1g representations have characters +1 under [Graphic 76] rotations, which transform Bz to −Bz. Thus the energies of HOMO and SOMO can not depend on the sign of Bz and must move at least quadratically with Bz. The two-dimensional eg representation on the other hand has zero character under [Graphic 77] rotations, which implies that the constituents of eg transform under such rotations either with different signs or into each other; indeed under a [Graphic 78] rotation LUMO+ is mapped onto LUMO− and vice versa.


Figure 3: (a) Dependence of the single particle orbital energies on the magnetic field strength. From this, the effective orbital moment of the LUMOs, here depicted in their complex representation, can be extracted as μorb = 33.7 μeVT−1. The energies of the SOMO and HOMO orbitals depend quadratically on the magnetic field and involve a much lower scale than the LUMOs, as seen in the close-up in panel (b).

Finally, the interaction of electronic spins with magnetic fields is represented by adding a Zeeman term [Graphic 79] to Equation 1,


where gS = 2 and S is the total spin operator on the molecule written in the frontier orbital basis.

Effective low-energy Hamiltonian

Putting everything together, an effective low-energy Hamiltonian including magnetic interaction terms for both orbital and spin degrees of freedom can thus be given. It reads


where [Graphic 80] is the Hamiltonian for the corresponding low-energy N-particle subblock as given by Equation 12 and Equation 24.

Dynamics and transport

Reduced density operator and current

The transport calculations for the molecule in an STM setup are done by using the formalism introduced in earlier works [18,26,27]. For the sake of clarity, in the following we briefly discuss the main steps to obtain the current through the molecule. The full system is described by the Hamiltonian


where [Graphic 81] describes the isolated molecule, see Equation 1. To incorporate image charge effects in our model, leading to renormalizations of the energies of the system’s charged states [28], we included a term [Graphic 82][11],


where [Graphic 83] is the particle number operator on the molecule. Electrostatic considerations regarding the geometry of the STM setup yielded δic ≈ 0.32 eV [11]. The Hamiltonians [Graphic 84] and [Graphic 85] corresponding to substrate (S) and tip (T), respectively, are describing noninteracting electronic leads. They read


where [Graphic 86] creates an electron in lead η with spin σ and momentum k. The tunneling Hamiltonian [Graphic 87] finally is given by


It contains the tunneling matrix elements [Graphic 88], which are obtained by calculating the overlap between the lead wavefunctions [Graphic 89] and the molecular orbitals [Graphic 90][26]. They yield the tunneling rates

[Graphic 91]

which are of the order of 10−6 eV and 10−9 eV for the substrate and the tip, respectively. Finally, the dynamics of the transport itself is calculated by evaluating the generalized master equation,


for the reduced density operator [26,29] ρred = TrS,T(ρ). The Liouvillian superoperator


contains the terms [Graphic 92] and [Graphic 93] describing tunneling from and to the substrate and the tip, respectively. To account for relaxation processes leading to de-excitation of molecular excited states, we included a relaxation term [Graphic 94], analogously to [30]:


It depends on the deviation of ρ from the thermal distribution ρth,N of the N-particle subblock, which is given by a Boltzmann distribution:


with β = (kBT)−1. Since [Graphic 95] acts separately on each N-particle subblock, it conserves the particle number on the molecule and thus does not contribute to transport directly. In this work, the relaxation factor 1/τ is around the same order of magnitude as the mean tip tunneling rate onto the molecule. In particular, we are interested in the stationary solution [Graphic 96] for which [Graphic 97]. Finally, the current through the system in the stationary limit can be evaluated as


yielding the current operator for lead η as [Graphic 98].

Transport characteristics

In this work, a tip–molecule distance of 5 Å was used and simulations were done at the temperature T = 1 K. We assumed a renormalization of the single particle energies δi = δ =1.83 eV (cf. Equation 3), an image-charge renormalization δic = 0.32 eV and a dielectric constant εr = 2.2 in order to fit our spectrum to the experiment of Swart et al. [14], which was carried out with CuPc on a NaCl(3 ML)/Cu(100) substrate with a workfunction of [Graphic 99] = 4.65 eV. With this, we find a triplet–singlet separation of the anionic ground state of 18 meV, which is in good agreement with experimental measurements of 21 meV [8]. Numerical results for the current and the differential conductance, according to Equation 37 and using the full Hamiltonian [Graphic 100] in Equation 29, are shown in Figure 4. Anionic (cationic) resonances at positive (negative) bias voltages are clearly seen.


Figure 4: Current and differential conductance curves exhibiting the anionic (cationic) resonance at positive (negative) bias voltage. Note that in contrast to all other results in this work, this curve is taken at a temperature of 60 K to emphasize the resonances in the dI/dV curve.

Notice that, in our model, the bias voltage at which a tip-mediated transition from the m-th neutral state to the n-th anionic state of the molecule is happening is


where e is the electron charge and αT accounts for the fact that in STM setups the bias voltage drops asymmetrically across the junction. Electrostatic considerations yielded αT = 0.59 for the tip and αS = −0.16 for the substrate [11]. If given without indices, Vres denotes the bias voltage corresponding to the groundstate-to-groundstate resonance.

The negative differential conductance at large negative bias in Figure 4 is caused by blocking due to population of excited states of the molecule. This has already been discussed in a previous work [27] and will not be of further interest here.

Transport simulations at finite magnetic fields

In Figure 5 we show the splitting of the anionic resonance with applied magnetic field in a dI/dV map. In the upper panel SOI is switched off, whereas in the lower panel it is switched on. One striking difference at first glance is the zero-field splitting for non-vanishing SOI, which is proportional to λ1 but enhanced by the bias drop, cf. Equation 38. For vanishing SOI, when Sz is a good quantum number, we can readily identify the corresponding transitions by using the effective spin Hamiltonian introduced in Equation 28. In the following, transitions from the neutral groundstate will be denoted by arabic numbers:

[Graphic 101]

while transitions from the neutral excited state will be denoted by Roman numerals:

[Graphic 102]

Other transitions are forbidden due to the selection rule for Sz, ΔSz = ±(1/2). The reason for the splitting into four lines observed in the upper panel of Figure 5 is that the orbital moment of the LUMO is not of the same size as the Bohr magneton.


Figure 5: Differential conductance maps as a function of the strength Bz of the magnetic field in z-direction. Upper (lower) panel: Spin–orbit interaction switched off (on). Solid and dashed lines depict the addition spectrum as calculated from the effective spin Hamiltonian, cf. Equation 28. Transitions starting from the neutral groundstate are denoted by solid lines, those from the neutral excited state by dashed lines.

For non-vanishing SOI, see lower panel of Figure 5, the definite assignment of transitions is not straightforward, at least for small magnetic fields. Since [Graphic 103] and [Graphic 104] are shifted downward by SOI, transition (2) now is the lowest lying transition, whereas transition (1) is shifted upward due to the positive contribution +λ1 to [Graphic 105]. Furthermore, transition (iv) is the only excited-state transition which can be definitely assigned to a line in the lower panel in Figure 5.

Figure 6 finally shows dI/dV maps as a function of the angle θ between the magnetic field and the z-axis. Panels (a), (b) and (c) show results obtained with vanishing SOI and panels (d), (e) and (f) are for finite SOI. Again, the results were fitted using the effective spin Hamiltonian introduced in Equation 28 with good agreement. The respective transitions can be identified by checking the assigned transitions in Figure 5 at the corresponding field strength.


Figure 6: Differential conductance maps vs the angle θ, formed by the applied magnetic field with the z-axis. Left (right) panels are without (with) SOI. Upper, middle and lower panels are calculated for a magnetic field strength of 1, 3 and 8 T, respectively. Solid and dashed lines depict the addition spectrum as calculated from the effective spin Hamiltonian, cf. Equation 28. Transitions starting from the neutral groundstate are denoted by solid lines, those from the neutral excited state by dashed lines.

Already at |B| = B = 1 T, cf. panels (a) and (d), the influence of SOI can be clearly seen. While for vanishing SOI any anisotropy of the dI/dV map is hidden beneath the temperature broadening, for finite SOI a slight θ-dependence can be observed. For B = 3 T, now also in the vanishing SOI case, Figure 6b, a slight anisotropy due to the orbital moment of the LUMOs can be observed, although still blurred by temperature. Again, at finite SOI in Figure 6e there is a much more pronounced dependence on θ. The high conductance areas at θ = 0° and θ = 180° for VbVres ≈ 0.8 meV correspond to the high conductance area in the middle of Figure 5 bottom, where many transitions are taking place at the same time. At B = 8 T, the magnetic field is dominating and a characteristic double cosine-like behaviour of the resonances can be observed, for both the case with no SOI, Figure 6c, and finite SOI, Figure 6f. For vanishing SOI, this behaviour is caused by the orbital moment of the LUMOs, since they interchange their positions when going from Bz to −Bz. The overall splitting between the double cosines, most evident at θ = 90°, is caused by the Zeeman term. The results for B = 8 T in Figure 6f at finite SOI are similar to those in Figure 6c, with the only difference that the cosine at large biases is more stretched, the one at low bias more compressed.


We established a model Hamiltonian for CuPc which accounts for electron–electron, spin–orbit and magnetic interactions in a minimal single particle basis represented by four frontier orbitals; the SOMO, the HOMO and two degenerate LUMOs. The distinct properties of these orbitals under rotations allowed us to deduce selection rules for matrix elements of the Coulomb interaction, which drastically reduce the number of nonvanishing terms and simplify the numerical diagonalization of the full many-body Hamiltonian. For the low-energy parts of the neutral and anionic blocks of the many-body spectrum we could further derive an effective spin Hamiltonian, capturing both SOI-induced splittings and magnetic anisotropy. Analogous Hamiltonians accounting for the effect of atomic SOI in molecular systems with orbital degeneracies have been derived for example in carbon nanotubes [31].

In order to study fingerprints of the SOI under realistic experimental conditions, we have studied the magnetotransport characteristics of a CuPc based junction in an STM setup. To this extent, a generalized master equation for the reduced density matrix associated to the full many-body Hamiltonian had to be solved in order to numerically obtain both the current and the differential conductance. Noticeably, by using the effective spin Hamiltonian, it was possible to reconstruct the nature of the many-body resonances observed in the numerical calculations.

In summary, we believe that our work significantly advances the present understanding of spin properties of CuPc. Moreover, the flexibility of our model Hamiltonian approach opens new perspectives for the investigation of other configurationally similar metallorganic compounds.

Supporting Information

Supporting Information File 1: Transformation from the atomic to the molecular orbital basis.
Format: PDF Size: 145.4 KB Download
Supporting Information File 2: Symmetries in the frontier orbitals basis.
Format: PDF Size: 185.8 KB Download
Supporting Information File 3: Details on the perturbative treatment of the SOI.
Format: PDF Size: 173.8 KB Download


The authors thank Thomas Niehaus, Jascha Repp and Dmitry Ryndyk for fruitful discussions. Financial support by the Deutsche Forschungsgemeinschaft within the research program SFB 689 is acknowledged.


  1. Gatteschi, D.; Sessoli, R.; Villain, J. Molecular Nanomagnets; Oxford University Press, 2006.
    Return to citation in text: [1] [2]
  2. Abragam, A.; Bleaney, B. Electron Paramagnetic Resonance of Transition Ions; Oxford University Press, 1970.
    Return to citation in text: [1]
  3. Mannini, M.; Pineider, F.; Danieli, C.; Totti, F.; Sorace, L.; Sainctavit, P.; Arrio, M.-A.; Otero, E.; Joly, L.; Cezar, J. C.; Cornia, A.; Sessoli, R. Nature 2010, 468, 417. doi:10.1038/nature09478
    Return to citation in text: [1] [2]
  4. Chiesa, A.; Carretta, S.; Santini, P.; Amoretti, G.; Pavarini, E. Phys. Rev. Lett. 2013, 110, 157204. doi:10.1103/PhysRevLett.110.157204
    Return to citation in text: [1] [2] [3]
  5. Misiorny, M.; Burzurí, E.; Gaudenzi, R.; Park, K.; Leijnse, M.; Wegewijs, M. R.; Paaske, J.; Cornia, A.; van der Zant, H. S. J. Phys. Rev. B 2015, 91, 035442. doi:10.1103/PhysRevB.91.035442
    Return to citation in text: [1] [2]
  6. Burzurí, E.; Gaudenzi, R.; van der Zant, H. S. J. J. Phys.: Condens. Matter 2015, 27, 113202. doi:10.1088/0953-8984/27/11/113202
    Return to citation in text: [1] [2] [3]
  7. Liao, M.-S.; Scheiner, S. J. Chem. Phys. 2001, 114, 9780. doi:10.1063/1.1367374
    Return to citation in text: [1] [2] [3]
  8. Mugarza, A.; Robles, R.; Krull, C.; Korytár, R.; Lorente, N.; Gambardella, P. Phys. Rev. B 2012, 85, 155437. doi:10.1103/PhysRevB.85.155437
    Return to citation in text: [1] [2] [3]
  9. Stepanow, S.; Mugarza, A.; Ceballos, G.; Moras, P.; Cezar, J. C.; Carbone, C.; Gambardella, P. Phys. Rev. B 2010, 82, 014405. doi:10.1103/PhysRevB.82.014405
    Return to citation in text: [1]
  10. Schulz, F.; Ijäs, M.; Drost, R.; Hämäläinen, S. K.; Harju, A.; Seitsonen, A. P.; Liljeroth, P. Nat. Phys. 2015, 11, 229. doi:10.1038/nphys3212
    Return to citation in text: [1]
  11. Siegert, B.; Donarini, A.; Grifoni, M. arXiv 2015, No. 1507.05504.
    Return to citation in text: [1] [2] [3] [4] [5]
  12. Ryndyk, D. A.; Donarini, A.; Grifoni, M.; Richter, K. Phys. Rev. B 2013, 88, 085404. doi:10.1103/PhysRevB.88.085404
    Return to citation in text: [1] [2]
  13. Repp, J.; Meyer, G.; Stojković, S. M.; Gourdon, A.; Joachim, C. Phys. Rev. Lett. 2005, 94, 026803. doi:10.1103/PhysRevLett.94.026803
    Return to citation in text: [1]
  14. Swart, I.; Sonnleitner, T.; Repp, J. Nano Lett. 2011, 11, 1580. doi:10.1021/nl104452x
    Return to citation in text: [1] [2] [3]
  15. Mann, J. B. Atomic Structure Calculations I. Hartree-Fock Energy Results for the Elements Hydrogen to Lawrencium; Los Alamos Scientific Laboratory of the University of California, 1967.
    Return to citation in text: [1] [2]
  16. Slater, J. C.; Koster, G. F. Phys. Rev. 1954, 94, 1498. doi:10.1103/PhysRev.94.1498
    Return to citation in text: [1]
  17. Froyen, S.; Harrison, W. A. Phys. Rev. B 1979, 20, 2420–2422. doi:10.1103/PhysRevB.20.2420
    Return to citation in text: [1]
  18. Siegert, B.; Donarini, A.; Grifoni, M. Phys. Status Solidi B 2013, 250, 2444. doi:10.1002/pssb.201350002
    Return to citation in text: [1] [2]
  19. Galassi, M.; Davies, J.; Theiler, J.; Gough, B.; Jungman, G.; Alken, P.; Booth, M.; Rossi, F. GNU Scientific Library Reference Manual, 3rd ed.; Network Theory Limited, 2009.
    Return to citation in text: [1]
  20. Bendix, J.; Brorson, M.; Schaffer, C. E. Inorg. Chem. 1993, 32, 2838–2849. doi:10.1021/ic00065a010
    Return to citation in text: [1]
  21. Yu, Z. G. Phys. Rev. B 2012, 85, 115201. doi:10.1103/PhysRevB.85.115201
    Return to citation in text: [1]
  22. Dzyaloshinsky, I. J. Phys. Chem. Solids 1958, 4, 241. doi:10.1016/0022-3697(58)90076-3
    Return to citation in text: [1]
  23. Moriya, T. Phys. Rev. Lett. 1960, 4, 228. doi:10.1103/PhysRevLett.4.228
    Return to citation in text: [1]
  24. Miyahara, S.; Fouet, J.-B.; Manmana, S. R.; Noack, R. M.; Mayaffre, H.; Sheikin, I.; Berthier, C.; Mila, F. Phys. Rev. B 2007, 75, 184402. doi:10.1103/PhysRevB.75.184402
    Return to citation in text: [1]
  25. Herzog, S.; Wegewijs, M. R. Nanotechnology 2010, 21, 274010. doi:10.1088/0957-4484/21/27/274010
    Return to citation in text: [1]
  26. Sobczyk, S.; Donarini, A.; Grifoni, M. Phys. Rev. B 2012, 85, 205408. doi:10.1103/PhysRevB.85.205408
    Return to citation in text: [1] [2] [3]
  27. Donarini, A.; Siegert, B.; Sobczyk, S.; Grifoni, M. Phys. Rev. B 2012, 86, 155451. doi:10.1103/PhysRevB.86.155451
    Return to citation in text: [1] [2]
  28. Kaasbjerg, K.; Flensberg, K. Phys. Rev. B 2011, 84, 115457. doi:10.1103/PhysRevB.84.115457
    Return to citation in text: [1]
  29. Darau, D.; Begemann, G.; Donarini, A.; Grifoni, M. Phys. Rev. B 2009, 79, 235404. doi:10.1103/PhysRevB.79.235404
    Return to citation in text: [1]
  30. Koch, J.; von Oppen, F. Phys. Rev. Lett. 2005, 94, 206804. doi:10.1103/PhysRevLett.94.206804
    Return to citation in text: [1]
  31. Laird, E. A.; Kuemmeth, F.; Steele, G. A.; Grove-Rasmussen, K.; Nygård, J.; Flensberg, K.; Kouwenhoven, L. P. Rev. Mod. Phys. 2015, 87, 703. doi:10.1103/RevModPhys.87.703
    Return to citation in text: [1]
Other Beilstein-Institut Open Science Activities