Article pubs.acs.org/JCTC
Optimizing Vibrational Coordinates To Modulate Intermode Coupling Paul M. Zimmerman*,† and Peter Smereka‡ †
Department of Chemistry and ‡Department of Mathematics, University of Michigan, Ann Arbor, Michigan 48109, United States S Supporting Information *
ABSTRACT: The choice of coordinate system strongly affects the convergence properties of vibrational structure computations. Two methods for efficient generation of improved vibrational coordinates are presented and justified by analysis of a model anharmonic two-mode Hessian and numerical computations on polyatomic molecules. To produce optimal coordinates, metrics which quantify off-diagonal couplings over a grid of Hessian matrices are minimized through unitary rotations of the vibrational basis. The first proposed metric minimizes the total squared off-diagonal coupling, and the second minimizes the total squared change in off-diagonal coupling. In this procedure certain anharmonic modes tend to localize, for example X−H stretches. The proposed methods do not rely on prior fitting of the potential energy, vibrational structure computations, or localization metrics, so they are unique from previous vibrational coordinate generation algorithms and are generally applicable to polyatomic molecules. Fitting the potential to the approximate n-mode representation in the optimized bases for all-trans polyenes shows that off-diagonal anharmonic couplings are substantially reduced by the new choices of coordinate system. Convergence of vibrational energies is examined in detail for ethylene, and it is shown that coupling-optimized modes converge in vibrational configuration interaction computations to within 1 cm−1 using only 3-mode couplings, where normal modes require 4-mode couplings for convergence. Comparison of the vibrational configuration interaction convergence with respect to excitation level for the two proposed metrics shows that minimization of the total off-diagonal coupling is most effective for low-cost vibrational structure computations.
■
INTRODUCTION Vibrations are fundamental motions of all molecules and materials and, thus, are intimately related to structure and bonding.1−6 To probe these properties, vibrational spectroscopy provides a highly sensitive measurement tool that is especially powerful when used in cooperation with computational vibrational structure theories. These theories, including vibrational self-consistent field7−11 (VSCF), configuration interaction12,13 (VCI), coupled cluster,14 and perturbation theory,15−22 have been developed to provide characterization of experimental signals that map onto detailed descriptions of molecular motion. Vibrational structure methodologies all rely on an underlying choice of vibrational coordinates for the solution of the Schrödinger equation and representation of the potential energy. For this reason, the choice of coordinate system has recently been highlighted as an important factor in the success of vibrational structure simulations.23−27 Toward such a goal and inspired by earlier works,28−41 especially that of Truhlar,23 Hirata et al.24 showed that vibrational coordinates could be optimized42−44 to provide more rapid convergence and lower variational energies than normal modes. This was done by iteratively performing VSCF and VCI computations, and mixing vibrational modes into a new set that minimized the energy of the targeted vibrational state. Interestingly, the © 2016 American Chemical Society
procedure tended to localize X-H stretches, where the anharmonic coupling within a single localized mode outweighed the effect of coupling to other modes.45−49 In another insightful work, Steele et al.27 showed that localized modes have advantages over normal modes in terms of spatial decay of coupling elements. Localization was shown to allow compact representations of vibrational coupling. This procedure did not require repeated vibrational structure calculations, making it tractable for larger systems. These studies showed that the local nature of chemical interactions50−53 is an important effect that should influence the choice of coordinate system. In Hirata’s study24 not all modes became localized, suggesting maximized locality does not necessarily provide an optimal vibrational coordinate system.54 Such considerations allow us to frame two key unanswered questions in vibrational mode representations: 1. What set(s) of vibrational coordinates produce minimal intermode harmonic and/or anharmonic coupling? 2. To what extent is locality an intrinsic property of the vibrational potential energy? To provide insight into these questions, we introduce direct minimization of intermode coupling as a means to create Received: December 9, 2015 Published: February 25, 2016 1883
DOI: 10.1021/acs.jctc.5b01168 J. Chem. Theory Comput. 2016, 12, 1883−1891
Article
Journal of Chemical Theory and Computation 1
improved vibrational coordinates. Two methods are presented which minimize coupling without use of explicit localization metrics, but still create localized vibrational modes. These methods are unique from previous schemes because they avoid too much locality, which can create large harmonic couplings, as well as repeated vibrational structure computations,24 which have high computational cost. It will be shown that direct minimization of off-diagonal coupling produces coordinate systems that can be fit to high accuracy with simplified potential energy representations, compared to normal modes. These coordinate systems can be constructed without prior potential fitting, so they can be created directly from ab initio Hessians. VCI computations demonstrate that optimized coordinate systems can allow more rapid convergence of vibrational energies than normal modes (in agreement with the work of Hirata24). The new procedures use unitary transformations of normal mode coordinates and operate without prior construction of a vibrational Hamiltonian. These properties allow the proposed techniques to be highly applicable to vibrational analysis of polyatomic molecules. The outline of this article is as follows. First, a model vibrational Hessian is considered, where coupling can be analytically quantified as a function of coordinate choice. This analysis will show that locality is sometimes preferred, and other times detrimental. Next, algorithms for minimizing anharmonic coupling in general polyatomic systems are demonstrated. These methods are applied to the shorter polyene series of molecules where the anharmonic couplings are analyzed. Finally, VCI computations demonstrate the convergence behavior of vibrational signatures when using these optimized coordinates.
F (θ ) =
(4)
with α=
β=
γ=
2 (axbx − ayby + 3(ωx 2 − ωy 2)δ) 3
and
F(θ ) = α − β cos 4θ α=
1 2 ζ + 2δ 2 3
β=
1 2 ζ − 2δ 2 3
When minimizing F, the magnitude of anharmonicity, ζ, competes with the harmonic coupling, δ. In the limit of ζ ≫ δ, the preferred θ = 0, and the modes that minimize F will be fully localized. For this case, any delocalization increases F because mixing anharmonic modes results in increased B122(x, y). For ζ ≪ δ, θ = π/4, and delocalized (normal) modes are preferred. Conversely to the strongly anharmonic case, localizing the vibrations increases the total coupling as measured by F. If we further specify that δ = 0 and ζ = 0, there is no anharmonicity and any two vectors (any choice of θ) are eigenvectors of A. In this case, choosing eigenvectors
(1)
We can define a basis for the vibrations by mixing the local coordinates, x and y, over angle θ using the unitary transformation
(2)
⎛ ⎜ ⎛1 ⎞ ⎛0⎞ ⎜ ⎜ ⎟ , ⎜ ⎟ or ⎝0⎠ ⎝1 ⎠ ⎜ ⎜ ⎝
In the new basis, A becomes B(θ ) = UTAU
1 ((ax − by)2 + (bx − ay)2 + 3(ωx 2 − ωy 2)2 6 − 4(bx 2 + by 2) − 12δ 2)
⎛ ω 2 δ ⎞ ⎛ ζx 0 ⎞ x ⎟+⎜ ⎟⎟ A(x , y) = ⎜⎜ ⎜ 2⎟ ⎝ δ ωx ⎠ ⎝ 0 ζy ⎠
THEORY Model Hessian Anharmonicity Analysis. We seek to show the relationship between locality and coupling in the Hessian matrix. First, while Hessian matrices are generally functions of all coordinates involved, ab initio computations generally only construct the Hessian by evaluation of the matrix at specific points. Our model Hessian is chosen to be a function of local coordinates, x and y, that can be expressed as
⎛ cos θ −sin θ ⎞ U=⎜ ⎟ ⎝ sin θ cos θ ⎠
1 ((ax − by)2 + (bx − ay)2 + 3(ωx 2 − ωy 2)2 6 + 4(bx 2 + by 2) + 12δ 2)
F(θ) therefore represents a measure of anharmonic coupling. This analysis shows that achieving minimal coupling between the modes (over a span of x and y) generally depends on all Hessian parameters. In general, modes that minimize F may be anywhere between fully delocalized and fully localized. For instance, fully decoupled, localized (θ = 0) anharmonic oscillators and fully harmonic oscillators (θ varies) both have minimal values of F(θ) = 0. An illustrative example is that of two degenerate, anharmonic oscillators with constant coupling such that ωx = ωy, ax = ay = ζ, and bx = by = 0. In this case, A becomes
■
⎛ ω 2 δ ⎞ ⎛a x + b y b x + b y ⎞ x x y x ⎟+⎜ x ⎟ A (x , y ) = ⎜ ⎜ ⎟ 2 ⎜ δ ω ⎟ b x b y b x a y + + y y y ⎠ ⎝ x y ⎠ ⎝
1
∫−1 ∫−1 B122(θ) dx dy = α − β cos 4θ − γ sin 4θ
(3)
At x = y = 0, the vibrational coordinate system consists of normal modes if θ is chosen such that B is diagonal. B is not, however, generally diagonal for x ≠ 0, y ≠ 0. We therefore introduce a function, F(θ), which measures the off-diagonal coupling integrated over the neighborhood of vibrational motion. In the new basis B,
1 ⎞ ⎛ 1 ⎟ ⎜− 2 2⎟ ⎜ , 1 ⎟ ⎜ 1 ⎟ ⎜ 2⎠ ⎝ 2
⎞ ⎟ ⎟ ⎟ ⎟ ⎠
corresponding to the fully localized and fully delocalized representations, respectively, makes no difference. Instead of minimizing the sum-squared off-diagonal coupling over space through F(θ), one could minimize the off-diagonal components of the Hessian’s total derivative. This will lead to a 1884
DOI: 10.1021/acs.jctc.5b01168 J. Chem. Theory Comput. 2016, 12, 1883−1891
Article
Journal of Chemical Theory and Computation
atomwise, projective approach analogous to Pipek−Mezey orbital localization55 maximizes
metric designed to minimize the anharmonic off-diagonal coupling. Given the total derivative, dA =
⎛bx by ⎞ ⎛ ax bx ⎞ δA δA ⎟ ⎟d x + ⎜ dx + dy = ⎜⎜ ⎟ ⎜ b a ⎟d y δx δy ⎝ bx by ⎠ ⎝ y y⎠
3N − 6 N
L=
i
⎤ 2 ⎡ T δA ⎤ 2 ⎡ T δA F ′(θ ) = ⎢U (θ ) U (θ )⎥ + ⎢U (θ ) U (θ )⎥ ⎣ ⎦12 δx δy ⎣ ⎦12
3N − 6
gives
L′ =
α′ =
1 [(ax − by)2 + (bx − ay)2 + 4(bx 2 + by 2)] 8
β′ =
1 ((ax − by)2 + (bx − ay)2 − 4(bx 2 + by 2)) 8
γ′ =
1 (axbx − ayby) 2
A
(V Aiα)2 )2 (6)
α=x ,y,z
in order to minimize the number of atoms involved in a vibrational mode. Here, ViAα refers to a Cartesian coordinate α of atom A involved in vibration i. A related technique, similar to Boys orbital localization,56,57 maximizes the distance between mode centers by maximizing
and summing over the squared off-diagonals of the rotated partial derivatives
F ′(θ ) = α′ − β′cos 4θ − γ ′sin 4θ
∑ ∑( ∑
N
∑ (∑ ∑ i
(5)
(V Aiα)2 RA)2
A α=x ,y,z
(7)
where RA is the position of atom A. Improved coordinate systems are likely to be generated by partially localizing modes, where locality is desirable for distant vibrations, but strong coupling results when L or L′ is maximized. This large coupling occurs because neither the Pipek−Mezey or Boys methods explicitly take into consideration any information about coupling.58 For this reason, locality metrics alone cannot be expected to yield low coupling vibrations. Given the results of the Theory section, we seek to create an improved set of vibrational coordinates without explicitly invoking localization metrics. The first set of vibrations is formed by creating a grid of M + 1 Hessians near the lowest energy geometry, and minimizing
No integration over x or y is performed because F′ is not a function of these variables. Interestingly, the minimizer of F′ closely resembles that of F, except the functional dependences of ω and δ have disappeared. By construction, minimizing F′ rotates the vibrational coordinates to maximize diagonal anharmonicity. F′ can always be set equal to zero if bx = by = 0. In this case θ = 0, the optimal vibrations are localized, and all anharmonic effects will appear on the Hessian diagonal. For nonzero values of bx and by, local modes are preferred unless there exists significant anharmonic off-diagonal coupling, regardless of the magnitude of harmonic coupling, δ. The two proposed coupling measures, F and F′, form plausible metrics based on energetic coupling which indicate whether or notand how muchto localize vibrations. Generally, localization is preferred when diagonal anharmonicity dominates over coupling (F) or when mode−mode anharmonicity is small (F′). Neither metric suggests all modes should be localized. Specific modes types, however, are likely to become localized: X−H stretches, for instance, likely fall into the categories of ζ ≫ δ and bx ≈ by ≈ 0. These results are independent of spatial locality metrics, suggesting locality is an emergent feature of the Hessian. An algorithm will now be presented to generalize this analysis to the many-mode case, and subsequent examples of optimal coupling modes will demonstrate their corresponding degree of locality. Algorithms for Local and Optimal Vibrations. These considerations allow us to formulate an algorithm to create improved vibrational coordinates. Consider the matrix, V′, composed of the 3N − 6 nonzero eigenvectors of the Hessian matrix, H, represented in mass-weighted Cartesian coordinates evaluated at a minimum energy geometry. Taking unitary transformations over the vector of all angles, θ, which rotate pairs of eigenvectors, a new basis, V, can be created
M
E=
3N − 6
∑ ∑ m=0
(H̃ ijm)2 (8)
i
