Spin molecular-orbit coupling and magnetic properties of the decorated honeycomb layers of Mo₃S₇(dmit)₃ crystals

Strong spin orbit coupling (SOC) in weakly interacting systems can lead to unconventional insulating states such as topological insulators. Strongly correlated phases such as the topological Mott insulator can emerge from the interplay of strong Coulomb repulsion and SOC which may be realized in Ir-based transition metal oxides.¹ For instance, in Na₂IrO₃ and Li₂IrO₃ materials, SOC removes the orbital degeneracy of 5d electrons leading to S = 1/2 pseudospins which interact through anisotropic and quantum compass exchange interactions on a honeycomb lattice. These materials are potential realizations of the Heisenberg-Kitaev model, closely related to the Kitaev model which sustains a spin liquid state.² 

Multinuclear coordinated organometallic complexes are strongly correlated systems in which spin-orbit coupling (SOC) can be relevant. The multinuclear complex, Mo₃S₇(dmit)₃, consists of honeycomb networks of Mo₃S₇(dmit)₃ molecules stacked on top of each other along the c-direction of the crystal. Since the Mo₃S₇(dmit)₃ molecules can be described by three Wannier orbitals,³ their packing on the honeycomb lattices within the layers lead to decorated honeycomb lattices, as shown in Fig. 1. Based on perturbative expansions including Coulomb repulsion and SMOC, the effective spin exchange model which describes the magnetic properties of these layers is a S = 1 XXZ+120⁰ quantum compass model on the honeycomb lattice.^4,5

In the present paper, we discuss the role played by the single-spin anisotropy induced by SMOC in the magnetic properties of the honeycomb layers of Mo₃S₇(dmit)₃ crystals shown in Fig. 1. We predict that isolated honeycomb layers of Mo₃S₇(dmit)₃ crystals are Néel ordered but increasing the magnetic frustration of the lattice can drive the system into a large-D phase.

The isolated triangular clusters can be described by a Hubbard-Heisenberg model in the presence of the spin molecular-orbit coupling (SMOC):

where:

with t_c the hopping between hybrid metal-ligand orbitals at sites in the triangular clusters and $aiσ†$ creates an electron at the ith Wannier orbital with spin σ.

An important question is how SOC affects the low energy states described by these Wannier orbitals. Generically the one-electron spin-orbit coupling takes the form H_SOC = ∑_αK_α ⋅s_α, where s_α is the spin of the α^th electron and K_α is a pseudovectorial operator that acts only on the spatial part of the wavefunction. For example, in the Pauli approximation,

where p is the single electron momentum operator and ϕ(r) is the potential through which the electron moves.

Spherical symmetry implies that for an isolated atom K ∝ L, the angular momentum operator. Thus, in ab initio calculations, it is common to write the SOC as a linear superposition of atomic contributions, for example, in the mean-field Breit-Pauli approximation^7,8

where $ZAeff$ is the effective charge of nucleus A which accounts for the screening effects of the rest of the electrons on the nuclear potential, l_αA = (r_α − R_A) × p_α, is the angular momentum of the α^th electron relative to the nucleus A at R_A. This atomic approach does not provide any insight into how SOC acts on the low-energy states described by the Wannier orbitals and requires detailed first principles calculations for each new material studied. A more chemically intuitive description should consider the coupling between the electronic spin and the molecular orbital degrees of freedom.

Symmetry dictates the general form of the SMOC.⁶ For cyclic molecules with C_N symmetry the coupling is⁶ 

where $cmσ†=i|m|∑jajσ†ei2πj/N/N$⁠. For a C₃ symmetric molecule, such as Mo₃S₇(dmit)₃, this simplifies to

where $Lαx,y,z$ are the Cartesian components of the one-electron molecular orbital “angular momentum” operator for the α^th electron, λ_xy, is the transverse SMOC and λ_z the longitudinal SMOC.

The SMOC hamiltonian (6) gives a far clearer understanding of the relevant low-energy physics than the atomistic description (4). This is because molecular orbitals, not atomic orbitals, are the natural basis for the discussion of molecular physics. As molecular orbitals for cyclic molecules are Bloch states of the Wanniers, H_SMOC describes the coupling of orbital currents running around the plane of the molecule to the electron’s spin. The low-energy model constructed from Wannier spinors from four component relativistic band structure calculations⁹ of Mo₃S₇(dmit)₃ yields single molecule SMOC exactly as predicted by Eq. (6).

Finally, the Hubbard-Heisenberg contribution describing the Coulomb and exchange energy in each triangular cluster is:

where U is the onsite Hubbard interaction, J_F is an intracluster exchange interaction, and $niσ=aiσ†aiσ$ the number operator. The direct ferromagnetic exchange, J_F < 0, plays a crucial role in generating magnetic anisotropies.⁴ 

The nearest-neighbor triangular clusters are connected through the hopping amplitude, t:

where ⟨lm⟩ denotes two nearest-neighbor triangular clusters and i a cluster site in a given molecule.

By straightforward diagonalization of the full hamiltonian (1), we have found that in the strong coupling limit, U ≫ t, isolated triangular clusters with four electrons effectively behave as pseudospin-1 localized moments.⁴ The many-body states of the cluster can be classified according to the z-component of the total angular momentum, J_z = L_z + S_z. The ground state of the isolated cluster with no SMOC is a triplet. As shown in Fig.2, SMOC splits the lowest energy triplet into a non-degenerate singlet (j = 0) which is the ground state and a doublet (j = ±1), where j denotes the eigenstates of J_z. Since we have an even number of electrons in the cluster, Kramers theorem does not apply and non-degenerate states are possible. SMOC can be effectively accounted for as a single-spin anisotropy contribution at each cluster. Hence, the effective spin model for the m-th Mo₃S₇(dmit)₃ molecule in the crystal is just:

where $Srm$ describes the effective pseudospin-1 localized at a triangular cluster. The single-spin anisotropy, D, depends strongly on SMOC as shown in Fig. 3. Using perturbation theory^4,5 we have obtained an effective spin exchange model describing the magnetic exchange coupling between the pseudospin-1 in the a − b planes of Mo₃S₇(dmit)₃. This model is an XXZ+120⁰ honeycomb quantum compass model⁴ leading to:

where J = (J_xx + J_yy)/2, Δ = J_zz/J and Q = (J_xx − J_yy)/2 and ϕ_i = 2π(i − 1)/3 where i labels the three bonds around each ▽ shown in Fig. 1. Thus, we see that the second term (proportional to J) is simply the XXZ model and the third term (proportional to Q) is the honeycomb 120° compass model.¹⁰ 

The exchange couplings and single-spin anisotropy entering model (10) have been obtained numerically from perturbation theory⁴ as well as analytically using a canonical transformation.⁵ In Fig. 3 we show the dependence of the model (10) parameters: J, Q, Δ and D with λ_xy for the longitudinal SMOC, λ_z = λ_xy/2, appropriate for Mo₃S₇(dmit)₃ crystals. With no SMOC, the exchange couplings are isotropic: J ≈ 0.04t_c, Q = 0, Δ = 1, J_xz = 0, as it should. As SMOC increases, spin exchange anisotropy, i. e., Q ≠ 0, Δ ≠ 1 (or equivalently J_xx ≠ J_yy ≠ J_zz) increases. On the other hand, Fig. 3 shows how the single-spin anisotropy is rapidly enhanced by SMOC.

One can see from the results of Fig. 3 that, for sufficiently large SMOC, the single-spin anisotropy overcomes the exchange couplings, D > J, ΔJ = J_zz. We can expect that in the limit, D ≫ J, J_zz, the model is dominated by single-spin anisotropy and the ground state of the model is just the tensor product of $Srℓz=0$ at each lattice site. Hence, it is important to first obtain an estimate of the critical D_c at which a transition to a large-D phase occurs and the magnitude of SMOC needed to reach such critical D. In order to do so, we explore below the magnetic properties of a slightly simplified version of the full model (10) in the limit of weak SMOC. This limit is relevant to the a − b planes of Mo₃S₇(dmit)₃, for which ab initio calculations⁹ give $λxyMo3S7(dmit)3=0.042tc$ with λ_z = λ_xy/2. We can see from Fig. 3 that for λ_xy ≲ 0.5 exchange coupling anisotropies are negligible which justifies using an isotropic Heisenberg model for describing the magnetic properties of the a − b planes of the crystal.

We have recently explored the Néel to large-D transition of the isotropic version of model (10) using a SU(3) description of the spins which introduces three Schwinger bosons to properly account for the three projections: $Siz=0,±1$ of the spin-1 at each lattice site.¹¹ We find a transition from the Néel ordered phase to a large-D phase for D_c ≈ 3J. From this condition and the dependence of D with λ_xy we can obtain a realistic estimate of SMOC needed for the large-D phase to occur in isolated a − b layers of Mo₃S₇(dmit)₃. In Fig. 3 we show the critical SMOC $(λxy/tc)critic≈0.38$ at which the D_c = 3J condition is satisfied. This critical value is about nine times larger than SMOC in Mo₃S₇(dmit)₃ crystals. Hence, we conclude that isolated a − b layers of Mo₃S₇(dmit)₃ should be in a magnetically ordered Néel state.

Magnetic frustration plays an important role in inducing disordered spin liquid states. We analyze the effect of a next-nearest neighbor exchange coupling, J′, on the magnetic properties of model (10) in the limit of weak SMOC i. e. with isotropic couplings. Hence, the model we study is a S = 1 J − J′ Heisenberg on the honeycomb lattice with single-spin anisotropy:

We first discuss the phase diagram of this model with no single-spin anisotropy, D = 0, based on recent SU(2) Schwinger boson mean-field theory (SBMFT) calculations.¹¹ For J′ = 0, the ground state of the model is a Néel state, as expected since the lattice is bipartite. A direct transition from Néel order to spiral order occurs for J′ ≈ 0.24.

We finally discuss the effect of J′ on the Néel to large-D transition of model (11) based on SU(3) SBMFT calculations of model (11). In Fig. 4 we show the critical D at which a spin gap opens up signaling the formation of the paramagnetic large-D phase. We find that the critical D_c = 3J obtained for J′ = 0 of the unfrustrated model, is strongly suppressed by J′. Hence, magnetic frustration favors the large-D phase. The kink around J′ = 0.12J in Fig. 4 is related to the expected poor performance of the SU(3) SBMFT approach as D_c is suppressed below J, D_c < J.¹¹ 

The magnetic properties of the a − b planes of Mo₃S₇(dmit)₃ crystals can be modeled through a S = 1, XXZ+120⁰ quantum compass model with single-spin anisotropy.⁴ For the weak SMOC present in Mo₃S₇(dmit)₃ crystals,⁹ quantum compass couplings and single-spin anisotropy play a minor role so that the effective model in this regime is just the isotropic S = 1 Heisenberg model on the honeycomb lattice. Hence, we predict that isolated honeycomb layers of Mo₃S₇(dmit)₃ are Néel ordered. However, we find that magnetic frustration induced by a next-nearest neighbor exchange coupling, J′, favors the large-D phase.

J. M. acknowledges financial support from (MAT2015-66128-R) MINECO/FEDER, UE, and A.R. from the French National Research Agency (Contract ORGANISO, ANR-15-CE09- 0017). Work at the University of Queensland was supported by the Australian Research Council (DP160100060).