Pages

Spin Densities with NWChem

Electron spin density around an isolated gadolinium atom in an electric field.

Lately, I've been interested in the chiral-induced spin selectivity (CISS) effect as a means of achieving enantioseparation and controlling chemistry with spin-polarized surfaces. Much of the pioneering work on this topic has been performed only in the past couple of years mainly by the groups of Ron Naaman and Yossi Paltiel at the Weizmann Institute in Israel.

The CISS effect relies partially on the mirror symmetry-breaking of chiral molecules. Enantiomers are typically defined as non-superposable mirror images of the same compound, however enantiomers are not precise mirror images - small structural deviations between the enantiomers cause imperfect mirroring. This, plus the quantum dynamical current of electrons that flow through the molecule as it becomes polarized, can force enantiomers with nearly identical electric dipole moments to have simultaneous opposite magnetic dipole interactions with ferromagnetic materials by which they can be distinguished. This description captures the spirit of the CISS effect, but is not strictly correct since the "magnetic dipole interaction" is actually due to a spin-exchange interaction between the chiral molecule and the localized spins in the ferromagnetic material.

For a discrimination via the exchange interaction, we should expect that two enantiomers will have strong differences in the spatial arrangement of their static or induced spin density - a net excess of spin-up or spin-down electrons. Modeling the dynamical case is difficult, but we can perform some simple calculations in NWChem to illustrate some of these concepts in ground state molecules. (See the ARTAIOS code for more advanced spin transport simulations).

NWChem has a convenient module which exports the spin density around the nuclei for us as .cube files which can be rendered in Avogadro or VESTA. Ideally, we would do a fully-relativistic calculation including spin-orbit coupling, however it seems that the data exporting tool included with NWChem is not yet compatible with the SODFT module so we are forced to settle with performing open-shell DFT with relativistic ZORA corrections. Here is an example plot of the spin density calculated around a gadolinium atom.

Spin density calculated for a Gd atom in NWChem. Blue regions correspond to negative spin density (spin down) and yellow regions correspond to positive spin density (spin up).

Below, the calculated spin densities for the enantiomers of CHBrClF at the B3LYP x2c-tzvppall level can be seen to have very small differences which illustrates the gist of the CISS effect. The electric dipole moments are seen to be exactly mirrored as we expect, but the different spin densities will be expected to cause different adsorption rates near a spin-polarized surface.

start atom

TITLE "R_Chirality"
CHARGE 0
GEOMETRY nucleus point
C       -0.44902151    -0.58263682    -0.42107461
Cl       0.70451160    -1.83047885     0.07158718
H       -0.56238473    -0.61588205    -1.49796314
Br       0.17221920     1.21315323     0.03592669
F       -1.63875135    -0.80322165     0.17213501               
END

BASIS spherical
   * library x2c-tzvppall
END

RELATIVISTIC
   zora on
   zora:cutoff 1d-30
END

DFT
   odft
   xc b3lyp
   convergence energy 1e-8
   grid fine
   iterations 1000
   print "final vectors analysis"
   vectors output R-CBrClFH.movecs
END
task dft optimize
task dft energy

PROPERTY
   vectors R-CBrClFH.movecs
   dipole
   quadrupole
   octupole
   electrondensity
   #efield
   grid pad 10.0 ngrid 100 output R-CBrClFH-elf-pad.cube
END
task dft property

Mirror images of the chiral molecule CHBrClF. Note the small differences in the distribution of spin density, particularly the extended lobe of negative spin density in the upper-left of the left-handed molecule.

Real molecules of biological interest can also be modeled on a laptop if they are small enough. As we know, life exclusively uses left-handed amino acids. Below are plots of the spin densities of the amino acid serine. The geometry optimization took a few hours on my laptop.

start atom

TITLE "L_Chirality"
CHARGE 0
GEOMETRY nucleus point
C       -0.24953497    -1.24818253     0.08197602
C        0.10675344     0.19922910     0.46348136
C       -0.69823886     1.13599301    -0.43132052
N        1.54334480     0.43702801     0.28568460
O        0.72750105    -1.89409115    -0.56459874
O       -1.31849955    -1.73983910     0.32047127
O       -0.37381591     2.47071840    -0.05036869
H       -0.22634707     0.35051766     1.49340254
H       -1.75961374     0.92216681    -0.29515100
H       -0.43371621     0.95929202    -1.48050680
H        1.70448825     1.40765616     0.04475728
H        2.05247547     0.25608986     1.14154027
H        1.47249430    -1.25597205    -0.61495504
H       -0.84595434     3.08360191    -0.61953790                                       
END

BASIS spherical
   * library x2c-tzvppall
END

RELATIVISTIC
   zora on
   zora:cutoff 1d-30
END

DFT
   odft
   xc b3lyp
   convergence energy 1e-8
   grid fine
   iterations 1000
   print "final vectors analysis"
   vectors output L-serine.movecs
END
task dft optimize
task dft energy

PROPERTY
   vectors L-serine.movecs
   dipole
   quadrupole
   octupole
   electrondensity
   #efield
   grid pad 10.0 ngrid 100 output L-serine-elf-pad.cube
END
task dft property

Mirror images of the amino acid serine. Note the large difference in spin density between the enantiomers around the carboxylate moieties.

We can also add in point charges to simulate the electric field gradient a particle will experience as it approaches a conductive metallic surface. The snapshots provided might give us some idea of the magnitude of the dynamical CISS effect via the amount of rearrangement of the spin density. A series of snapshots is shown for an isolated Gd atom.

start atom

TITLE "Gd"
CHARGE 0
GEOMETRY nucleus point
   Gd 0.0 0.0 0.0              
END

BASIS spherical
   * library x2c-tzvpall
END

RELATIVISTIC
   zora on
   zora:cutoff 1d-30
END

BQ  
  0.000  -20.000   0.000   -35.0  
  0.000   20.000   0.000    35.0   
END

DFT
   odft
   mult 9
   xc b3lyp
   convergence energy 1e-8
   grid fine
   iterations 1000
   print "final vectors analysis"
   vectors output Gd.movecs
END
task dft energy

PROPERTY
   vectors Gd.movecs
   dipole
   quadrupole
   octupole
   electrondensity
   #efield
   grid pad 10.0 ngrid 100 output Gd-elf-pad.cube
END
task dft property

Animation of the influence of an increasing electric field strength across an isolated gadolinium atom. The electron density drawn to the positive charge has a net "down" spin orientation. The CISS effect is believed to work in a similar way - a chiral molecule approaching a conductive surface becomes polarized and the electron density has different spin polarization depending on the particular enantiomer.

Further Reading

Naaman, R.; Paltiel, Y.; Waldeck, D. H. Chiral Molecules and the Electron Spin. Nat Rev Chem 20193 (4), 250–260. https://doi.org/10.1038/s41570-019-0087-1.