Novel equations to predict vibrational spectroscopic and electrodynamics properties of molecules

Mathematical relations linking electric and magnetic field-dependent physical properties of molecules have been unveiled. These relations are analogous to Maxwell relations in thermodynamics and are derived from mixed third-order partial derivatives of every alternative Legendre representation of the energy of molecules with respect to the electric or magnetic field and normal coordinates. Some of these novel physical relationships have practical applications in the low computational cost calculation of parameters commonly used in vibrational spectroscopy like the Stark and Zeeman Tuning Rates. Furthermore, other equalities have shown connections and alternative ways of computing physical properties used in electrodynamics as permanent dipolar moments and polarizabilities.


Introduction
In thermodynamics [1], the relation that gives the internal energy (U) as a function of the extensive parameters like entropy (S), volume (V ), and number of particles (N) is known as the fundamental relation of a thermodynamical system (equation (1)). If the fundamental equation of a system is known, every thermodynamic attribute is completely and precisely determined. The derivative of the internal energy as a function of the extensive parameters gives rise to the intensive parameters, which are important physical properties such as the temperature (T), pressure (P) and electrochemical potential (μ). Therefore, the differential equation of the fundamental relation of the internal energy is expressed as a function of the intensive parameters and the differentials of the extensive parameters (equation (2)).
A mathematical formalism known as Legendre transformations [1] is used in thermodynamics to obtain other fundamental equations from the fundamental relation of the internal energy. In such reformulations, the intensive parameters replace extensive parameters as mathematically independent variables. All these fundamental equations are known as thermodynamic potentials. Gibbs free energy (G), Enthalpy (H), Helmholtz free energy (F) and Grand Canonical potential (C) are some of the most commonly used thermodynamic potentials and their corresponding differential equations (3).   A formal aspect of thermodynamics is the set of Maxwell relations [1]. These relations arise from the equality of the mixed partial second derivatives of the fundamental relation expressed in any of the various possible alternative representations. One example of a Maxwell relation is shown in equation (4).
The second derivatives of fundamental equations are also descriptive of material properties of direct physical interest like the coefficient of thermal expansion (α), the isothermal or adiabatic compressibility [2], and the molar heat capacity at constant volume or pressure [1]. Using Jacobian operations and Maxwell relations, the isothermal compressibility (κ T ) can be related to adiabatic compressibility (κ S ), and the heat capacity at constant volume (c V ) to the heat capacity at constant pressure (c P ) (equation (5)) [1]. An analogous mathematical approach used in thermodynamics to develop equation (5) will be used later on. In this work, all these equivalent mathematical relations will be derived for a molecular system where the fundamental equation is the potential energy of the molecule and the independent variables are the nuclear normal coordinates and the elements of an external uniform electric (or magnetic) field vector. Similar equalities are formulated from mixed third-order partial derivatives of the energy and other Legendre representations. The obtained relations carry practical uses in alternative and faster ways of computing properties used in vibrational spectroscopy and electrodynamics such as the Stark and Zeeman Tuning Rates, permanent dipolar moments and polarizabilities. In most cases, these parameters can be obtained by means of these novel relations without the application of an external electric or magnetic field during their computation, thus reducing the cost of the calculation. Furthermore, relations obtained by Jacobian operations show other connections between these electrodynamics parameters.
The vibrational Stark effect [3][4][5][6][7][8][9][10][11] is the shifting of vibrational frequencies due to the presence of an external electric field. These frequencies can be measured using Infrared or Raman Spectroscopy [12,13]. The dependence of the wavenumber (ñ ) of a particular mode m with an electric field strength (F) can be described using a Taylor expansion as: The vibrational Zeeman effect [14] is the same phenomenon caused by an external magnetic field. The Stark Tuning Rate [15][16][17] (STR) and the Zeeman Tuning rate [18] (ZTR) are defined as the derivative of the frequency with respect to the electric field (F) or magnetic field (B) at the limit of zero field. These parameters are broadly used for estimating the electric or magnetic field magnitude experimentally by monitoring the shift in the frequency of molecular probes. It is a common practice to calculate theoretically the STR and ZTR due to the experimental difficulties in measuring them.
In this paper, a theoretical development of these new relationships is presented alongside some particular practical applications. Subsequently, as a proof of concept, the Stark Tuning Rate is numerically computed employing one of these novel equations for different chemical compounds typically used as molecular probes. Moreover, for the purpose of comparison, the STR is also calculated in a traditional manner with the computational approach normally adopted in the literature.

Theoretical development
2.1. Mixed second and third-order partial derivative equalities The potential energy (E) of a molecular system with fixed nuclear positions is a function of all atom coordinates. When a molecule is under the presence of an external uniform electric field (EF), the energy also depends on the vector element values of this new variable (F x , F y , F z ). For convenience, the normal coordinates of a molecule calculated in the absence electric field (Q 1 , K, Q n ; where n is the total number of normal modes) are used instead of nuclear Cartesian coordinates [12,13]. Therefore, the potential energy of a molecule under the effect of an EF can be expressed as a function of the vector elements of the EF and the normal coordinates of the molecule as shown in equation (7).
The total differential of E is given by: F F Q  Q   n   , , , ,  , , , ,  , , , ,   1  , , , , ,   1 , , , , , y z n x z n x y n x y z n x y z n 1 1 1 It can be noticed that the first derivative of the energy with respect to the EF vector elements yields the vector elements of the electric dipolar moment ( , , x y z m m m ) as [13]: , , , , , , , , In addition, the first derivative of the energy with respect to the normal mode coordinates produces the normal forces ( f 1 , K, f n ) as [13]: Now, equation (8) can be summarised as: The second derivatives of the energy yield other known parameters such as the polarizability tensor elements ( i j , a where i, j, k=x, y, z) and the normal force constants (λ m where m n 1   ) [12,13]: A Legendre transformation [1] can be applied to the energy to define other fundamental expressions of the energy, so-called energy potentials, depending on f 1 , K, f n and μ x , μ y , μ z in addition to Q 1 , K, Q n and F x , F y , F z 1 R F F F f f , , , , , 14 x y z n The Legendre transformation implies that: where i, j, k=x, y, z and m p n 1 ,   , and that: When the molecular geometry is optimised, the forces are zero, then R equals E. When the electric field is absent, T equals E. If both conditions are satisfied, Y equals E. Thus, the same applies to the total differentials of such energy potentials.
The differentials of the new fundamental expressions of the energy (R, T and Y) can be obtained similarly to the one for E: Based on the symmetry of partial derivatives, the equality of mixed second-order partial derivatives of E, R, T, and Y produces the relations shown in section A.1. From all equations only the following two will be used in further derivations.
The second and third derivatives of E and T in sections A.1 and A.2 are defined as restrained since the molecular normal coordinates are constant during derivation, while the derivatives of R and Y are defined as unrestrained since the nuclear coordinates are variable and the forces are constant. In the following text, an upperscore will denote unrestrained parameters.
Exact equivalent equations relating the magnetic field (B), the magnetic dipole moment (m ) and the magnetic polarizability (a  ) instead of their electric counterparts are obtained for the case of vibrational Zeeman effects (equations not shown).

Vibrational Stark effect
The vibrational Stark effect provides to experimentalists the possibility of determining the local electric field of molecules at diverse environments by measuring the shift of their vibrational frequencies. As shown in equation (6), knowing in advance the lower order derivatives of the Taylor expansion provides the necessary information for precisely calculating the magnitude of such field. The first-order Stark Tuning Rate (STR or σ) is defined as the derivative of the wavenumber 2 of a particular mode m ( m ñ ) with respect to the electric field at F=0 [3]. The second derivative is known as second-order Stark Tuning Rate.

Restrained first-order Stark Tuning Rate
The restrained Stark Tuning Rate (rSTR or σ) has been previously used by Brewer et al [5]. The rSTR is defined as the first derivative of the wavenumber with respect to the field holding Q constant (the nuclear configuration coordinates are constant): A shortcut relation to determine the rSTR avoiding the application of an electric field using the mixed thirdorder derivative relation is stated in equation (29). Based on this equation, the derivative of the normal force 2 It can also be the frequency (ν).
constant λ m with respect to the electric field at F i (right term) at fixed nuclear positions (Q p constants, p " ) can be determined as the second derivative of the electrical dipolar moment over the i axis with respect to the normal coordinate m.
The equation (32) [12] relates m ñ with λ m : To obtain the rSTR, equations (32) and (31) can be used to express σ as a function of the normal force constant: Subsequently, by performing a replacement using equation (29), the desired relation (equation (34)) is derived.
shows how rSTR can be computed alternatively to equation (31). Starting from an optimised molecular geometry in vacuum, a frequency analysis can be performed to obtain all normal force constants, thus, their vibrational frequencies. For a selected frequency, the corresponding normal coordinate is obtained. A normal coordinate is a displacement vector for all atom in the Cartesian space. This vector is used to displace the coordinates of each atom by multiplying the vector by small displacement values and adding this values to the relaxed structure. For each distorted geometry, an electrical dipolar moment is calculated. Finally, the rSTR is obtained from the second derivative of the dipolar moment along the i direction versus the displacement values applied to the normal mode m. This should be done for i = x, y and z to obtain the rSTR vector and module. The electric dipolar moment of a molecule can be computed, by means of quantum mechanics, by applying the electric dipolar operator to the wavefunction of the molecule instead of using the equation (17). This way of computing the rSTR provides a computational alternative that does not require the application of an electric field, compared to the traditional method in which the normal force constants are numerically computed at different electric fields. The determination of rSTR by the traditional method requires fixing the nuclear positions while applying an electric field to the structure optimised in vacuum. This fact impies that the molecular geometry is not in an energetically relaxed state. This is a drawback of the traditional method since such approach is subject to computational problems during the calculation of the frequencies, given that, when an electric field is applied to a fixed unrelaxed molecular structure, the molecular forces are not zero, and the frequencies analysis may not converge.

Unrestrained first-order Stark Tuning Rate
The unrestrained equation (35) implies that all the normal mode forces are constants 3 .
In the particular case when all forces are zero (f p 0, p = " ), the molecular geometry is in a minimum (relaxed state) or maximum (transition state) of energy. Combining the equation (32) together with the equation (35), it leads to equation (36).
From the mixed third-order partial derivative relations, an equality for the derivative in equation (36) can be obtained.
By using the demonstration in section A.3 (equation (83)), it can be proven that equation (13) is equivalent to equation (37).
Based on equation (30) and replacing by the inverse of equation (37) to obtain: and finally: The combination of equations (36) and (39) yields the final relation for the alternative calculation of the unrestrained first-order Stark Tuning Rate (uSTR): It is demonstrated in section A.4 of the appendix (equation (86)) that equation (40) is equivalent to equation (41):

Connection between restrained and unrestrained first-order Stark Tuning Rates
Presently, it is feasible to relate the restrained Stark Tuning Rate (σ) with the unrestrained Stark Tuning Rate (s ) using equation (98) in section A.6 where l=f or Q, indistinctly: Now replacing equation (42) by equations (35) and (31) leads to equation (43): For the case l=f, the last derivative in equation (43) can be replaced by equation (27) and then using equation (88) it is proven that the equation (44) is valid for l=f or Q indistinctly.

Unrestrained second-order Stark Tuning Rate
There is also a relation that links the second-order Stark Tuning Rate [20] with the polarizability. The second derivative of equation (32) with respect to the electric field provides under unrestrained configuration (constant forces) the following: The second derivative of λ m is obtained by taking the derivative of equation (39): where: Using equation (46) in conjunction with equations (47) and (39) affords: (45) using equations (48) and (39) yields:

Replacing in equation
With relation (49) it is possible to calculate the unrestrained second coefficient of the Taylor expansion in equation (6) by determining the normal force constant as well as the second derivative of the permanent dipolar moment and the polarizability with respect to the normal force.

Restrained second-order Stark Tuning Rate
Another relation links the restrained second-order Stark Tuning Rate [20] with the polarizability. The second derivative of equation (32) with respect to the electric field at restrained nuclear configuration provides: The second derivative of λ m is obtained by taking the derivative of equation (29): Using equation (51) in conjunction with equation (52) affords: With relation (54) it is possible to calculate the restrained second coefficient of the Taylor expansion of equation (6) by determining the normal force constant as well as the second derivative of the permanent dipolar moment and the polarizability with respect to the normal force along the normal coordinate m.

Dipolar moments and polarizabilities 2.3.1. Restrained and unrestrained dipolar moments
So far, the permanent dipolar moment was defined as the derivative of the energy with respect to the electric field at constant normal coordinates (μ; equation (55)). In quantum chemistry, this restrained dipolar moment μ is usually computed by applying the dipolar moment operator to the wavefunction of a molecule with fixed nuclear geometry. Alternatively, an unrestrained dipolar moment can also be defined ( ; m equation (55)). Under this definition, when forces are constant, the nuclear configuration has no geometric restrain. At constant null forces, the system is in a relaxed state (minimum of energy) or a transition state (maximum of energy). The restrained dipolar moment implies that the nuclear configuration is fixed and the electronic configuration is unrestrained. While for the unrestrained dipolar moment, the nuclear and electronic configuration are unrestrained. The numerical calculation of m by using the equation (55) requires a geometry optimisation for each electric field value which increase the computational cost compared to μ.
It is feasible to link both dipolar moment definitions by combining equations (55) and (98) derived from Jacobian properties, in section A.6: where l=f or Q, indistinctly. By replacing For a relaxed molecular structure (at f m =0), the unrestrained and restrained dipolar moments are equal provide information regarding the molecular geometry distortion given by an electric field. However, as it was done for the derivation of equation (44), these derivatives can be replaced by their mixed second order derivative counterparts (equation (27)) and yield the final expression: With equation (57), it is possible to compute the unrestrained dipolar moment when the molecular structure is away from a minimum energetic state.

Restrained and unrestrained polarizabilities
Similarly to the permanent dipolar moment, the polarizability can be defined as restrained (α) (equation (58)) and unrestrained (ā ) (equation (59)) where, for α the nuclear configuration is restrained and the electronic configuration is unrestrained, while for ā the nuclear and electronic configurations are unrestrained.
By renaming the term α as the electronic polarizability and ā as the total polarizability, the equation (63) states that the total polarizability (ā ) is equivalent to the sum of the electronic (α) and the nuclear polarizabilities (ã ). The necessity for separation of the total polarizability in nuclear and electronic contributions has been previously discussed [21].
The total (unrestrained) polarizability can be computed numerically by means of its definition (equation (59)) or by computing the electronic and nuclear polarizability terms and equation (63). The electronic polarizability is a well-known property usually computed by default in several quantum chemistry programs [22][23][24]. The nuclear polarizability can be determined either by using the equation (62) or by using in conjunction equations (61) and (56). In this way, the determination of the nuclear polarizability term at null electric field can be done avoiding the application of an electric field.

Vibrational Zeeman effect
The vibrational Zeeman effect [14] is the shifting of vibrational frequencies due to the presence of an external magnetic field (B). Analogously to the STR, the first and second derivatives of the wavenumber with respect to the magnetic field yield the firstand second-order Zeeman Tuning Rates [18] (ZTR), respectively. The exact equivalent equations found above can be obtained relating the magnetic field (B), the magnetic dipole moment (m ) and the magnetic polarizability (a  ) instead of their electric counterparts: Restrained Zeeman Tuning Rate: Unrestrained second-order Zeeman Tuning Rate:

Proof of concept
As a proof of principle, in this section, the unrestrained Stark Tuning Rate (uSTR) was numerically computed for a set of small molecules in two ways as a matter of comparing the traditional method and the proposed method in this paper. The traditional way uses the equation by definition (equation (35)) and the newly developed method uses the singular equation (41).

Methods
Quantum calculations were performed using methodologies based on Density Functional Theory (DFT). The hybrid three parameters exchange functional from Becke combined with the correlation functional from Lee, Yang and Parr (B3LYP) was chosen to compute the electric field effects on the wavenumber modes [25][26][27]. All calculations were performed with Gaussian 03 [24] using the basis set TZV. The molecules selected for this study are commonly used molecular probes: CO, NO, HCN, water, methane and 4-Chloro-Benzonitrile. For the HCN and 4-Chloro-Benzonitrile, the STR was computed for the vibrational mode corresponding to C-N vibration. For water, the vibrational mode of 1567 cm −1 was selected. For methane, the vibrational mode with highest vibrational frequency was chosen. The electric field vector was pointed in the direction of the vibration of the selected mode. A geometry optimisation was performed for all structures in vacuum with an extremely tight converge criteria 4 . A frequency analysis for all optimised structures in vacuum was performed to obtain the normal mode vectors for all the molecules. Unrestrained first-order Stark Tuning Rate was then computed using the traditional method and the novel method presented in this work.

Traditional method for computing uSTR
From the frequency analysis, vibrational frequencies were obtained for all molecules at the electric field modules of F=±0.01, ±0.005, ±0.002, ±0.001, ±0.0005, ±0.0002, 0.0 au. Geometries were previously tightly optimised at each field. STR was numerically computed by fitting a 6th degree polynomial to the frequency versus electric field module tables and by taking the derivative at F=0 according to equation (35).

Results
The obtained uSTR values from both methods were compared for precision and computational cost. Table 1 shows that the values of the computed uSTR are close and well within the numerical error. The percentage differences 5 in the uSTR for the standard and the new method seem to be higher in polar molecules (H 2 O, HCN and 4-Chloro-Benzonitrile) and with greater number of atoms. One possible explanation is that at higher number of atoms the numerical errors in the determination of the derivatives increase. The results from the traditional method cannot be considered as exact, neither from the new method. Nevertheless, the standard method may be prompt to higher numerical errors due to during the application of the electric field and optimisation, the molecular geometry changes, allowing the molecule to slightly rotate. This difference could be reduced if the electric field vector was applied in the direction of the net dipolar moment. This may explain the reason for a higher difference in polar molecules considering that the higher the dipolar moment is, the greater the rotation. The frequency analysis could be also another source of imprecision. However, the relation found in this paper that connects both methods are mathematically equivalent and further analysis should be done to identify the source of inaccuracy in each method and to determine which method is more precise. The new method may yield more accurate uSTR values as it does not require the application of an electric field and requires just a one time frequency analysis computation. This reduces the mentioned sources of error.
A comparison of the computation time is shown in table 2. It can be seen that for bimolecular species the time spent for computing STR using the standard method is about 1.5 times longer than using the new method. For three-atom molecules, the ratio almost doubles and for molecules of more than 5 atoms it is ∼3 fold. It is expected that the ratio would increase for higher number of atoms. The numerical computation time for the new method can be considerably reduced by using less number of points. In this example, to obtain a precise uSTR value with the new method, more than enough number of points were used (27). However, this number can be reduced at least by 2 or more to yield a faster calculation. The main goal of this section is to show that the applicability of equation (41) and it is not focused specifically on the timing performance.

Discussion and conclusions
The intention of this work is to present the developed relations from sections A.1 and A.2, their uses to obtain physically meaningful relations uncovered in sections 2.2, 2.3 and 2.4, and finally showing a simple application example with the unrestrained Stark Tuning Rate in section 3. According to the traditional way of computing the uSTR, the geometry of several structures at diverse electric field strengths ought to be optimised prior to the frequency analysis yielding the wavenumber. This procedure is computationally very expensive and is increasingly hindered by the number of atoms in the molecule. The relation obtained in equation (41) significantly reduces the computation time of the uSTR and may increase the precision in the calculation of its value. Additionally, the application of an electric field is not required to compute the unrestricted and restricted STRs. Incorporation of this method for computing the restrained and unrestrained first order and second order uSTR in traditional quantum computing packages is encouraged. With little extra computational cost, all these extra parameters can be obtained after the frequency analysis calculation. A more focused study to optimise the performance (speed and precision) of the implementation of the new STR method should be done. The study should determine and compare the speed and precision of this method with the traditional one with more diverse molecules. The novel mathematical relationships developed here can be used as shortcut equations to compute restrained/unrestrained first/second-order Stark Tuning Rate (STR) and Zeeman Tuning Rate (ZTR). Furthermore, an equation that connects the uSTR and rSTR has been found, in addition to other important electric or magnetic field-dependent properties of molecules relevant for vibrational spectroscopy. For electrodynamics, novel definitions for unrestrained dipolar moment and nuclear polarizability were introduced and an important relation between electronic, nuclear and total polarizabilities was developed. These relationships are obtained from mixed third order partial derivatives of the Legendre transformation of the energy which depend on the normal coordinates and the external uniform electric or magnetic field vector. A direct application of these relationships is the faster and efficient computation of the STR and ZTR. Additionally, this work can be a starting point for further future theoretical developments in this field. Whether this novel framework will prove to be useful for other practical applications is yet to be discovered. Q F Mixed third-order partial derivatives of E, R, T and Y Analogously, the equality of the mixed third-order partial derivatives [19] of E, R, T, and Y provides the following useful relations: 6

A.3. Specific demonstration 1 (normal force constant equality)
The normal force constant is defined as the second derivative of the energy with respect to a normal coordinate while the rest of the normal coordinates and all the elements of the electric field vector are kept constant.
Considering a system with n normal vibrational modes: The equality of both normal force constants is proven based on Jacobian properties [1] as shown in section A. 5: The following is the demonstration of equation (87). This demonstration can be extended for equations (88), (89) and (90).
The derivative of the left term of equation (87) can be rewritten using Jacobian's notation as shown in equation (91).
x Q The first variable to change, in this case Q 1 , should be in the second position. The first step is to multiply equation (91) by a unity Jacobian which is constructed by using the denominator of equation (91) with the first variable to change replaced: