Calculation of nuclear magnetic resonance shieldings using frozen-density embedding

We have extended the frozen-density embedding (FDE) scheme within density-functional theory [T. A. Wesolowski and A. Warshel, J. Phys. Chem. 97, 8050 (1993)] to include external magnetic fields and applied this extension to the nonrelativistic calculation of nuclear magnetic resonance (NMR) shieldings. This leads to a formulation in which the electron density and the induced current are calculated separately for the individual subsystems. If the current dependence of the exchange-correlation functional and of the nonadditive kinetic-energy functional are neglected, the induced currents in the subsystems are not coupled and each of them can be determined without knowledge of the induced current in the other subsystem. This allows the calculation of the NMR shielding as a sum of contributions of the individual subsystems. As a test application, we have calculated the solvent shifts of the nitrogen shielding of acetonitrile for different solvents using small geometry-optimized clusters consisting of acetonitrile and one solvent molecule. By comparing to the solvent shifts obtained from supermolecular calculations we assess the accuracy of the solvent shifts obtained from FDE calculations. We find a good agreement between supermolecular and FDE calculations for different solvents. In most cases it is possible to neglect the contribution of the induced current in the solvent subsystem to the NMR shielding, but it has to be considered for aromatic solvents. We demonstrate that FDE can describe the effect of induced currents in the environment accurately.     

NMR solvent shifts of acetonitrile from frozen density embedding calculations

We present a density functional theory (DFT) study of solvent effects on nuclear magnetic shielding parameters. As a test example we have focused on the sensitive nitrogen shift of acetonitrile immersed in a selected set of solvents, namely water, chloroform, and cyclohexane. To include the effect of the solvent environment in an accurate and efficient manner, we employed the frozen-density embedding (FDE) scheme. We have included up to 500 solvent molecules in the NMR computations and obtained the cluster geometries from a large set of conformations generated with molecular dynamics. For small solute-solvent clusters comparison of the FDE results with conventional supermolecular DFT calculations shows close agreement. For the large solute-solvent clusters the solvent shift values are compared with experimental data and with values obtained using continuum solvent models. For the water --> cyclohexane shift the obtained value is in very good agreement with experiments. For the water --> chloroform NMR solvent shift the classical force field used in the molecular dynamics simulations is found to introduce an error. This error can be largely avoided by using geometries taken from Car-Parrinello molecular dynamics simulations.     

Is it worthwhile to go beyond the local-density approximation in subsystem density functional theory?

Frozen density embedding (FDE) theory is one of the major techniques aiming to bring modeling of extended chemical systems into the realm of high accuracy calculations. To improve its accuracy it is of interest to develop kinetic energy density functional approximations specifically for FDE applications. In the study reported here we focused on optimizing parameters of a generalized gradient approximation-like kinetic energy functional with the purpose of better describing electron excitation energies. We found that our optimized parametrizations, named excPBE and excPBE-3 (as these are derived from a Perdew-Burke-Ernzerhof-like parametrization), could not yield improvements over available functionals when applied on a test set of systems designed to probe solvatochromic shifts. Moreover, as several different functionals yielded very similar errors to the simple local-density approximation (LDA), it is questionable whether it is worthwhile to go beyond the LDA in this context.     

Potential-energy surfaces of local excited states from subsystem- and selective Kohn-Sham-TDDFT

Calculating excited-state potential-energy surfaces for systems with a large number of close-lying excited states requires the identification of the relevant electronic transitions for several geometric structures. Time-dependent density functional theory (TDDFT) is very efficient in such calculations, but the assignment of local excited states of the active molecule can be difficult. We compare the results of the frozen-density embedding (FDE) method with those of standard Kohn-Sham density-functional theory (KS-DFT) and simpler QM/MM-type methods. The FDE results are found to be more accurate for the geometry dependence of excitation energies than classical models. We also discuss how selective iterative diagonalization schemes can be exploited to directly target specific excitations for different structures. Problems due to strongly interacting orbital transitions and possible solutions are discussed. Finally, we apply FDE and the selective KS-TDDFT to investigate the potential energy surface of a high-lying π → π^∗ excitation in a pyridine molecule approaching a silver cluster.     

Topological analysis of electron densities from Kohn-Sham and subsystem density functional theory

In this study, we compare the electron densities for a set of hydrogen-bonded complexes obtained with either conventional Kohn-Sham density functional theory (DFT) calculations or with the frozen-density embedding (FDE) method, which is a subsystem approach to DFT. For a detailed analysis of the differences between these two methods, we compare the topology of the electron densities obtained from Kohn-Sham DFT and FDE in terms of deformation densities, bond critical points, and the negative Laplacian of the electron density. Different kinetic-energy functionals as needed for the frozen-density embedding method are tested and compared to a purely electrostatic embedding. It is shown that FDE is able to reproduce the characteristics of the density in the bonding region even in systems such as the F-H-F(-) molecule, which contains one of the strongest hydrogen bonds. Basis functions on the frozen system are usually required to accurately reproduce the electron densities of supermolecular calculations. However, it is shown here that it is in general sufficient to provide just a few basis functions in the boundary region between the two subsystems so that the use of the full supermolecular basis set can be avoided. It also turns out that electron-density deformations upon bonding predicted by FDE lack directionality with currently available functionals for the nonadditive kinetic-energy contribution.     

Wave‐function frozen‐density embedding: Approximate analytical nuclear ground‐state gradients

We report the derivation of approximate analytical nuclear ground‐state uncoupled frozen density embedding (FDEu) gradients for the resolution of identity (RI) variant of the second‐order approximate coupled cluster singles and doubles (RICC2) as well as density functional theory (DFT), and an efficient implementation thereof in the KOALA program. In order to guarantee a computationally efficient treatment, those gradient terms are neglected which would require the exchange of orbital information. This approach allows for geometry optimizations of single molecules surrounded by numerous molecules with fixed nuclei at RICC2‐in‐RICC2, RICC2‐in‐DFT, and DFT‐in‐DFT FDE level of theory using a dispersion correction, required due to the DFT‐based treatment of the interaction in FDE theory. Accuracy and applicability are assessed by the example of two case studies: (a) the Watson‐Crick pair adenine‐thymine, for which the optimized structures exhibit a maximum error of about 0.08 Å for our best scheme compared to supermolecular reference calculations, (b) carbon monoxide on a magnesium oxide surface model, for which the error amount up to 0.1 Å for our best scheme. Efficiency is demonstrated by successively including environment molecules and comparing to an optimized conventional supermolecular implementation, showing that the method is able to outperform conventional RICC2 schemes already with a rather small number of environment molecules, gaining significant speed up in computation time. © 2016 Wiley Periodicals, Inc.     

多参考态二级微扰方法的改进

本文提出一个新的扩大模型空间的方案用于改进多参考态二级微扰理论(MRPT2)计算. 新方案保持了原方案中扩大模型空间之前的简单程序结构, 理论上完全可以避免势能面计算中入侵态的出现, 并在一系列比较计算中得到证实. 新MRPT2程序是研究分子激发态和电子光谱的有用工具.     

A flexible implementation of frozen-density embedding for use in multilevel simulations

A new implementation of frozen-density embedding (FDE) in the Amsterdam Density Functional (ADF) program package is presented. FDE is based on a subsystem formulation of density-functional theory (DFT), in which a large system is assembled from an arbitrary number of subsystems, which are coupled by an effective embedding potential. The new implementation allows both an optimization of all subsystems as a linear-scaling alternative to a conventional DFT treatment, the calculation of one active fragment in the presence of a frozen environment, and intermediate setups, in which individual subsystems are fully optimized, partially optimized, or completely frozen. It is shown how this flexible setup can facilitate the application of FDE in multilevel simulations.     

Atoms-in-molecules analysis for planewave DFT calculations--a numerical approach on a successively interpolated charge density grid

We used a successive charge interpolation scheme and Ridders method for differentiation, to acquire accurate charge densities and their higher derivatives in electronic structure calculations. This enables us to search bond critical points using arbitrary charge density grids. We applied the planewave-DFT code, VASP, to generate the charge density of selected benchmark molecules. The properties of bond critical points are in good agreement with those obtained by complementary implementations. We validated our GRID implementation by performing electronic structure calculations using the Gaussian 03 program package and various tools for analysis of the charge density provided by the AIMPAC package. In particular, we carefully investigate the influence of effective core potentials on the location of bond critical points, especially for a short chemical bond, which is crucial in the present pseudopotential-based planewave DFT calculations. We expect our generic implementation will not only be useful for the analysis of chemical bonding in molecules, but will particularly provide a microscopic understanding of extended systems including periodic boundary conditions.     

A technique for orbital exponent optimization in ab initio HF?SCF?LCAO?MO calculations

A technique for Slater orbital exponent optimization in an HF? SCF? LCAO? MO calculation is proposed in which orbital exponent variation is incorporated into the SCF scheme. This is accomplished by rewriting Slater's rules so that the shielding terms depend on the molecular charge distribution through the elements of the population matrix. The SCF scheme then includes a calculation of a new set of orbital exponents from the coefficients of self-consistent molecular orbitals obtained from the previous set of exponents. The process is iterated until the energy attains its lowest value. The technique is illustrated by minimal basis calculations on LiH, BH, and HF. Near optimization is obtained with considerably less effort than is necessary for other reported techniques. Aside from interesting properties, the technique can be important for extended basis calculations where exponent optimization is a difficult task.     

A density-functional approach to polarizable models: a Kim-Gordon response density interaction potential for molecular simulations

A combined linear-response-frozen electron-density model has been implemented in a molecular-dynamics scheme derived from an extended Lagrangian formalism. This approach is based on a partition of the electronic charge distribution into a frozen region described by Kim-Gordon theory [J. Chem. Phys. 56, 3122 (1972); J. Chem. Phys. 60, 1842 (1974)] and a response contribution determined by the instantaneous ionic configuration of the system. The method is free from empirical pair potentials and the parametrization protocol involves only calculations on properly chosen subsystems. We apply this method to a series of alkali halides in different physical phases and are able to reproduce experimental structural and thermodynamic properties with an accuracy comparable to Kohn-Sham density-functional calculations.     

A modified multi-reference second order perturbation theory

A new scheme with extended model space is proposed to improve the calculation of multi-reference second order perturbation theory (MRPT2). The new scheme preserves the concise code structure of the original program, and avoids intruder states in constructions of the potential energy surface, which is confirmed by a series of comparable calculations. The new MRPT2 program is an available tool for the research of molecular excited states and electronic spectrum.     

Basis set superposition error in N-body clusters

We present a new scheme for calculating the basis set superposition error (BSSE) in N-body clusters. It is based on the assumption that each n-body term can be expressed as a sum of only two-body contributions. The conventional Boys–Bernardi method can be used thus for calculating BSSE-corrected energy terms. The scheme is illustrated by some calculations on the hydrogen fluoride trimers and tetramers. The results are compared to the ones obtained with the site–site function counterpoise (SSFC), the hierarchical (Valiron–Mayer) function counterpoise (VMFC), the pairwise additive function counterpoise (PAFC), and the successive reaction counterpoise (SRCP) schemes.     

Computation of methodology-independent single-ion solvation properties from molecular simulations. III. Correction terms for the solvation free energies, enthalpies, entropies, heat capacities, volumes, compressibilities, and expansivities of solvated ions

The raw single-ion solvation free energies computed from atomistic (explicit-solvent) simulations are extremely sensitive to the boundary conditions (finite or periodic system, system or box size) and treatment of electrostatic interactions (Coulombic, lattice-sum, or cutoff-based) used during these simulations. However, as shown by Kastenholz and Hu?nenberger [J. Chem. Phys. 124, 224501 (2006)], correction terms can be derived for the effects of: (A) an incorrect solvent polarization around the ion and an incomplete or/and inexact interaction of the ion with the polarized solvent due to the use of an approximate (not strictly Coulombic) electrostatic scheme; (B) the finite-size or artificial periodicity of the simulated system; (C) an improper summation scheme to evaluate the potential at the ion site, and the possible presence of a polarized air-liquid interface or of a constraint of vanishing average electrostatic potential in the simulated system; and (D) an inaccurate dielectric permittivity of the employed solvent model. Comparison with standard experimental data also requires the inclusion of appropriate cavity-formation and standard-state correction terms. In the present study, this correction scheme is extended by: (i) providing simple approximate analytical expressions (empirically-fitted) for the correction terms that were evaluated numerically in the above scheme (continuum-electrostatics calculations); (ii) providing correction terms for derivative thermodynamic single-ion solvation properties (and corresponding partial molar variables in solution), namely, the enthalpy, entropy, isobaric heat capacity, volume, isothermal compressibility, and isobaric expansivity (including appropriate standard-state correction terms). The ability of the correction scheme to produce methodology-independent single-ion solvation free energies based on atomistic simulations is tested in the case of Na(+) hydration, and the nature and magnitude of the correction terms for derivative thermodynamic properties is assessed numerically.     

An auto-coherent variational scheme for Yukawa potential

The variational principle is used to obtain solutions to Schrödinger's equation for a particle in the Yukawa potential. A Laguerre basis set extended by an extra function is employed in the calculations. A special parameter used in the extra function and its relation with the systems energy results in utilizing an auto-coherent scheme. Considerable improvement seems to be achieved especially in the critical region where the screening parameter approaches its threshold value.     

Approaching the bulk limit with finite cluster calculations using local increments: the case of LiH

Finite-cluster calculations employing high-level wavefunction-based ab initio methods and extended atomic-orbital basis sets are used to determine local energy increments for bulk LiH. It is shown that these increments can be converged with respect to cluster size and point-charge embedding so as to yield bulk cohesive energies with an accuracy of better than 1 mE(h), both at the Hartree-Fock and at correlated levels. Instrumental for the efficiency of the scheme is the introduction of non-orthogonal orbitals, at an intermediate stage.     

SCF calculations of the interactions of alkali and halide ions with the mercury surface

From extensive ab initio calculations on the interactions between mercury clusters and alkali and halide ions we have derived analytical pair-potential functions for the interaction between the ion and an extended mercury (111) surface. A novel correction scheme is proposed in order to reduce the shortcomings of cluster model. A preferred adsorption above the twofold bridge site was found for Li⁺ and Na⁺ and above the threefold hollow site for all other ions. The ab initio results have been fitted to analytical functions that can be used in computer simulations.     

A quartic force field coordinate substitution scheme using hyperbolic sine coordinates

Quartic force fields (QFF) are currently the most cost‐effective method for the approximation of potential energy surfaces for the calculation of anharmonic vibrational energies. It is known, although, that its performance can be less than satisfactory due to limitations related to slow convergence of the series. In this article, we present a coordinate substitution scheme using a combination of Morse and sinh coordinates, well adapted for its use with cartesian normal coordinates. We derive expressions for analytical integrals for use in VSCF and VCI calculations and show that the simultaneous substitution of symmetric and antisymmetric normal coordinates by Morse and sinh coordinates, respectively, significantly improves the vibrational transition frequencies for these modes in a well‐balanced fashion. The accuracy of this substitution scheme is demonstrated by comparing one and two‐dimensional sections of substituted and unsubstituted QFF with ab initio potential energy grids, as well as with vibrational energy calculations using as test cases two well‐studied benchmark molecules: water and formaldehyde. We conclude that the coordinate substitution scheme presented constitutes a very attractive alternative to simple QFFs in the context of cartesian normal coordinates.