## Publications

126 results found

Li M, Cucinotta C, Horsfield A, 2024, A computational model for a molecular chemical sensor, *Nanoscale*, Vol: 16, Pages: 5334-5342, ISSN: 2040-3364

In this study, we propose that a molecular junction with a sharp Negative Differential Resistance (NDR) current peak could improve the selectivity, thereby functioning as a potential molecular sensor for molecule recognition. Using DFT–NEGF simulations, we investigate the connection between molecule–molecule coupling, molecule–electrode coupling and the corresponding NDR peak shape. Based on this analysis we propose three design rules to control the sensitivity of a sensor and determine that one mechanism for NDR is for a localised molecular orbital involved in resonant tunneling to enter and leave the bias window. Our findings provide useful insight into the development of single molecule sensors for molecule recognition.

Wells T, Foulkes W, Dudarev S,
et al., 2023, The Einstein-de Haas effect in an Fe₁₅ cluster, *Journal of Physics: Condensed Matter*, Vol: 35, Pages: 1-16, ISSN: 0953-8984

Classical models of spin-lattice coupling are at present unable to accurately reproduce results for numerous properties of ferromagnetic materials, such as heat transport coefficients or the sudden collapse of the magnetic moment in hcp-Fe under pressure. This failure has been attributed to the absence of a proper treatment of effects that are inherently quantum mechanical in nature, notably spin-orbit coupling. This paper introduces a time-dependent, non-collinear tight binding model, complete with spin-orbit coupling and vector Stoner exchange terms, that is capable of simulating the Einstein-de Haas effect in a ferromagnetic Fe15 cluster. The tight binding model is used to investigate the adiabaticity timescales that determine the response of the orbital and spin angular momenta to a rotating, externally applied Β field, and we show that the qualitative behaviours of our simulations can be extrapolated to realistic timescales by use of the adiabatic theorem. An analysis of the trends in the torque contributions with respect to the field strength demonstrates that SOC is necessary to observe a transfer of angular momentum from the electrons to the nuclei at experimentally realistic Β fields.The simulations presented in this paper demonstrate the Einstein-de Haas effect from first principles using a Fe cluster.

Mario Z, Horsfield A, Lischner J, 2023, Accelerating GW calculations through machine learned dielectric matrices, *npj Computational Materials*, Vol: 9, ISSN: 2057-3960

The GW approach produces highly accurate quasiparticle energies, but its application to large systems is computationally challenging due to the difficulty in computing the inverse dielectric matrix. To address this challenge, we develop a machine learning approach to efficiently predict density–density response functions (DDRF) in materials. An atomic decomposition of the DDRF is introduced, as well as the neighborhood density–matrix descriptor, both of which transform in the same way under rotations. The resulting DDRFs are then used to evaluate quasiparticle energies via the GW approach. To assess the accuracy of this method, we apply it to hydrogenated silicon clusters and find that it reliably reproduces HOMO–LUMO gaps and quasiparticle energy levels. The accuracy of the predictions deteriorates when the approach is applied to larger clusters than those in the training set. These advances pave the way for GW calculations of complex systems, such as disordered materials, liquids, interfaces, and nanoparticles.

Bustamante C, Gadea E, Todorov T,
et al., 2023, Fluorescence in quantum dynamics: accurate spectra require post-mean-field approaches, *Journal of Chemical Physics*, Vol: 158, Pages: 1-11, ISSN: 0021-9606

Real time modeling of fluorescence with vibronic resolution entails the representation of the light–matter interaction coupled to a quantum-mechanical description of the phonons and is therefore a challenging problem. In this work, taking advantage of the difference in timescales characterizing internal conversion and radiative relaxation—which allows us to decouple these two phenomena by sequentially modeling one after the other—we simulate the electron dynamics of fluorescence through a master equation derived from the Redfield formalism. Moreover, we explore the use of a recent semiclassical dissipative equation of motion [C. M. Bustamante et al., Phys. Rev. Lett. 126, 087401 (2021)], termed coherent electron electric-field dynamics (CEED), to describe the radiative stage. By comparing the results with those from the full quantum-electrodynamics treatment, we find that the semiclassical model does not reproduce the right amplitudes in the emission spectra when the radiative process involves the de-excitation to a manifold of closely lying states. We argue that this flaw is inherent to any mean-field approach and is the case with CEED. This effect is critical for the study of light–matter interaction, and this work is, to our knowledge, the first one to report this problem. We note that CEED reproduces the correct frequencies in agreement with quantum electrodynamics. This is a major asset of the semiclassical model, since the emission peak positions will be predicted correctly without any prior assumption about the nature of the molecular Hamiltonian. This is not so for the quantum electrodynamics approach, where access to the spectral information relies on knowledge of the Hamiltonian eigenvalues.

Li B, Xiao C, Harrison N,
et al., 2023, Role of electron localisation in H adsorption and hydride formation in the Mg basal plane under aqueous corrosion: a first-principles study, *Physical Chemistry Chemical Physics*, Vol: 25, Pages: 5989-6001, ISSN: 1463-9076

Understanding hydrogen-metal interactions is important in various fields of surface science, including the aqueous corrosion of metals. The interaction between atomic H and a Mg surface is a key process for the formation of sub-surface Mg hydride, which may play an important role in Mg aqueous corrosion. In the present work, we performed first-principles Density Functional Theory (DFT) calculations to study the mechanisms for hydrogen adsorption and crystalline Mg hydride formation under aqueous conditions. The Electron Localisation Function (ELF) is found to be a promising indicator for predicting stable H adsorption in the Mg surface. It is found that H adsorption and hydride layer formation is dominated by high ELF adsorption sites. Our calculations suggest that the on-surface adsorption of atomic H, OH radicals and atomic O could enhance the electron localisation at specific sites in the sub-surface region, thus forming effective H traps locally. This is predicted to result in the formation of a thermodynamically stable sub-surface hydride layer, which is a potential precursor of the crucial hydride corrosion product of magnesium.

Smutna J, Wenman MR, Horsfield AP,
et al., 2023, The bonding of H in Zr under strain, *Journal of Nuclear Materials*, Vol: 573, Pages: 1-12, ISSN: 0022-3115

Accurate computer simulation is important for understanding the role of irradiation-induced defects in zirconium alloys found in nuclear reactors. Of particular interest is the distribution and trapping of hydrogen, and the formation of zirconium hydride. These simulations require an accurate representation of Zr-H bonding in order to predict the behaviour of H around atomic-scale defects, dislocation lines, and dislocation loops. Here we explore the bonding of H in Zr under strain, how well it is represented by state-of-the-art Embedded Atom Method (EAM) potentials, and what physics is needed for an accurate representation in a Linear Combination of Atomic Orbitals (LCAO) Density Functional Theory (DFT) framework. For H in dilute solution under hydrostatic strain in the range -10% to +10%, solution energies and Zr-H bond lengths computed using EAM potentials are found to be in poor agreement with plane-wave DFT results. We note that the bond lengths are in a poor agreement even in equilibrium. LCAO basis sets are used to explore the importance of electron distribution around H atoms, and the transfer of electrons between H and Zr. The electron distribution around H atoms is found to be important to the explanation of the difference between octahedral and tetrahedral interstitial sites for H, with H in a tetrahedral site having very similar bonding to H in zirconium hydrides. The interatomic electron transfer has a smaller impact but is needed for maximum accuracy.

Li M, Cucinotta CS, Horsfield AP, 2022, The influence of surface Fe on the corrosion of Mg, *Journal of Physics and Chemistry of Solids*, Vol: 170, ISSN: 0022-3697

Iron is a common impurity in magnesium alloys, and is acknowledged to accelerate Mg corrosion, contributing to Mg’s poor corrosion resistance. However, an atomistic understanding of this acceleration effect is still incomplete. Here we use Density Functional Theory simulations performed with the Quantum Espresso package to investigate several Fe/Mg models, calculating the associated work functions, atomic charges, and H and Fe absorption energies. Compared with a pure Mg slab, we find that Fe’s existence increases the work function and decreases the H adsorption energy. Furthermore, a general trend is observed that the Fe absorption energy decreases with increasing interaction between Fe atoms on the Mg substrate. Based on these results, a mechanism based on charge redistribution is put forward to explain how Fe accelerates the corrosion of Mg. Our findings provide insight into Mg’s corrosion process at the atomic level that might inform future measures to prevent corrosion.

Fogarty R, Horsfield A, 2022, Molecular dynamics study of structure and reactions at the hydroxylated Mg(0001)/bulk water interface, *The Journal of Chemical Physics*, Vol: 157, ISSN: 0021-9606

A molecular level understanding of the aqueous Mg corrosion mechanism will be essential in developing improved alloys for battery electrodes, automobile parts, and biomedical implants. The structure and reactivity of the hydroxylated surface is expected to be key to the overall mechanism because (i) it is predicted to be the metastable surface state (rather than the bare surface) under a range of conditions and (ii) it provides a reasonable model for the outer corrosion film/water interface. We investigate the structure, interactions, and reactivity at the hydroxylated Mg(0001)/water interface using a combination of static Density Functional Theory calculations and second-generation Car–Parrinello ab initio molecular dynamics. We carry out detailed structural analyses into, among other properties, near-surface water orientations, favored adsorption sites, and near-surface hydrogen bonding behavior. Despite the short timescale (tens of ps) of our molecular dynamics run, we observe a cathodic water splitting event; the rapid timescale for this reaction is explained in terms of near-surface water structuring lowering the reaction barrier. Furthermore, we observe oxidation of an Mg surface atom to effectively generate a univalent Mg species (Mg+). Results are discussed in the context of understanding the Mg corrosion mechanism: For example, our results provide an explanation for the catalytic nature of the Mg corrosion film toward water splitting and a feasible mechanism for the generation of the univalent Mg species often proposed as a key intermediate.

Fogarty R, Li B, Harrison N,
et al., 2022, Structure and interactions at the Mg(0001)/water interface: An ab initio study, *Journal of Chemical Physics*, Vol: 156, ISSN: 0021-9606

A molecular level understanding of metal/bulk water interface structure is key for a wide range of processes including aqueous corrosion, our focus, but their buried nature makes experimental investigation difficult and means we must mainly rely on simulations. We investigate the Mg(0001)/water interface using second generation Car-Parrinello molecular dynamics (MD) to gain structural information, combined with static density functional theory calculations to probe the atomic interactions and electronic structure (e.g calculating the potential of zero charge). By performing detailed structural analyses of both metal-surface atoms and the near-surface water we find, amongst other insights: i) water adsorption causes significant surface roughening, ii) strongly adsorbed water covers only one quarter of available surface sites and iii) adsorbed water avoids clustering on the surface. Static calculations are used to gain a deeper understanding of the structuring observed in MD. For example, we use an energy decomposition analysis combined with calculated atomic charges to show adsorbate clustering is unfavorable due to Coulombic repulsion between adsorption site surface atoms. Results are discussed in the context of previous simulations of metal/water interfaces. The largest differences for the Mg(0001)/water system appear to be the high degree of surface distortion and minimal difference between the metal work function and metal/water potential of zero charge. The structural information in this paper is important for understanding aqueous Mg corrosion, as the Mg(0001)/water interface is the starting point for key reactions. Furthermore, our focus on understanding the driving forces behind this structuring leads to important insights for general metal/water interfaces.

Bertoni A, Fogarty R, Sanchez C,
et al., 2022, QM/MM optimization with quantum coupling: Host–guest interactions in a pentacene-doped p-terphenyl crystal, *The Journal of Chemical Physics*, Vol: 156, ISSN: 0021-9606

In this work, we present a novel force-based scheme to perform hybrid quantum mechanics/molecular mechanics (QM/MM) computations. The proposed scheme becomes especially relevant for the simulation of host–guest molecular systems, where the description of the explicit electronic interactions between a guest molecule and a classically described host is of key importance. To illustrate its advantages, we utilize the presented scheme in the geometry optimization of a technologically important host–guest molecular system: a pentacene-doped p-terphenyl crystal, a core component of a room-temperature MASER device. We show that, as opposed to the simpler and widely used hybrid scheme ONIOM, our Quantum-Coupling QM/MM scheme was able to reproduce explicit interactions in the minimum energy configuration for the host–guest complex. We also show that, as a result of these explicit interactions, the host–guest complex exhibits an oriented net electric dipole moment that is responsible for red-shifting the energy of the first singlet–singlet electronic excitation of pentacene.

Zauchner MG, Dal Forno S, Cśanyi G,
et al., 2021, Predicting polarizabilities of silicon clusters using local chemical environments, *Machine Learning: Science and Technology*, Vol: 2, Pages: 1-16, ISSN: 2632-2153

Calculating polarizabilities of large clusters with first-principles techniques is challenging because of the unfavorable scaling of computational cost with cluster size. To address this challenge, we demonstrate that polarizabilities of large hydrogenated silicon clusters containing thousands of atoms can be efficiently calculated with machine learning methods. Specifically, we construct machine learning models based on the smooth overlap of atomic positions (SOAP) descriptor and train the models using a database of calculated random-phase approximation polarizabilities for clusters containing up to 110 silicon atoms. We first demonstrate the ability of the machine learning models to fit the data and then assess their ability to predict cluster polarizabilities using k-fold cross validation. Finally, we study the machine learning predictions for clusters that are too large for explicit first-principles calculations and find that they accurately describe the dependence of the polarizabilities on the ratio of hydrogen to silicon atoms and also predict a bulk limit that is in good agreement with previous studies.

Bustamante C, Gadea E, Horsfield A,
et al., 2021, Dissipative equation of motion for electromagnetic radiation in quantum dynamics, *Physical Review Letters*, Vol: 126, ISSN: 0031-9007

The dynamical description of the radiative decay of an electronically excited state in realistic many-particle systems is an unresolved challenge. In the present investigation electromagnetic radiation of the charge density is approximated as the power dissipated by a classical dipole, to cast the emission in closed form as a unitary single-electron theory. This results in a formalism of unprecedented efficiency, critical for ab initio modeling, which exhibits at the same time remarkable properties: it quantitatively predicts decay rates, natural broadening, and absorption intensities. Exquisitely accurate excitation lifetimes are obtained from time-dependent DFT simulations for C2+, B+, and Be, of 0.565, 0.831, and 1.97 ns, respectively, in accord with experimental values of 0.57±0.02, 0.86±0.07, and 1.77–2.5 ns. Hence, the present development expands the frontiers of quantum dynamics, bringing within reach first-principles simulations of a wealth of photophysical phenomena, from fluorescence to time-resolved spectroscopies.

Bustamante C, Todorov T, Sanchez C,
et al., 2020, A simple approximation to the electron-phonon interaction in population dynamics, *Journal of Chemical Physics*, Vol: 153, ISSN: 0021-9606

The modeling of coupled electron–ion dynamics including a quantum description of the nuclear degrees of freedom has remained a costly and technically difficult practice. The kinetic model for electron–phonon interaction provides an efficient approach to this problem, for systems evolving with low amplitude fluctuations, in a quasi-stationary state. In this work, we propose an extension of the kinetic model to include the effect of coherences, which are absent in the original approach. The new scheme, referred to as Liouville–von Neumann + Kinetic Equation (or LvN + KE), is implemented here in the context of a tight-binding Hamiltonian and employed to model the broadening, caused by the nuclear vibrations, of the electronic absorption bands of an atomic wire. The results, which show close agreement with the predictions given by Fermi’s golden rule (FGR), serve as a validation of the methodology. Thereafter, the method is applied to the electron–phonon interaction in transport simulations, adopting to this end the driven Liouville–von Neumann equation to model open quantum boundaries. In this case, the LvN + KE model qualitatively captures the Joule heating effect and Ohm’s law. It, however, exhibits numerical discrepancies with respect to the results based on FGR, attributable to the fact that the quasi-stationary state is defined taking into consideration the eigenstates of the closed system rather than those of the open boundary system. The simplicity and numerical efficiency of this approach and its ability to capture the essential physics of the electron–phonon coupling make it an attractive route to first-principles electron–ion dynamics.

Fogarty R, Smutna J, Wenman M,
et al., 2020, Beyond two-center tight binding: Models for Mg and Zr, *Physical Review Materials*, Vol: 4, ISSN: 2475-9953

We describe a systematic approach to building ab initio tight-binding models and apply this to hexagonal metals Mg and Zr. Our models contain three approximations to plane-wave density functional theory (DFT): (i) we use a small basis set, (ii) we approximate self-consistency, and (iii) we approximate many-center exchange and correlation effects. We test a range of approximations for many-center exchange-correlation and self-consistency to gauge the accuracy of each in isolation. This systematic approach also allows us to combine multiple approximations in the optimal manner for our final tight-binding models. Furthermore, the breakdown of errors into those from individual approximations is expected to be a useful guide for which approximations to include in other tight-binding models. We attempt to correct any remaining errors in our models by fitting a pair potential. Our final tight-binding model for Mg shows excellent agreement with plane-wave results for a wide range of properties (e.g., errors below 10% for self-interstitial formation energies and below 3% for equilibrium volumes) and is expected to be highly transferable due to the minimal amount of fitting. Calculations with our Zr model also show good agreement with plane-wave results (e.g., errors below 2% for equilibrium volumes) except for properties where self-consistency is important, such as self-interstitial formation energies. However, for these properties we are able to generate a tight-binding model which shows excellent agreement with non-self-consistent DFT with a small basis set (i.e., many-center effects are captured accurately by our approximations). As we understand the source of remaining errors in our Zr model we are able to outline the methods required to build upon it to describe the remaining properties with greater accuracy.

Smutna J, Fogarty RM, Wenman MR,
et al., 2020, Systematic development of ab initio tight-binding models for hexagonal metals, *Physical Review Materials*, Vol: 4, Pages: 043801-1-043801-18, ISSN: 2475-9953

A systematic method for building an extensible tight-binding model from ab initio calculations has been developed and tested on two hexagonal metals: Zr and Mg. The errors introduced at each level of approximation are discussed and quantified. For bulk materials, using a limited basis set of spd orbitals is shown to be sufficient to reproduce with high accuracy bulk energy versus volume curves for fcc, bcc, and hcp lattice structures, as well as the electronic density of states. However, the two-center approximation introduces errors of several tenths of eV in the pair potential, crystal-field terms, and hopping integrals. Environmentally dependent corrections to the former two have been implemented, significantly improving the accuracy. Two-center hopping integrals were corrected by taking many-center hopping integrals for a set of structures of interest, rotating them into the bond reference frame, and then fitting a smooth function through these values. Finally, a pair potential was fitted to correct remaining errors. However, this procedure is not sufficient to ensure transferability of the model, especially when point defects are introduced. In particular, it is shown to be problematic when interstitial elements are added to the model, as demonstrated in the case of octahedral self-interstitial atoms.

Cui Y, King DJM, Horsfield AP,
et al., 2020, Solidification orientation relationships between Al3Ti and TiB2, *Acta Materialia*, Vol: 186, Pages: 149-161, ISSN: 1359-6454

Orientation relationships (ORs) can form during solidification by a variety of mechanisms that are often difficult to distinguish after solidification. Here we study three ORs formed by the nucleation of Al3Ti on TiB2, and by the pushing and engulfment of TiB2 by growing Al3Ti facets in hyperperitectic Al-rich melts. The nucleation OR is identified by growing a relatively large TiB2 crystal, solidifying multiple small Al3Ti crystals on one (0001) facet of TiB2, and measuring the resulting OR by electron backscatter diffraction (EBSD). Pushing and engulfment ORs are investigated by statistical analysis of EBSD measurements, density functional theory (DFT) calculations of interface energies, and imaging of cross-sections of TiB2 particles being pushed and engulfed by Al3Ti facets. It is shown that the lowest energy OR is formed by nucleation as well as by pushing/engulfment. The higher energy ORs, formed by pushing and engulfment, correspond to local interfacial energy minima and can be explained by rotation of TiB2 particles on Al3Ti facets during pushing.

Wells T, Horsfield A, Foulkes WMC,
et al., 2019, The microscopic Einstein-de Haas effect, *Journal of Chemical Physics*, Vol: 150, ISSN: 0021-9606

The Einstein-de Haas (EdH) effect, where the spin angular momentum of electrons is transferred to the mechanical angular momentum of atoms, was established experimentally in 1915. While a semiclassical explanation of the effect exists, modern electronic structure methods have not yet been applied to model the phenomenon. In this paper, we investigate its microscopic origins by means of a noncollinear tight-binding model of an O2 dimer, which includes the effects of spin-orbit coupling, coupling to an external magnetic field, and vector Stoner exchange. By varying an external magnetic field in the presence of spin-orbit coupling, a torque can be generated on the dimer, validating the presence of the EdH effect. The avoided energy level crossings and the rate of change of magnetic field determine the evolution of the spin. We also find that the torque exerted on the nuclei by the electrons in a time-varying B field is not only due to the EdH effect. The other contributions arise from field-induced changes in the electronic orbital angular momentum and from the direct action of the Faraday electric field associated with the time-varying magnetic field.

Todorov TN, Horsfield AP, 2019, Multiple-probe electronic open boundaries with bad contacts, *Physical Review B*, Vol: 99, ISSN: 2469-9950

We revisit the weak-coupling limit of the hairy probes method for electronic open boundaries. In this limit, the electronic density matrix is approximately stationary. We exploit this fact to combine hairy probes with electron-phonon scattering at the level of low-order transitions between eigenstates of the electronic Hamiltonian. This provides a method that is computationally very efficient, and whose results can be interpreted quite straightforwardly. The resultant time-dependent hybrid method is illustrated through the numerical calibration of a thermoelectric nanoscale thermometer.

Coury MEA, Dudarev SL, Foulkes WMC,
et al., 2018, Erratum: Hubbard-like Hamiltonians for interacting electrons in s, p, and d orbitals (vol 93, 075101, 2016), *Physical Review B*, Vol: 98, ISSN: 2469-9950

Thong A, Shaffer M, Horsfield AP, 2018, Rectification and negative differential resistance via orbital level pinning, *Scientific Reports*, Vol: 8, ISSN: 2045-2322

A donor-acceptor system, 4-thiophenyl-azafulleroid (4TPA-C60), is investigated at the point of HOMO/LUMO resonance and beyond to understand how negative differential resistance (NDR) features may be observed in such systems. Our previous investigation showed that charge transfer between the occupied and unoccupied states at resonance hindered crossing of the HOMO and LUMO levels, thus preventing the formation of an NDR feature. In this work, it is shown that the negative differential resistance feature of 4TPA-C60 can be tailored based on the couplings at the metal/molecule interface. Ab initio calculations show that limited charge extraction from atomically sharp contacts results in a HOMO-LUMO pinning effect which delays the onset of the NDR feature. Subsequent unpinning of the states can only occur when additional charge extraction channels enter the bias window, highlighting an important role which non-frontier states play in charge transport. The proposed charge transfer mechanism is then exploited by introducing a fluorine atom into the C60 cage to tune the energies of the acceptor, and narrow the width of the current peak. These findings not only demonstrate the importance of the metal/molecule interface in the design of molecular electronic architectures but also serve to inform future design of molecular diodes and RTDs.

Charlton RJ, Fogarty R, Bogatko S,
et al., 2018, Implicit and explicit host effects on excitons in pentacene derivatives, *Journal of Chemical Physics*, Vol: 148, ISSN: 0021-9606

Anab initiostudy of the effects of implicit and explicit hosts on the excited state properties ofpentacene and its nitrogen-based derivatives has been performed using ground state density func-tional theory (DFT), time-dependent DFT and ∆SCF. We observe a significant solvatochromicredshift in the excitation energy of the lowest singlet state (S1) of pentacene from inclusion inap-terphenyl host compared to vacuum; for an explicit host consisting of six nearest neighbourp-terphenyls, we obtain a redshift of 65 meV while a conductor-like polarisable continuum model(CPCM) yields a 78 meV redshift. Comparison is made between the excitonic properties of pen-tacene and four of its nitrogen-based analogues, 1,8-, 2,9-, 5,12-, and 6,13-diazapentacene with thelatter found to be the most distinct due to local distortions in the ground state electronic struc-ture. We observe that a CPCM is insufficient to fully understand the impact of the host due tothe presence of a mild charge-transfer (CT) coupling between the chromophore and neighbouringp-terphenyls, a phenomenon which can only be captured using an explicit model. The strengthof this CT interaction increases as the nitrogens are brought closer to the central acene ring ofpentacene.

Zauchner M, Horsfield AP, todorov T, 2018, Efficient electron open boundaries for simulating electrochemical cells, *Physical Review B: Condensed Matter and Materials Physics*, Vol: 97, Pages: 1-1, ISSN: 1098-0121

Nonequilibrium electrochemistry raises new challenges for atomistic simulation: we need to perform molecular dynamics for the nuclear degrees of freedom with an explicit description of the electrons, which in turn must be free to enter and leave the computational cell. Here we present a limiting form for electron open boundaries that we expect to apply when the magnitude of the electric current is determined by the drift and diffusion of ions in a solution and which is sufficiently computationally efficient to be used with molecular dynamics. We present tight-binding simulations of a parallel-plate capacitor with nothing, a dimer, or an atomic wire situated in the space between the plates. These simulations demonstrate that this scheme can be used to perform molecular dynamics simulations when there is an applied bias between two metal plates with, at most, weak electronic coupling between them. This simple system captures some of the essential features of an electrochemical cell, suggesting this approach might be suitable for simulations of electrochemical cells out of equilibrium.

Xu W, Horsfield AP, Wearing D,
et al., 2017, Classical and quantum calculations of the temperature dependence of the free energy of argon, *Computational Materials Science*, Vol: 144, Pages: 36-41, ISSN: 0927-0256

The free energy is central to statistical mechanics and thermodynamics, and its accurate calculation via. computational modelling is important for a large number of applications, especially when its experimental value is hard to obtain. Several established and general methods for calculating the Helmholtz free energy across different length scales, including continuum, atomistic and quantum mechanical, are compared and analyzed. A computational approach is then proposed to calculate the temperature dependences of internal energy and absolute Helmholtz free energy for solid and liquid phases with the coupling of thermodynamic integration (TI) and harmonic approximation calculations from both classical molecular dynamics (MD) and density functional theory (DFT). We use the Lennard-Jones system as an example (i.e. argon) for the demonstration of the approach. It is observed that the free energy transits smoothly from being describable by the harmonic approximation to including anharmonic effects at a transition temperature around 0.56 Tm; below this temperature, the quantum behavior of atoms is important. At higher temperatures (T > 0.56 Tm), the TI and harmonic approximation results for the Helmholtz free energy functions become increasingly divergent with the increase of temperature. This work demonstrates that a multiscale approach employing TI, MD, and DFT can provide accurate calculations of the temperature dependence of absolute Helmholtz free energy for both solid and liquid phases.

Horsfield AP, Haase A, Turin L, 2017, Molecular recognition in olfaction, *Advances in Physics: X*, Vol: 2, Pages: 937-977, ISSN: 2374-6149

The mechanism by which the chemical identity of odourants is established by olfactory receptors is a matter of intense debate. Here we present an overview of recent ideas and data with a view to summarising what is known, and what has yet to be determined. We outline the competing theories, and summarise experimental results employing isotopes obtained for mammals, insects, and individual receptors that enable us to judge the relative correctness of the theories.

boleininger M, Horsfield AP, 2017, Efficient local-orbitals based method for Ultrafast Dynamics, *Journal of Chemical Physics*, Vol: 147, ISSN: 1089-7690

Computer simulations are invaluable for the study of ultrafast phenomena, as they allow us to directly access the electron dynamics. We present an efficient method for simulating the evolution of electrons in molecules under the influence of time-dependent electric fields, based on the Gaussian tight binding model. This model improves upon standard self-charge-consistent tight binding by the inclusion of polarizable orbitals and a self-consistent description of charge multipoles. Using the examples of bithiophene, terthiophene, and tetrathiophene, we show that this model produces electrostatic, electrodynamic, and explicitly time-dependent properties in strong agreement with density-functional theory, but at a small fraction of the cost.

Thong AZ, Shaffer MS, Horsfield AP, 2017, HOMO-LUMO coupling: the fourth rule for highly effective molecular rectifiers, *Nanoscale*, Vol: 9, Pages: 8119-8125, ISSN: 2040-3372

Three rules for creating highly effective unimolecular rectifiers that utilize asymmetric anchoring groups have been proposed by Van Dyck and Ratner [Ratner et al., Nano Lett., 2015, 15, 1577–1584]. This study investigates their proposed rectification mechanism in a functionalised azafullerene system (4TPA–C60) and identifies a fourth rule. NEGF-DFT shows that 4TPA–C60 fulfills the three design rules and finds that a saturated bridge is not required to fulfil the third rule, contrary to previous belief. Instead a twisted-π bridge decouples the donor and acceptor states whilst still providing a high conductance pathway. The molecular junction has a calculated rectification ratio of 145 at a bias of ±1 V and the U-type rectification mechanism is driven by the pinning of the HOMO to the LUMO when the device is forward biased, but not when reverse biased. The switching behaviour is a result of a charge dipole forming at different interfaces for different bias directions. An additional design rule is thus proposed: charge transport should allow bias dependent coupling of filled to unfilled states. The findings in this work not only help in understanding charge transport in molecular rectifiers, but also have wider implications for the design of molecular resonant tunneling devices.

Todorov TN, Cunningham B, Dundas D,
et al., 2017, Non-conservative forces in bulk systems, *Materials Science and Technology*, Vol: 33, Pages: 1442-1446, ISSN: 1316-2012

The ability of electrons and phonons in solids to equilibrate is fundamental to the thermodynamic andtransport properties of these systems. Nonetheless in mesoscopic systems sufficiently high current densitiescan lead to situations where the phonon subsystem can no longer reach a steady state. In previous workthis phenomenon was connected with the demonstration that interatomic forces under current are non-conservative, with the ability to do net work around closed paths in the nuclear subspace. Here we considerthese effects in bulk systems. We arrive at a critical current density beyond which current flow results in theuncompensated stimulated emission of a blast of forward-travelling phonons. The resultant atomic motion isillustrated with model simulations of long atomic wires under current. While the critical current density forthese effects is very high compared with those in electroplasticity phenomena, it is hoped that the paper willstimulate further research into non-conservative dynamics in extended conductors and the possible relevanceof the effect to electropulsing.

Boleininger M, Guilbert A, Horsfield AP, 2016, Gaussian polarizable-ion tight binding, *Journal of Chemical Physics*, Vol: 145, ISSN: 1089-7690

To interpret ultrafast dynamics experiments on large molecules, computer simulation is required due to the complex response to the laser field. We present a method capable of efficiently computing the static electronic response of large systems to external electric fields. This is achieved by extending the density-functional tight binding method to include larger basis sets and by multipole expansion of the charge density into electrostatically interacting Gaussian distributions. Polarizabilities for a range of hydrocarbon molecules are computed for a multipole expansion up to quadrupole order, giving excellent agreement with experimental values, with average errors similar to those from density functional theory, but at a small fraction of the cost. We apply the model in conjunction with the polarizable-point-dipoles model to estimate the internal fields in amorphous poly(3-hexylthiophene-2,5-diyl).

Horsfield AP, Boleininger M, D'Agosta R,
et al., 2016, Efficient simulations with electronic open boundaries, *Physical Review B*, Vol: 94, ISSN: 1550-235X

We present a reformulation of the Hairy Probe method for introducing electronic open boundariesthat is appropriate for steady state calculations involving non-orthogonal atomic basis sets. As acheck on the correctness of the method we investigate a perfect atomic wire of Cu atoms, and aperfect non-orthogonal chain of H atoms. For both atom chains we find that the conductance hasa value of exactly one quantum unit, and that this is rather insensitive to the strength of couplingof the probes to the system, provided values of the coupling are of the same order as the meaninter-level spacing of the system without probes. For the Cu atom chain we find in addition thataway from the regions with probes attached, the potential in the wire is uniform, while withinthem it follows a predicted exponential variation with position. We then apply the method to aninitial investigation of the suitability of graphene as a contact material for molecular electronics.We perform calculations on a carbon nanoribbon to determine the correct coupling strength of theprobes to the graphene, and obtain a conductance of about two quantum units corresponding totwo bands crossing the Fermi surface. We then compute the current through a benzene moleculeattached to two graphene contacts and find only a very weak current because of the disruption ofthe π-conjugation by the covalent bond between the benzene and the graphene. In all cases we findthat very strong or weak probe couplings suppress the current.

Xu W, Horsfield AP, Wearing D,
et al., 2016, Diversification of MgO//Mg interfacial crystal orientations during oxidation: A density functional theory study, *Journal of Alloys and Compounds*, Vol: 688, Pages: 1233-1240, ISSN: 1873-4669

In this work we use computer simulations to explain the variety of crystal orientations observed at interfaces between MgO and Mg when Mg single crystals are oxidized. Using first-principles density functional theory simulations we investigate the interfacial stability of MgO//Mg interfaces, and find that a combination of interfacial chemical bonding energy and epitaxial strain stored in the oxide layers can change the relative stability of competing MgO//Mg interfaces. We propose that a combination of the oxygen chemical potential at the interface plane and the epitaxial strain energy stored in the oxide layers is responsible for the differences in observed interfacial crystal orientations–a key insight for the design and development of Mg alloys reinforced by MgO particles.

This data is extracted from the Web of Science and reproduced under a licence from Thomson Reuters. You may not copy or re-distribute this data in whole or in part without the written consent of the Science business of Thomson Reuters.