2
$\begingroup$

I am working with lithium sulfate and am wondering how to figure out the ionic radius of the Li atoms. I did a DFT/B3LYP/G-311G ++ (d,p) calculation.

$\endgroup$
1
  • 4
    $\begingroup$ Ionic radius is not a well-defined quantum-mechanical quantity, so estimating it properly in an unbiased can be quite involved. $\endgroup$
    – Hayden S
    Jul 10, 2023 at 21:44

1 Answer 1

2
$\begingroup$

Does it really need to be done in Gaussian?

I found a paper(1) titled "Atomic and Ionic Radii of Elements 1-96" where they calculated it for a large fraction of the periodic table:

Abstract: Atomic and cationic radii have been calculated for the first 96 elements, together with selected anionic radii. The metric adopted is the average distance from the nucleus where the electron density falls to 0.001 electrons per bohr3, following earlier work by Boyd. Our radii are derived using relativistic all-electron density functional theory calculations, close to the basis set limit. They offer a systematic quantitative measure of the sizes of non-interacting atoms, commonly invoked in the rationalization of chemical bonding, structure, and different properties. Remarkably, the atomic radii as defined in this way correlate well with van der Waals radii derived from crystal structures. A rationalization for trends and exceptions in those correlations is provided.

To do it, they used the non-empirical hybrid-exchange correlation functional PBE0, together with a very large and uncontracted atomic natural orbital-relativistic correlation consistent(ANO-RCC) basis set and the Douglas–Kroll–Hess second-order scalar relativistic Hamiltonian. Since I dropped out of my master degree back in the pandemic I don't have access to a Gaussian anymore, but I posted here previously a question regarding volumes of surfaces with uniform electronic density with some examples. Now I adapted the examples to try reproduce these calculations of Rahm et al., at least for Li you are interested in, and also for the next three elements in the alkali metal group, Na, K and Rb.

My attempt was done using NWChem(2), a major free software computational chemistry package(3) available in the repositories from many linux distributions. The particular version I used was the one currently available from Fedora 38 repositories:

[user@fedora ~]$ dnf list installed | grep 'nwchem'
nwchem.x86_64                                        7.0.2-12.fc38                       @fedora                
nwchem-common.noarch                                 7.0.2-12.fc38                       @fedora                
nwchem-openmpi.x86_64                                7.0.2-12.fc38                       @fedora                
[user@fedora ~]$ 

So I created input files for single point energy calculations, using PBE0 and ANO-RCC, where each ion is centered at x=0, y=0, z=0 coordinates, as shown below. The largest value for radius calculated in the paper we are trying reproduce is a little under 3Å, so the input files also direct NWChem to output density cubes measuring 6Åx6Åx6Å, that should be enough to avoid any part of the volume with electronic density above 0.001 electrons per bohr³ to be left outside the cube. The grid spacings were also chosen so our voxel volume (21.4μų) is close to the one used in the paper (18.5μų):

[user@fedora ionic_radii2]$ more *.nw | cat
::::::::::::::
input_lithium.nw
::::::::::::::
echo

start molecule

title "Single Point | pbe0/ANO-RCC"
charge 1

memory total 2 gb

geometry units angstroms print xyz autosym
Li      0.000000    0.000000    0.000000
end

basis
  Li library ANO-RCC rel
end

relativistic
  douglas-kroll
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE lithium
 LimitXYZ
 -3.0 3.0 216
 -3.0 3.0 216
 -3.0 3.0 216
 spin total
 gaussian
 output lithium_density.cube
end

task dplot
::::::::::::::
input.nw
::::::::::::::
echo

start molecule

title "Single Point | pbe0/aug-pcS-1/pcS-0"
charge 0

memory total 2 gb

geometry units angstroms print xyz autosym
C      0.913139    0.913139    0.000000
H      0.676070    1.981751    0.000000
C      0.000000    0.000000    0.780217
C      0.000000    0.000000   -0.780217
H      1.378212   -1.576370    0.000000
C      0.334232   -1.247371    0.000000
H      1.981751    0.676070    0.000000
C     -1.247371    0.334232    0.000000
H     -0.405381   -2.054282    0.000000
H     -2.054282   -0.405381    0.000000
H     -1.576370    1.378212    0.000000
end

basis
  H library pcS-0
  C library aug-pcS-1
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE propellane
 LimitXYZ
 -4.0 4.0 10
 -4.0 4.0 10
 -2.0 2.0 10
 spin total
 gaussian
 output propel_density.cube
end

task dplot::::::::::::::
input_potassium.nw
::::::::::::::
echo

start molecule

title "Single Point | pbe0/ANO-RCC"
charge 1

memory total 2 gb

geometry units angstroms print xyz autosym
K      0.000000    0.000000    0.000000
end

basis
  K library ANO-RCC rel
end

relativistic
  douglas-kroll
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE Potassium
 LimitXYZ
 -3.0 3.0 216
 -3.0 3.0 216
 -3.0 3.0 216
 spin total
 gaussian
 output potassium_density.cube
end

task dplot
::::::::::::::
input_rubidium.nw
::::::::::::::
echo

start molecule

title "Single Point | pbe0/ANO-RCC"
charge 1

memory total 2 gb

geometry units angstroms print xyz autosym
Rb      0.000000    0.000000    0.000000
end

basis
  Rb library ANO-RCC rel
end

relativistic
  douglas-kroll
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE Rubidium
 LimitXYZ
 -3.0 3.0 216
 -3.0 3.0 216
 -3.0 3.0 216
 spin total
 gaussian
 output rubidium_density.cube
end

task dplot
::::::::::::::
input_sodium.nw
::::::::::::::
echo

start molecule

title "Single Point | pbe0/ANO-RCC"
charge 1

memory total 2 gb

geometry units angstroms print xyz autosym
Na      0.000000    0.000000    0.000000
end

basis
  Na library ANO-RCC rel
end

relativistic
  douglas-kroll
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE Sodium
 LimitXYZ
 -3.0 3.0 216
 -3.0 3.0 216
 -3.0 3.0 216
 spin total
 gaussian
 output sodium_density.cube
end

task dplot
[user@fedora ionic_radii2]$ 

This page has instructions on how to install and run NWChem under Fedora. To install it there, you can run this command in the terminal:

dnf install nwchem-openmpi

Once installed, I did a parallel run on each input, using 8 threads on a AMD Ryzen 7 1700 machine, as follows:

module load mpi/openmpi
mpirun -np 8 nwchem_openmpi input_lithium.nw > output_Li.out
mpirun -np 8 nwchem_openmpi input_sodium.nw > output_Na.out
mpirun -np 8 nwchem_openmpi input_potassium.nw > output_K.out
mpirun -np 8 nwchem_openmpi input_rubidium.nw > output_Rb.out

Once finished, I checked the outputs produced no error and the cube files were properly recorded. The runtimes are also shown, by printing the last line of each output file:

[user@fedora ionic_radii2]$ grep -i 'error' *.out
[user@fedora ionic_radii2]$ tail -n 1 *.out
==> output_K.out <==
 Total times  cpu:     2009.8s     wall:     2116.0s

==> output_Li.out <==
 Total times  cpu:       49.0s     wall:      102.6s

==> output_Na.out <==
 Total times  cpu:      394.4s     wall:      474.6s

==> output_Rb.out <==
 Total times  cpu:    11640.9s     wall:    11949.9s
[user@fedora ionic_radii2]$ ls -l -s *.cube
131428 -rw-r--r--. 1 user user 134580624 jul 16 12:09 lithium_density.cube
131428 -rw-r--r--. 1 user user 134580626 jul 16 13:42 potassium_density.cube
131428 -rw-r--r--. 1 user user 134580625 jul 16 17:02 rubidium_density.cube
131428 -rw-r--r--. 1 user user 134580623 jul 16 13:04 sodium_density.cube
[user@fedora ionic_radii2]$

The timings show computational cost scales quickly with the size of the system. While the calculation for Li+ takes a bit over 100s wall time, the one for Rb+ takes over 3h. My initial idea was to try reproduce the results for all the alkali metals, but when I saw how quick processing time increased I decided to stop at rubidium.

Now the remaining problem is to do volume integration on the electron densities data present in the cube files. In the original study they used the partitioning method decribed in this paper(4). I sticked to the less sophisticated method I used in my past question, of just counting points above the given threshold of 0.001 electrons per bohr³, and then calculating the ion volume as a ratio of the box volume. In the original study, that covers 96 elements, they raise the issue of most of them except those in groups 1,2 and 18 have non-spherical grids, what complicates estimating average radii. But as all elements we are dealing with are from group 1, we can estimate the radii from the respective volumes by just assuming them spherical. The following Python 3 code implements it (must be saved to a .py file and run in the same folder the cube files are located):

from math import pi
from os import listdir
import re

def main():
    cubes = [file for file in listdir() if 'density.cube' in file]
    for file in cubes:
        with open(file, 'r') as f:
            text = f.read()

        ##Matches data in format the cube file uses to store
        ##density values, like 0.77410E-07 or 0.38588E+01.
        p = re.compile(r"0\.\d{5}E[-+]\d{2}")
        
        ## From the LimitXYZ section in the input files. If you change the cube volume, this must be updated.
        cube_volume = 6*6*6
        isosurface_threshold = 0.001

        data = p.findall(text)
        numeric_data = [float(value) for value in data]
        fraction_inside = len([density
                               for density in numeric_data
                               if density >= isosurface_threshold]) / len(numeric_data)

        bound_volume = cube_volume * fraction_inside
        ## From the sphere volume formula.
        spheric_radius = pow((3*bound_volume)/(4*pi), 1/3)
        print('For',
              file,
              'the calculated average radius is',
              round(spheric_radius, 2),
              'angstroms')

main()

It then outputs:

Python 3.9.16 (15ea147a6c65, May 30 2023, 13:28:11)
[PyPy 7.3.11 with GCC 13.1.1 20230511 (Red Hat 13.1.1-2)] on linux
Type "help", "copyright", "credits" or "license()" for more information.
>>> 
===== RESTART: /home/user/Downloads/OneDrive_1_20-07-2023/volume_measure.py ====
For lithium_density.cube the calculated average radius is 0.97 angstroms
For potassium_density.cube the calculated average radius is 1.74 angstroms
For rubidium_density.cube the calculated average radius is 1.91 angstroms
For sodium_density.cube the calculated average radius is 1.32 angstroms
>>> 

The results I got were very close to those in the original work:

Ion Results from Rahm et al. (Å) My results (Å)
Li+ 0.98 0.97
Na+ 1.33 1.32
K+ 1.75 1.74
Rb+ 1.91 1.91

On how these results compare with values derived from empirical data, they don't discuss the case for monopositive cations as in depth as is done with the neutral atoms. For the neutral atoms, they plot the difference between their calculated radii against crystallographic radii, as shown bellow. On average their radii are smaller by 0.01 Å, and the standard deviation corresponds to 0.211Å, with most alkali metals acting as outliers:

enter image description here

The determination of radii also depend on other factors, like coordination numbers. Regarding the ionic species, bellow I put some values compiled by Shannon(5) from empirical data, for comparison:

Ion coordination number (CN) crystal radius (CR) effective ionic radius (IR)
Li+ IV 0.730 0.590
Li+ VI 0.90 0.76
Na+ IV 1.13 0.99
Na+ V 1.14 1.00
Na+ VI 1.16 1.02
Na+ VII 1.26 1.12
Na+ VIII 1.32 1.18
Na+ XII 1.53 1.39
K+ IV 1.51 1.37
K+ VI 1.52 1.38
K+ VII 1.60 1.46
K+ VIII 1.65 1.51
K+ IX 1.69 1.55
K+ X 1.73 1.59
K+ XII 1.78 1.64
Rb+ VI 1.66 1.52
Rb+ VII 1.70 1.56
Rb+ VIII 1.75 1.61
Rb+ X 1.80 1.66
Rb+ XI 1.83 1.69
Rb+ XII 1.86 1.72
Rb+ XIV 1.97 1.83

References:

(1). Rahm, Martin, et al. “Atomic and Ionic Radii of Elements 1 –96.” Chemistry – A European Journal, vol. 22, no. 41, Oct. 2016, pp. 14625–32. DOI.org (Crossref), https://doi.org/10.1002/chem.201602949.

(2). Aprà, E., et al. “NWChem: Past, Present, and Future.” The Journal of Chemical Physics, vol. 152, no. 18, May 2020, p. 184102. DOI.org (Crossref), https://doi.org/10.1063/5.0004997.

(3). Lehtola, Susi, and Antti J. Karttunen. “Free and Open Source Software for Computational Chemistry Education.” WIREs Computational Molecular Science, vol. 12, no. 5, Sept. 2022. DOI.org (Crossref), https://doi.org/10.1002/wcms.1610.

(4). Tang, W., et al. “A Grid-Based Bader Analysis Algorithm without Lattice Bias.” Journal of Physics: Condensed Matter, vol. 21, no. 8, Feb. 2009, p. 084204. DOI.org (Crossref), https://doi.org/10.1088/0953-8984/21/8/084204.

(5). Shannon, R. D. “Revised Effective Ionic Radii and Systematic Studies of Interatomic Distances in Halides and Chalcogenides.” Acta Crystallographica Section A: Crystal Physics, Diffraction, Theoretical and General Crystallography, vol. 32, no. 5, Sept. 1976, pp. 751–67. scripts.iucr.org, https://doi.org/10.1107/S0567739476001551.

$\endgroup$
1
  • 1
    $\begingroup$ I gave my +1 long ago, but I want to say that I'm very surprised that this only has 2 upvotes so far! It was a very high-quality and high-effort answer! $\endgroup$ Sep 30, 2023 at 21:00

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .