Vibrational Study and Force Field of the Citric Acid Dimer Based on the SQMMethodology

1 Cátedra de Quı́mica General, Instituto de Quı́mica Inorgánica, Facultad de Bioquı́mica, Quı́mica y Farmacia, Universidad Nacional de Tucumán, Ayacucho 471, 4000 S. M. de Tucumán, Argentina 2 Cátedra de Fisicoquı́mica I, Instituto de Quı́mica Fı́sica, Facultad de Bioquı́mica, Quı́mica y Farmacia, Universidad Nacional de Tucumán, San Lorenzo 456, T4000CAN S, M. de Tucumán., Argentina 3 Centro de Quı́mica Inorgánica, CEQUINOR/CCT, Departamento de Quı́mica, Facultad de Ciencias Exactas, Universidad Nacional de La Plata, 47 y 115. CC. 962 (B1900AVV), 1900 La Plata, Argentina 4 Laboratorio de Fisicoquı́mica, Departamento de Fı́sica, Facultad de Ciencias Exactas y Tecnoloǵıa, Universidad Nacional de Tucumán, Avenida Independencia 1800, T4000CAN S. M. de Tucumán., Argentina


Introduction
The structural and vibrational studies of citric acid are of chemical, biochemical, and technological interest [1][2][3][4][5][6][7][8].The crystal and molecular structure of anhydrous citric acid was determined by Nordman et al. [6] and Glusker et al. [7] by means of X-ray crystal analyses.Both authors demonstrated that two molecules of the acid are linked through hydrogen bonds between their carboxylic groups originating the centre-symmetrical cyclic structure of the dimeric species.Many previous studies in relation to the chemical and spectroscopic properties of this acid as monomer have been reported [9][10][11][12][13][14] while an ab-initio/HF force field study for the vibrations of citric acid and citrate trianion was published by Tarakeshwar and Manogaran [15].
In the latter study, the calculated geometries and vibrational frequencies only for the monomer citric acid were used to assign the vibrational spectra.A complete vibrational characterization of this acid is very important to obtain information about chemical and physical components in food, biological materials, and pharmaceutical products at the molecular level [14,16].So far, there are no publications about experimental or high-level theoretical studies on the geometries and force field for the dimer citric acid.Hence, obtaining reliable parameters by DFT theoretical methods is an appealing alternative to have a precise knowledge of the normal vibration modes.The aim of this paper is to report a DFT theoretical and experimental study of dimer citric acid by infrared spectroscopy in order to effect a complete assignment of every observed bands in the vibrational Advances in Physical Chemistry spectra by taking into account the presence of the dimeric species in the solid phase.For this purpose, the optimized geometries and frequencies for the normal vibration modes, considering the dimer with two units of the free acid, were calculated.For a complete and reliable assignment of the compound, the DFT calculations were combined with the SQMFF methodology [17][18][19] in order to fit the theoretical frequencies values to the experimental ones.In addition, the force constants were calculated and compared with those reported for the monomer by using the HF/4-21G method [15].Here, we demonstrated that a DFT molecular force field for the dimer calculated by using the B3LYP/6-31G * combination is well represented.DFT normal mode assignments, in terms of the potential energy distribution, are in agreement with those obtained from the normal coordinate analysis.The total energy calculated with the 6-31G * basis set for citric acid dimer was corrected for basis set superposition error (BSSE) by the standard Boys-Bernardi counterpoise method [20].Besides, the dimer electronic properties were evaluated by NBO [21][22][23][24] and Atoms in Molecules (AIMs) [25,26] studies in order to analyze the nature and magnitude of the intermolecular interactions.

Experimental and Theoretical Calculations
A pure anhydrous Mallinckrodt commercial sample of citric acid (C 6 H 8 O 7 ) (CA) was used.The FTIR spectrum of the solid substance in KBr pellets was recorded with a Bruker IFS 66 FTIR-spectrophotometer from 4000 to 400 cm −1 .FT-Raman spectrum was also collected on a Bruker IFS 66 FT-IR spectrophotometer provided with a NIR Raman attachment equipped with an Nd:YAG laser at 1064 nm laser.
The starting point for the dimer geometry optimization was a structure with bond lengths taken from the Xray crystal structure [6,7].These initial geometries were optimized by using the hybrid B3LYP method [27,28] with the 6-31G * and 6-311++G * * basis sets, as implemented in the GAUSSIAN 03 program [29].The atom structure and labelling can be seen in Figure 1.The electronic charge density topological analysis was performed by using the Atoms in Molecules (AIMs) methodology [25], by the AIM2000 program package [26] while the NBO calculations were performed by using the NBO 3.1 program [24], as implemented in the GAUSSIAN 03 package [29].The harmonic wavenumbers and the valence force fields in Cartesian coordinates were calculated at the B3LYP/6-31G * and 6-311++G * * levels of approximation.The resulting force fields were transformed into "natural" internal coordinates by using the MOLVIB program [30,31].The natural internal coordinates are presented in Table S1 of the Supplementary Material at doi/0.1155/2011/347072 and have been defined as proposed by Fogarasi et al. [32].In this case, some coordinates are slightly different from those used in the previous study [15].The six inter-monomer coordinates for this dimer are similar to those defined by other authors for other dimeric acids [33][34][35][36], which can be seen in Supplementary Table S2.Following the SQMFF procedure [37,38] the harmonic force fields were evaluated at B3LYP/6-31G * level.In the scaling procedure the diagonal force constants associated with a particular symmetry coordinate are multiplied by a small number of scale factors ( f i , f j , etc.).As a result, each off-diagonal constant is automatically scaled by the geometrical mean of factors ( f i × f j ) 1/2 for the corresponding diagonal constants.Finally, the scale factors are adjusted by an iterative method until they reproduce the experimental wavenumbers as possible.In this case, we use the Rauhut and Pulay scale factors [17,18].The potential energy distribution (PEDs) components higher than or equal to 10% are subsequently calculated with the resulting SQM.The nature of every vibration modes was carried out by means of the GaussView program [39].S3 shows the total energies and the corresponding dipole moment values for the dimer structures with the B3LYP method by using 6-31G * and 6-311++G * * basis sets compared with the ones corresponding to two monomer units.The energy difference between the cyclic dimer and the two citric acid isolated monomers with the B3LYP/6-31G * method is 82.86 kJ mol −1 while by using the 6-311++G * * basis set the value at 68.64 kJmol −1 decreases.Reva and Stepanian [40] predict an energy difference of 86.00 kJ mol −1 between the cyclic dimer and two isolated monomers of benzoic acid (ab initio/4-31G) while the calculated energy difference between their monomers is 40.12 kJ mol −1 .In this case, the total energy for the citric acid dimer is smaller than the sum of the energies of each monomer (Supplementary Table S3) by using both basis sets.Such observation might mean that the dimer presence even in the gas phase is preferable to the isolated monomer, as observed in the 4-hydroxybenzoic acid dimer [36].This result is very important to carry out the assignment of the vibrational spectra of this acid.In general, our theoretical values are in agreement with the experimental values for the citric acid [6,7].The dimer stability was investigated by means of NBO [21][22][23] and AIM calculations [25].In the NBO analysis, the second-order perturbation energies E (2) (donor → acceptor) involving the most important delocalization for the dimer clearly reveal the energies and occupancies of the main contributions: LP(2)O19 → σ * O22-H21 (115.2 kJ mol −1 ), LP(2)O27 → σ * O19-H20 (38.4 kJ mol −1 ) and LP(2)O27 → σ * C3-H6 (2.4 kJ mol −1 ) to the delocalization energy due to the OH-O bonds between both structures.In Bader's topological analysis [25] the localization of the Bond Critical Point (BCP) and Ring Critical Point (RCP) in the ρ(r) and the Laplacian values ∇ 2 ρ(r) at these points is important to characterize molecular electronic structure.The critical point has the typical properties of the closed-shell interaction when the value of ρ(r) is relatively low, the relationship |λ 1 |/λ 3 is <1, and ∇ 2 ρ(r) is positive indicating that the interaction is dominated by the charge contraction away from the interatomic surface toward each nucleus.The topological properties for the citric acid dimer can be seen in Supplementary Table S4.The values obtained clearly reveal the formation of hydrogen bonds between two COOH groups of each monomer.In addition, the seven RCPs calculated by means AIM analysis for the dimer corroborate those contributions found in the NBO study.

Vibrational Analysis
The dimer's structure has C 1 symmetry and 120 vibrational normal modes where all vibrations are IR and Raman active.The recorded infrared spectra for the compound in the solid phase can be seen in Figure 2. The assignment of the experimental bands to the 120 expected normal vibration modes was made on the basis of the Potential Energy Distribution (PED) in terms of symmetry coordinates by using the B3LYP/6-31G * data and taking into account the previous assignment for the monomer [15] and for related molecules [33][34][35][36]40].Table 1 shows the observed frequencies and the assignment of the vibrational normal modes while the PED based on the 6-31G * basis set is shown in Supplementary Table S5.The theoretical calculations reproduce the normal wavenumbers by using the 6-31G * basis set, with an initial rmsd value of 11.5 cm −1 while when the SQMFF method is applied by using Pulay's scaling factors, the final value is 1.5 cm −1 .The best results are obtained with a B3LYP/6-31G * calculation because the used transferable scale factors are only defined for this method.The SQM force fields for the compound dimer can be obtained at request.The broad band observed between 3400 and 2200 cm −1 in the experimental IR spectra of the compound in the solid phase can be attributed to the O-H hydrogen bond formed by the spatial arrangement of molecules in the lattice crystal [41,42].Note that the theoretical infrared spectra for the monomer structure demonstrate a little agreement with the experimental spectrum (Figure 2), especially in the high wavenumbers cm −1 region.On the other hand, in the theoretical Raman spectrum shown in Figure 3, in relation to the experimental one, a higher intensity in the bands attributed to the OH stretching is observed.The existence of molecules in the solid phase as dimer justifies such a difference because this species restricts the motion of the COOH groups and is strongly affected by the double hydrogen bond, as was also observed in the 4-hydroxybenzoic acid in solid state [36].Thus, the wavenumbers in the dimer are shifted with respect to the monomer and the splitting of the antisymmetric and symmetric modes is particularly large.In this compound, both monomer and dimer structures are present in the solid phase and, for this reason, a comparison between the average calculated infrared spectra from B3LYP/6-31G * level for the monomer and dimer by using average wavenumbers and intensities with the corresponding experimental ones demonstrates a good correlation, especially in the high wavenumbers region.The discussion of the assignment of the most important groups for the compound (see Table 1 and Supplementary Table S5) is presented as follows.

Band Assignments
OH Modes.In the previous assignment for the monomer the theoretical HF/4-21G calculations [15] predict that the behaviour of the central COOH group is different from the terminal groups, as in our case by using the B3LYP/6-31G *  IR bands at 1340, 1308, and 1081 cm −1 , respectively, as seen in Table 1.In the region of lower wavenumbers, the OH torsion modes are expected and these modes are also affected by the presence of H bonding. Previously, those modes in the monomer were assigned, according to HF/4-21G calculations [15], to the bands at 645, 598, 556, and 383 cm −1 while by using the B3LYP/6-31G * level they are calculated between 655 and 387 cm −1 ; for this reason, the OH torsion modes are associated with the Raman bands at 642, 627, 492, and 415 cm −1 .
CH 2 Modes.These modes in the monomer were previously observed in the 3000-2900 cm −1 region [15].In this case, the out-of-phase and in-phase antisymmetric stretching modes of these groups can be assigned to the shoulder and IR band at 3035 and 2994 cm −1 , as indicated in Table 1, while the outof-phase and in-phase symmetric modes are assigned to the weak IR band at 2950 cm −1 and to the medium intensity Raman band at 2982 cm −1 .In the monomer, the scissoring modes were assigned at 1427 cm −1 while for the dimer those modes were also clearly calculated in this region; thus, they are assigned to the bands at 1430, 1424, 1420, and 1413 cm −1 .
In the above-mentioned study [15], the wagging modes were associated with the very weak IR band at 1325 cm −1 .Here, as predicted by calculation, the bands at 1389, 1365, 1308, 1292, and 1036 cm −1 are assigned to the wagging modes.Previously, in the monomer the twisting and rocking modes were assigned at 1350 and 966 cm −1 , respectively [15].In the dimer, the IR bands located at 1358, 1292, and 1214 cm −1 and the Raman band at 1285 cm −1 are associated with the rocking modes, while the expected twisting modes are associated with the IR bands at 904 and 881 cm −1 and with the shoulder in the Raman spectrum at 1036, 656 and 514 cm −1 .
COO Modes.In the previous studies for the monomer [15] the band at 1706 cm −1 and the shoulder at 1695 cm −1 were assigned to the C=O stretching modes.Note that the band at 1706 cm −1 (Table 1) was related to the strong Raman band at 1735 cm −1 and assigned to the same mode in spite of the difference between them (19 cm −1 ).In the dimer, the six expected C=O stretching modes are clearly predicted by SQM/B3LYP/6-31G * calculation at 1786, 1783, 1775, 1771, 1726, and 1677 cm −1 .On the other hand, the IR spectrum clearly shows the presence of more than one band in this region; thus, two very strong bands at 1756 and 1708 cm −1 can be seen, each band presenting a shoulder at 1745 and 1698 cm −1 , respectively.Therefore, the IR and Raman bands between 1756 and 1698 cm −1 are clearly assigned to those modes.In this way, those shoulders should be also associated with other C=O stretching modes, and so this fact would confirm the presence of the dimer in the solid state.In both HF/4-21G [15] and B3LYP/6-31G * calculations, the C-O stretching modes are coupled with the C-OH scissoring modes.In the monomer [15], those modes were assigned to the IR bands at 1195, 1141, and 1108 cm −1 , while in the dimer, the IR bands at 1242, 1174, and 1140 cm −1 and the medium intensity Raman band at 904 cm −1 are assigned to these stretching modes.Previously, the medium intensity band in the IR spectrum at 1242 cm −1 was not observed and therefore not assigned.Taking into consideration their relative position, intensities, and the splitting predicted by calculation, the expected out-of-plane deformation modes (γCOO), in-phase and out-of-phase modes, were predicted as pure modes between 792 and 485 cm −1 and coupled with other modes at 655 and 499 cm −1 ; for this reason, they were assigned to the bands at 807, 795, 686, 574, and 492 cm −1 , as observed in Table 1.In the previous study [15], only one deformation mode (δCOO) was assigned to the band at 420 cm −1 .Here, as predicted by calculation, those expected deformation modes for the dimer are assigned to the IR band and shoulder at 781 and 789 cm −1 , respectively, and to the Raman bands at 700, 627, 574, and 541 cm −1 .
The rocking modes (ρCOO) for the monomer were not previously assigned and, in this case, the ρCOO modes were easily assigned to the IR bands at 599, 550, and 438 cm −1 and to the Raman bands at 367 and 348 cm −1 .Finally, the twisting (τwCOO) modes for the monomer were only predicted by the HF/4-21G level at 88, 68, 44, and 39 cm −1 and, in our case, for the dimer they are predicted strongly coupled between them and in the low wavenumbers region; this way, they could not be assigned because there are not observed bands in the vibrational spectra in those regions.
Skeletal Modes.In the previous assignment [15], the bands in the infrared spectrum at 1090, 1053, 900, 850, and 786 cm −1 were attributed to the C-C stretching modes, while the CCC deformation modes were not assigned because they appear at very low frequencies.Here, the bands in the IR spectrum at 1430, 1051 and 945 cm −1 and in the Raman spectrum at 1469, 1053, 966, 850, and 842 cm −1 are assigned to the C-C stretching modes.The calculations predicted the CCC deformation modes in the lower wavenumbers region; for this reason, they are assigned to the Raman bands between 323 and 212 cm −1 while the two modes predicted at 108 and 92 cm −1 were not assigned.
Intermonomer Modes.The intermonomer vibrational modes of the dimer usually appear in the 270-10 cm −1 region [33][34][35][36]40] and they are related to restricted translations or rotations of one molecule against the other (see Supplementary Table S2).Here, the six intermonomer modes appear strongly coupled with other modes and they are observed in the Raman spectrum with the exception of the γ a (OH-O) mode that was predicted by calculations with a higher contribution PED (42%) at 827 cm −1 ; for this reason, it is assigned to the very weak Raman band at 842 cm −1 .Table 1 shows the theoretical assignments for the dimer intermonomer vibrational modes.

Force Field
The force constants for the dimer were estimated by using Pulay et al. [37,38] scaling procedure as mentioned above.
The force constants expressed in terms of simple valence internal coordinates were calculated from the corresponding scaled force fields by using the MOLVIB program [30,31].It is interesting to compare the principal force constants, which were collected in Table 2 and calculated at the B3LYP/6-31G * level for some vibrations, with those reported for the monomer of this acid by using the HF/4-21G method [15].The calculated force constants values for the monomer by using the ab-initio method are slightly different from the values for the dimer by using the B3LYP/6-31G * level, as expected.The HF method overestimates the calculated force constants values due to the fact that it neglects the electron correlation and unharmony.On the other hand, the little variation in the f(C=O) force constant value in the monomer with reference to the dimer is justified because in the latter species two molecules of the acid are linked by hydrogen bonds between their carboxylic groups; therefore, the involved C=O bonds decrease the double bonds character.For this reason, a lower value in the f(νO-H) force constant can be also seen in the dimer in relation to the monomer.Obviously, the force constants related to the CH 2 groups in the monomer have lower values than the dimer because in the latter species the number of those groups is higher.

Conclusions
We conclude the following.
(i) The assignment previously made for the monomer [15] was corrected and completed in accordance with the present theoretical B3LYP/6-31G * results.Also, the assignment of the intermonomer vibrational modes for the dimer is carried out and the assignments of the 120 normal vibration modes corresponding to citric acid dimer are reported.
(iii) Numerous bands of different intensities observed in the vibracional spectra not previously assigned, in this work were assigned to the citric acid dimer.
(iv) The NBO analysis reveals that in the citric acid dimer there are strong intramolecular interactions of charge transfer from O lone pairs to the σ * (O-H) antibonds.
(v) The AIM studies of the charge density for the dimer satisfy the criteria of hydrogen bond interactions.
(vi) The complete B3LYP/6-31G * force field for the citric acid dimer was determined as well as the principal force constants for stretching and deformation modes.

Figure 1 :
Figure 1: Theoretical structure of the citric acid dimer.The H-bonding is indicated by dashed lines.

Figure 3 :
Figure 3: Comparison between (a) the calculated Raman spectrum from B3LYP/6-31G * level for the monomer, (b) the calculated Raman spectrum obtained from B3LYP/6-31G * level using average frequencies and intensities for both monomer and dimer species and (c) the experimental Raman spectrum in solid state.

Table 1 :
Observed and calculated wavenumbers (cm −1 ) and assignments for citric acid dimmer.
[15]tretching; δ: scissoring; wag: wagging; γ: out-of plane deformation; ρ: rocking; τ: torsion, τw: twisting; a: antisymmetric; s: symmetric; s: strong; m: medium; w: weak; v: very; sh: shoulder; br: broad.aThiswork.bCalculatedbyHF/4-21G method from[15].cTheoretical values from SQM/B3LYP/6-31G * calculations for citric dimer acid.#Intermonomer coordinates.leveland, as a consequence, the central OH stretching mode belonging to that central group is calculated at higher wavenumbers than the other ones.In the dimer, those modes are assigned according to the SQM calculations; for this reason, the OH stretching modes are associated with the Raman band at 3535 cm −1 , the IR bands at 3498, 3446, 3291, and 2933 cm −1 , and the shoulder located at 3350 cm −1 .The C-OH scissoring mode belonging to the COOH groups involved in the H-bonding formation is calculated at higher wavenumbers than the other ones (1483 and 1402 cm −1 ) probably because they are coupled with the corresponding C-O stretching modes; for this reason, both modes are assigned to the bands at 1493 and 1389 cm −1 .The remaining C-OH scissoring modes are assigned to the shoulder and * level using average frequencies and intensities for both monomer and dimer species.

Table 2 :
Scaled force constants for the citric acid dimer.