Computation of the gauge-including magnetically induced current

Historically, computations of magnetic properties led to a problem that arises from the choice of the gauge origin. In a homogeneous medium, a magnetic field B can be defined through a vector potential A(r) (B = ∇ × A(r)) but the opposite does not hold due to the arbitrary scalar function found in B = ∇ × (A(r) + ∇ Φ(r))¹⁰. In standard DFT computations the use of finite basis sets introduces a gauge dependence in the magnetic calculations that can be overcome only in the limit of a complete basis set. To avoid this problem, the idea to use magnetic field-dependent basis functions goes back to the calculation of ring currents by the seminal work of London²⁴. In this way the gauge dependence on the origin can be eliminated. These functions are known as gauge—including atomic orbitals (GIAO) defined in the Coulomb gauge (∇ ⋅ A(r) = 0). In large-scale computations small basis sets can be used to achieve rapid convergence^5,10.

The shielding tensor σ_C,λη for the nucleus C is defined as the second derivative of the energy with respect to the external magnetic field B and the nuclear magnetic moment μ_C in the limit of zero magnetic field and zero nuclear magnetic moment (λ, η = x, y, z)²⁵.

An alternative expression to compute the nuclear magnetic shielding tensor via the Biot–Savart law^1,26–28 is

Atomic and electronic structure of the studied superatoms

We studied two previously reported phosphine-protected gold-based clusters, [HAu₉(PPh₃)₈]²⁺ (1, ref. ^21,30) and [HPtAu₈(PPh₃)₈]⁺ (2, ref. ^22,23). The position of the hydrogen atom in the metal core has been previously deduced by NMR analysis in clusters 1²¹ and 2^22,23. Their DFT-optimized structures are shown in Fig. 1. Cluster 1 shows a C_4v symmetry in the core (formally (HAu₉)²⁺) and a C₄ symmetry including the ligand layer. Cluster 2 shows a lower C₂ symmetry both in the metal core (formally (HPtAu₈)⁺) and in the ligand layer. The two “crown-like” metal cores have an atom (Au/Pt) at the center and can be characterized by four distances d1–d4 and their respective distortions d1_a,b–d3_a,b (Fig. 1).

Fig. 1

Structure of the studied clusters.

a [HAu₉(PPh₃)₈]²⁺ (cluster 1 in text) and (b) [HPtAu₈(PPh₃)₈]⁺ (cluster 2). Characteristic distances in metal cores are labeled (dN and dN_a,b, N = 1–4). Au: yellow, P: orange, C: grey, core-hydrogen: white. Some PPh₃ ligands and protons in the phenyl rings are omitted for clarity.

For the cluster 1 the distance d3 (2.77 Å) is longer than d4 (2.70 Å) both preserving a perfect square on top and bottom from the central Au atom (square edge being 3.57 and 3.40 Å respectively). In cluster 2 there are three characteristic distances from the central Pt atom, a shorter one to the top square (2.66 Å) and two longer (2.68 and 2.71 Å) to the bottom parallelogram. Both parallelograms are distorted from the perfect square showing longer distances on top (edges 3.65 and 3.58 Å) than on bottom (edges 3.10 and 3.15 Å). This indicates the effect of the central atom (Au/Pt) on the geometry. For cluster 1 the distance between H and the central Au is 1.72 Å and for cluster 2 the corresponding distance to central Pt atom is shorter, 1.63 Å.

The electronic structure of 1 and 2 was analyzed by projecting the frontier orbitals to cluster-centered spherical harmonics Y_lm. The projected density of states (Supplementary Note 1 and Supplementary Fig. 1) shows that both clusters are 8-electron superatoms with a closed-shell configuration and significant HOMO-LUMO energy gaps (1.88 eV for 1 and 1.93 eV for 2) as expected from their chemical formula and overall charge¹¹. The Bader atomic charge analysis³¹ showed for both hydrogens nearly similar, vanishing charges, −0.09 e in 1 and −0.06 e in 2. This, combined with analysis of the superatomic orbitals (Supplementary Fig. 1) indicates that both hydrogens can be thought of as neutral atoms, contributing their s-electron to the superatomic electron system. Frontier orbital analysis of 1 done in ref. ²¹ led to similar conclusions.

NMR chemical shifts

First, we calculated proton shifts for the ortho-, meta-, and para-positions in PPh₃ ligands for clusters 1 and 2 and compared with experimental data (Table 1 and Supplementary Table 1). This served also as an internal test for the DFT approximations (xc functional and basis set) in the calculations.

Table 1

Experimental and DFT ¹H NMR chemical shift δ and standard deviation σ in ppm of ligand protons in clusters 1 and 2 (PBE xc functional and SDD pseudopotential).

		PBE/SDD
1	Expa	δ	σ
o-proton	7.25	7.44	0.89
m-proton	6.62	7.27	0.39
p-proton	7.05	7.53	0.22
2	Expb	δ	σ
o-proton	7.29	7.68	1.07
m-proton	6.59	6.98	0.60
p-proton	7.02	7.24	0.22

In cluster 1 the measured signals²¹ fragmented into three groups. A broad signal at 7.25 ppm was assigned to the o-proton. This smeared out triplet peak was explained by the spin–spin coupling with the neighboring phosphorous and the m-proton. The 6.62 ppm triplet signal was assigned to the m-proton and the 7.05 ppm triplet was assigned to the p-proton. The DFT results showed a 7.26–7.44 ppm chemical shift for the o-proton depending on the level of theory employed. This proton showed the highest standard deviation (0.78–0.89 ppm) independently of the level of theory employed within the three o-, m-, and p-proton analyzed, in good agreement with the chemical shifts observed. Computed values for the m-proton were in the range of 7.16–7.27 ppm with standard deviation of 0.38–0.39 ppm, and for the p-proton in the range of 7.44–7.53 ppm with standard deviation of 0.22–0.26 ppm.

Since the ligand layer of cluster 2 is similar to that of cluster 1, the measured o-, m-, and p-proton shifts showed almost identical values to the ones measured for 1. The computed values were within 0.4 ppm from the experimental values for o- and m-protons and within 0.2 ppm for p-protons (Table 1 and Supplementary Table 1). We note that all our computed values for the proton shifts are based on a single structure of the cluster and a full exploration of the distribution of the proton shifts due to ligand dynamics and potential solvent effects is not possible at the moment.

Next, we turn to the discussion of the NMR shift of the hydrogen in the metal cores of 1 and 2. The hydrogen peak of 1 was observed at 15.1 ppm²¹ in comparison with the previously reported 5.4 ppm for 2^22,23 (Table 2 and Supplementary Table 2). This anomalously large difference (Δ = 9.7 ppm) is rather surprising in view of similar atomic geometries of 1 and 2. Our computed values are notably shifted upfield (Table 2 and Supplementary Table 2) but show consistently a qualitatively similar difference in the range of Δ = 4.3–4.5 ppm, independent of the level of DFT approximations. The systematic upfield shift and smaller Δ obtained in the computations could be explained by the limitations in the description of the electron density and the scalar-relativistic approximation used in the effective core potential^32–36. These approximations preclude the calculations of the shielding tensor on the true all-electron density around the Pt and Au atoms including also the full spin-orbit effects.

Table 2

Experimental and DFT ¹H NMR chemical shift δ in ppm of the core-hydrogen in clusters 1 and 2 (PBE xc functional and SDD pseudopotential).

	Exp	PBE/SDD
1	15.1a	3.88
2	5.4b	−0.55
Δ	9.7	4.43

The anomalous difference of the hydrogen shifts in 1 and 2 prompted the authors of ref. ²¹ to speculate the reasons arising either from different hydrogen charge or different magnetic shielding around the hydrogen in these superatoms. However, our charge analysis indicated essentially identical hydrogen charge in both systems as discussed above. Hence, we next seek the explanation by examining the computed magnetically induced currents in 1 and 2. We chose to use the PBE xc functional and the SDD pseudopotential in these computations.

Magnetically induced currents

The nine-component current density tensor $$ {\mathcal{J}}_{\tau }^{\left(\lambda \right)}$$ (Eq. (3)) is a complex object and difficult to visualize. But, given a particular direction of the external magnetic field, one can analyze the current density as a vector function J(r). We chose an external magnetic field B = 1 T in a direction pointing from the central Au/Pt atom of clusters 1 and 2 to the hydrogen site and analysed J(r) in a cube containing 20 × 20 × 20 grid points centered at the central atom.

3D visualizations of J(r) (Supplementary Fig. 2 and Supplementary Movies 1 and 2) reveal a surprising qualitative indirect effect of the central Pt atom in cluster 2 on MICs around the core-hydrogen site. It is seen that Pt-doping of the gold changes significantly the diatropic current around the core-hydrogen as compared to cluster 1. To understand this interesting effect, a more quantitative analysis of this phenomenon is presented in Figures 2 and 3 and in Table 3. Figure 2 shows the projected vectorial field of J(r) as well as the diatropic and paratropic contributions in the x-y planes containing either the core-hydrogen or the central Au/Pt atoms in 1 and 2. Pt-doping in 2 has a tiny paratropic current basically at the metal site but induces a large area of diatropic current in the metal core of 2. In cluster 1, the corresponding gold site induces a deep maximum in the diatropic current immediately around the metal site but the paratropic contribution comes into play in the vicinity, neutralizing the diatropic effect. This interplay is seen quantitatively in Fig. 3 that shows circularly integrated currents in the x-y planes containing the central metal site. Radial integrals of the contributions up to 4 Å radius show diatropic/paratropic components of (−2.619/0.715) ŲnA/T in 2 as compared to (−1.626/1.184) ŲnA/T in 1 (table 3). This causes the total integrated current to be −1.905 ŲnA/T in 2 vs. −0.442 ŲnA/T in 1 indicating an enhanced aromaticity in that plane in cluster 2.

Fig. 2

Magnetically induced currents.

Paratropic (red) and diatropic (blue) vectorial field of the magnetically induced current density J(r) for (a–c) [HAu₉(PPh₃)₈]²⁺ (cluster 1) and (d–f) [HPtAu₈(PPh₃)₈]⁺ (cluster 2). Vectors projected on an x–y plane located at the z–coordinate of the core-hydrogen and the central atom (Au/Pt) and signed modulus of J(r) are shown. The external magnetic field B is oriented perpendicular to the x–y plane and pointing towards the reader.

Fig. 3

Integrated currents.

Circularly integrated current profile of the signed modulus of J(r) for (a,b) [HAu₉(PPh₃)₈]²⁺ (cluster 1, top row) and (c,d) [HPtAu₈(PPh₃)₈]⁺ (cluster 2, bottom row). The integrals are done in the x-y planes shown in Fig. 2, having either core-hydrogen or Au/Pt at the origin. Diatropic (blue), paratropic (red), and total (black) contributions are indicated.

Table 3

Radial integrals of the circularly integrated current strength (Fig. 3) (in Ų nA/T) at H and Au or Pt x-y planes at radius r from (0,0,z) respective coordinate in clusters 1 and 2.

1	H z–coordinate			Au z–coordinate
r (Å)	Tot	Para.	Dia.	Tot	Para.	Dia.
1	−0.045	0.001	−0.046	−0.715	0.000	−0.715
2	0.137	0.190	−0.052	−0.644	0.250	−0.893
3	0.225	0.892	−0.667	−0.422	1.013	−1.435
4	0.218	1.253	−1.034	−0.442	1.184	−1.626
2	H z–coordinate			Pt z–coordinate
r (Å)	Tot	Para.	Dia.	Tot	Para.	Dia.
1	−0.11	0	−0.11	−0.362	0.013	−0.374
2	−0.191	0.108	−0.299	−1.173	0.109	−1.282
3	−0.295	0.51	−0.805	−1.341	0.686	−2.027
4	−0.63	0.667	−1.297	−1.905	0.715	−2.619

The enhanced aromaticity in the Pt-doped cluster 2 has an indirect effect also in the plane containing the core-hydrogen, revealing a key mechanism to understand the previously reported proton NMR shifts.^21–23 Analogously to the central x-y plane discussed above, the Pt-doped cluster has also a larger diatropic environment around the core-hydrogen as compared to the all-gold cluster. Radially integrated diatropic/paratropic currents up to 4 Å in the plane containing the core-hydrogen are (−1.297/0.667) ŲnA/T for 2 and (−1.034/1.253) ŲnA/T for 1. As a result, the total integrated current around the core-hydrogen in 1 is paratropic (0.218 ŲnA/T) while in 2 it is diatropic (−0.63 ŲnA/T). This qualitative and quantitative difference in the magnetic shielding between 1 and 2 provides the explanation for the reported²¹ downfield shift of the core-hydrogen in 1.

Ground state DFT calculations

The structures of clusters 1 and 2 were optimized using the real-space grid based GPAW^37,38 program with a uniform grid spacing of 0.2 Å. Pt, Au, C, and P atoms were described by valence of 16, 11, 4, and 5, respectively, and the frozen core approximation was used for inner electrons. The PAW setups for Pt and Au include scalar-relativistic effects. The PBE³⁹ functional was employed for the exchange-correlation interaction. Structure optimization was continued until the force on each atom was less than 0.05 eV/Å. To identify the superatom symmetries of the frontier orbitals, the Y_lm analysis¹¹ was performed for all the clusters with a cutoff radius of 4.0 Å. A Bader charge analysis³¹ was performed for the studied systems using the density obtained in a single point computation in GPAW.

NMR chemical shifts and magnetic currents

The magnetic shielding tensors and gauge-including magnetically induced currents were computed in deMon2k²⁹ code employing the Auxiliary Density Functional Theory (ADFT) approach. For the ¹H NMR chemical shift and the MIC computations the geometries optimized in GPAW were used with PBE³⁹ and BP86^40,41 level for exchange-correlation effects. The Stuttgart–Dresden (SDD) pseudopotential⁴² (Au and Pt are described by a valence of 19 and 18 electrons respectively) along with the DZVP⁴³ basis set (P, C, and H) and the effective core potentials LANL2DZ⁴⁴ along with D95⁴⁵ basis set (P, C, and H) were employed. All the computations were performed in combination with the GEN-A2* auxiliary function set. Each ¹H chemical shift was referenced to the ¹H TMS chemical shift computed in the respective level of theory.