## Abstract

The monotonic and non-monotonic variations of atomic properties within and between the rows of the periodic table underlie our understanding of chemistry and accounting for these variations has been a signature strength of quantum mechanics (QM). However, the computational burden of QM motivates the development of more efficient means of describing electrons and reactivity. The recently developed LEWIS^{•} model incorporates lessons learnt from QM into a force field that includes electrons as explicit pseudo-classical particles. Here, we extend LEWIS^{•} across the 2*p* and 3*p* elements, and show that it is capable of reproducing both monotonic and non-monotonic variations of chemically important atomic properties in a cost-effective manner. An indicator of the strength of the construct is the ability of pairwise potentials trained on ionization energies and the order of spin configurations to predict atomic polarizabilities. In this manner, some insights of QM are uncoupled from its onerous computational burden.

## 1. Introduction

The challenge of describing the complex behaviour of bound electrons has spawned a spectrum of approaches, ranging from explicit to implicit treatments. In molecular mechanics, electrons are typically represented implicitly by energies that depend not only on the partial charge carried by each atom, but also on bond lengths, bond angles, dihedral angles and, most recently, polarization, all of which are calculated on the fly. Such force fields include OPLS [1], CHARMM [2] and AMBER [3]. For reactive force fields, in which bonds can be made and broken [4–6], electrons are also implicit in the dependence of energies on bond types, which are also calculated on the fly. Examples include ReaxFF [7], REBO [8,9] and EVB [10,11]. In essence, the effects of the electrons are captured empirically by many-body interactions among the atoms.

At the other end of the spectrum, electrons are described explicitly and rigorously by quantum mechanics (QM). An early triumph was its explanation of the variations in the properties of the chemical elements within and between rows of the periodic table. As we teach introductory chemistry students, (i) our understanding of the chemical behaviour of the elements relies heavily on the variation of the electron affinity (EA) and ionization energy (IE), and (ii) while nuclear charge increases monotonically, the EA and IE vary non-monotonically according to the half- and complete-filling of subshells of atomic orbitals that are (a) defined by wave mechanics and (b) admit only one electron of each spin. Not only is this a complete and satisfying description of the fundamental properties underlying the behaviour of the elements, but QM has since gone from strength to strength, finding ever broader applications.

It remains, however, that the complexity and cost of QM calculations seriously limit its applicability. For systems with numerous electrons, reliability is compromised by the shortcuts that must be taken in basis sets and configurations. Even in its more efficient density functional (DFT) incarnation, QM studies of reactions generally require restricting degrees of freedom by hypothesizing reaction pathways. However, especially in condensed phases, the relevance of numerous degrees of freedom is often a central issue. Given this uncharted territory, a significant contribution can be made by an approach capable of efficiently sampling and defining reaction pathways in complex media.

An alternative approach is to retain explicit electrons, but lean hard towards the more intuitive particle side of their wave–particle duality, while incorporating features that arise from their wave behaviour. The effect is to render the charge and spin densities of DFT in pointillist fashion (to borrow a term from the art world). The question is whether such a coarse description can do useful chemical work. To assess the possibilities, we need to start where QM started, i.e. with an account of the atomic properties that are important for ultimately describing molecular behaviour. Recently, we developed potentials for repulsions between pseudo-particles representing like-spin and unlike-spin electrons. After designing them on the basis of data for two elements, we then showed that they could also describe aufbau for three other elements scattered among the first three rows of the periodic table [12]. As the next step, we now wish to fill in the gaps in order to explore the potential for describing trends in periodic properties. Of particular interest, in addition to the ionization and spin excitation energies that contribute to aufbau, are the atomic polarizabilities. While a great deal of effort has gone into adding polarizabilities to force fields with implicit electrons, polarizabilities are natural features of models with explicit electrons, and we take advantage of that here. The combination of ionization energies (IEs), EAs and polarizabilities is a prerequisite for a reliable description of molecules because they respectively reflect the depth and width of the electron potential well.

## 2. LEWIS^{•}

### (a) Particles

LEWIS^{•} unpacks atoms into fully charged valence electrons (*e*) and atomic kernels (*X*) that include the core electrons. The independent mobility of these particles naturally builds in reactivity (i.e. the breaking and making of bonds, whether in homolytic or heterolytic fashion), as well as polarizability. In addition to Cartesian coordinates, each electron particle has one of two possible spins, which enables the description of different spin configurations, and a variable cloud diameter (i.e. uncertainty of position), with an associated quantum kinetic energy,
*m*_{e} the electron mass and *d*_{e} the electron cloud diameter.

### (b) Interactions and potential forms

The interactions between the particles are strictly pairwise to maximize efficiency, and smoothly differentiable to ensure energy conservation in simulations. All are coulombic at long distances, to conform to classical electrostatics, and finite at short distances on account of the diffuse distribution of the electrons [12]. Of course, numerous functions have these features and, as for the functionals in DFT, the quality of a given set of potentials must be evaluated empirically. The search for suitable potential forms is the most demanding aspect of developing the model. The following set of forms is the most effective that we have discovered among the many combinations tested so far.

The dependence of the pair interactions on interparticle distance, *r*, is expected to depend on how diffuse the charge associated with each particle is. Assuming a simple scaling relationship, we define *r*_{eff}=*r*/*d*_{eff}, where the *d*_{eff} is a monotonically increasing function of the characteristic diameters of the two charge distributions. For the energy of attraction between an electron and an atomic kernel of element *X*
*κ*_{eX} are adjustable parameters specific for the element *X* and the best form found for *d*_{eff} is
*σ*_{eX} scaling the diameter, *d*_{e}, of the valence electron and *ρ*_{eX} corresponding to the cloud diameter of the implicit core electrons. *ρ*_{eX} is set to zero for hydrogen and is found to correlate with the atomic radii of other elements (see below).

For the energy of repulsion between a pair of unlike-spin electrons, the form is similar
*d*_{eff} is
*d*_{e,i} and *d*_{e,j} are the cloud diameters of the two electrons. For repulsions between like-spin electrons, the consequences of the fermionic requirement of wave function antisymmetry are accounted for by adding an ‘exchange term’ for which the best form found so far is
*χ*_{ee}, *τ*_{ee1} and *τ*_{ee2} determine the amplitude, zero crossing and damping, respectively. As discussed below, an important feature of equation (2.6) is the change of sign with distance.

## 3. Training set and optimization

### (a) Training set

To describe a system, optimized potentials are required for the interactions between all pairs of particles. Each monoatomic species (containing just a single kernel) requires *U*_{ee(l)}, *U*_{ee(u)} and *U*_{eX}. However, to be suitably transferable, *U*_{ee(l)} and *U*_{ee(u)} must be independent of environment. Here, we use the *U*_{ee(l)} and *U*_{ee(u)}, obtained earlier by training on properties of hydrogen and chlorine atoms [12]. Only, the parameter values for *U*_{eX} are optimized for each newly introduced element, *X*, as we did earlier for three scattered elements, boron, fluorine and aluminium.

Each *U*_{eX} has five parameters (three in equation (2.2) and two in equation (2.3)). To optimize these parameters, we built training sets that include experimental data for the atoms in all degrees of ionization of their valence electrons: in particular, (i) the differences between their ground-state energies and (ii) the minimum energies required for each to achieve other spin configurations. We chose these data for training because they are important for reactivity and not addressed by conventional classical models. They are also highly demanding owing to the non-monotonic dependence of the energies on total electron spin for a given ionization state or on atomic number within a row of the periodic table.

Table 1 shows all the experimentally documented configurations that LEWIS^{•} can describe, with the order of energy following Hund's rules. In the case of sulfur, this makes a total of 14 data points in the training set: seven for the IEs and EA, and seven for the variations in total electron spin. The relatively small training sets of carbon and silicon each have a total of nine data points. Elements of groups 1A and 2A are not included because they lack sufficient ionization and spin excitation data for this approach to parameter optimization.

### (b) Parameter optimization

Our objective function is a weighted sum of the absolute value of the percent deviation from the experimental energy values in the training set. Because the energy depends on structure, evaluation of the objective function requires an additional layer of optimization: for each set of potential parameters, the structural parameters (coordinates of all particles and cloud sizes of all electrons) are varied to find the minimum energy for each species (as described in the electronic supplementary material).

The weights in the objective function were assigned according to the anticipated importance of each training variable for future simulations of reactions. The greatest weights were given to the EAs and first IEs. In addition, while low weights were given to the energies of the spin excitations, high priority was given to obtaining the correct sequence of spin states.

## 4. Results

### (a) Replication of training data

Table 2 summarizes the optimized values of the parameters for all the pair potentials, with the new *U*_{eX}s shown in italics, and figure 1 shows *U*_{ee(l)}, *U*_{ee(u)} and *U*_{eS}, as an example of *U*_{eX}. As illustrated in figure 1, all three interactions are coulombic at long distances, to conform to classical electrostatics. However, they are finite at short distances, on account of quantum uncertainty. Exchange effects are reflected in the difference between the repulsions between electrons of like-spin (dashed line) and unlike-spin (solid line). Here, the sign change (observed as the crossover between the dashed and solid black lines) is an important feature that provides for Pauli exclusion at short distances and Hund's rules (table 1) at longer distances.

The quality of the fit to the training data is illustrated in figure 2 for sulfur. All IEs, especially the EA, are fit very accurately. Energies for changes in total spin are well replicated except that the second spin excitation for the neutral atom (with six valence electrons) is significantly underestimated. This excitation promotes an electron to a new shell which is challenging for LEWIS^{•}. It is possible that this is due to the restriction to pair-wise interactions or the coarse graining of the density inherent in the pointillist description. However, it remains that the spin orders are always correct and that is our main interest for future applications.

Similar accuracies are obtained for all the other elements. Figure 3*a* (as well as tables S1 and S2 in the electronic supplementary material) shows the ability of LEWIS^{•} to simultaneously fit all the EAs and IEs in the training set, under the constraint that the order of the spin states is also reproduced, using the single set of electron–electron repulsions shown in the top panel of figure 1. The extra weights given to the EA and the first ionization energy are reflected in the quality of the fit.

Surveying the trends in figure 3*a*, it is notable that, even though a lower priority was given to higher-order IEs, the values for all the oxidation states are not only in good agreement with the individual experimental values, but also follow the non-monotonic experimental variations, including those within rows of the periodic table (figure 3*a*). The one exception is the experimentally observed dip in first ionization energy from nitrogen to oxygen, although the corresponding dip in the third row is successfully replicated, as are all the corresponding dips in the higher-order IEs. Getting the trends right is important, not only for internal consistency, but also ultimately for simulating chemistry. It should be noted that by fitting the EAs and first ionization energies (IE1), we automatically replicate electronegativity on the Mulliken scale, calculated as (EA+IE1)/2.

### (b) Predictions outside the training set

By way of validation, we now ask whether the LEWIS^{•} potentials can predict trends in atomic properties outside the training set. In particular, we focus on atomic polarizabilities, as a growing appreciation of their importance for molecular behaviour has motivated extensive efforts to add them to conventional (non-reactive) atomistic force fields [17–26]. In LEWIS^{•}, polarizabilities occur naturally with field-induced redistribution of the explicit, independently mobile, valence electrons (as described in the electronic supplementary material). Although the LEWIS^{•} interaction potentials are trained on data that vary non-monotonically within rows of the periodic table (figure 3*a*), LEWIS^{•} correctly predicts monotonically varying polarizabilities, with values that compare well with experiment, as shown in figure 3*b* (and electronic supplementary material, table S3). The trends within rows and the jumps between rows are well matched, and the quantitative agreement is good for all but aluminium (where the calculated value is double the experimental value, and not shown in the graph). The better quantitative agreement at the right end of the rows may reflect the availability of larger training sets for these elements. Of course, polarizability can be added to smaller training sets to obtain potentials that better describe this property for elements of interest.

### (c) Description of structure

In its pointillist rendering of the electron charge and spin density, LEWIS^{•} provides ready insights into atomic structure. Figure 4*a* shows ground-state electronic structures for the elements considered, including a metal (to the left of the green line) and two metalloids (between the green and blue lines), whereas figure 4*b* shows ground- and excited-state electronic structures for the valence ions of sulfur. Throughout, we see that electrons of like-spin avoid each other, with distributions around the kernel that are reminiscent of the major lobes of hybrid orbitals: two three and four electrons of like-spin form linear, triangular and tetrahedral arrangements, respectively. This is most obvious in panel (*a*) for spin *α* (magenta) in the structures of boron, carbon, nitrogen and oxygen. In some cases, the symmetries of the two spin groups allow pairing between electrons of opposite spin. Examples include the *αβ* (S^{+4} ground state) and 4*α*3*β* (ground states of F, Cl and S^{−}) spin configurations. The energy-minimized symmetries are identical for all species with the same spin configuration. For example, the symmetries of nitrogen and phosphorus in panel (*a*) are not only identical to each other, but also to the ground state of S^{+}.

Quantitative details of the structures are shown in figure 3*c*. It is satisfying that the kernel radius is always larger for the 3*p* elements than for any of the 2*p* elements, just as is the case for covalent atomic radii [16]. The radii of valence electron clouds and the radial distances of the valence electrons from the kernel are the values that minimize the ground-state energy, as the electrons compromise between avoiding each other while approaching the kernel. Like experimental radii, they shrink monotonically across a row and expand going down a family. It is notable that the monotonic behaviour of radii within rows occurs in spite of the non-monotonic behaviour of the properties on which the potentials were trained (figure 3*a*).

## 5. Conclusion

The successes of LEWIS^{•} are significant on several levels. Philosophically, LEWIS^{•} is most closely related to the original, orbital-free conception of DFT. However, good orbital-free functionals have proven elusive and the Kohn–Sham orbital approach has become dominant based on its straightforward treatment of the kinetic energy. The present pointillist representation of the electron charge and spin density has a similarly straightforward evaluation of the kinetic energy, while also approximating other contributions to the energy with efficient and intuitive pairwise potentials. In these respects, it is similar to eFF [27]. However, because eFF assumes that the electrons occupy Gaussian orbitals, its potentials are restricted to forms that do not serve as well as the heuristic forms of LEWIS^{•}. Of course, as for functionals in DFT, future exploration in LEWIS^{•} can be expected to provide still better potentials than the ones that we have discovered so far.

Regarding the present progress, it is simultaneously surprising and encouraging that a single set of pseudo-classical pairwise potentials for repulsions between valence electrons is transferable across three rows of the periodic table to closely mimic the non-monotonic variation of EA and ionization energy (and therefore also electronegativity), while correctly ordering the energies of different spin states. Also non-trivial is the ability of LEWIS^{•} so trained to naturally predict polarizabilities and trends in atomic radii. We know of no other method that can produce this information with pairwise interactions.

Although one can strive for greater accuracy with improved potentials, the present trends in atomic radii, electronegativities and polarizabilities already suggest that the present approach can reproduce bond lengths, bond energies and bond polarizations well enough to justify the effort of developing kernel–kernel potentials. Optimism regarding kernel–kernel potentials is supported by the effectiveness of the kernel–kernel potentials developed in the earlier LEWIS model of water [28,29] which allowed for an excellent description of the solvation and dynamics of proton defects in bulk water [30]. Limitations of that construct were its modelling of the valence electrons in pairs with a fixed cloud radius. Electron pairs cannot describe free radicals or fractional-bond order, and a fixed electron cloud radius compromises transferability from oxygen to other non-metals. With LEWIS^{•}, both of these limitations are overcome.

## Data accessibility

The datasets supporting this article have been uploaded as part of the electronic supplementary material.

## Authors' contributions

S.E. carried out all the coding and computations, J.H. supervised the work and the manuscript was written jointly.

## Competing interests

We declare we have no competing interests.

## Funding

This work was supported in part by NIH grant no. EB001035 and NSF grant no. 1305713. Additional computational support was provided by the Brandeis HPC.

## Acknowledgements

We appreciate Chen Bai's input regarding the effectiveness of differential evolution for parameter optimization.

- Received June 3, 2015.
- Accepted July 28, 2015.

- © 2015 The Author(s)

Published by the Royal Society. All rights reserved.