Xun Zhang

Although terahertz waves are susceptible to water of crystallization (bound water), we found that the terahertz absorption spectra of cysteine hydrochloride (LCH) and its monohydrate (LCHM) are highly similar. To explain this particular phenomenon, density functional theory (DFT) and the independent gradient model (IGM) were used to obtain the vibration mode and intermolecular interaction of LCH and LCHM. The molecular polarities of LCH and LCHM were then obtained by calculating their molecular polarity index (MPI). The characteristic peak positions in the terahertz spectra of LCH and LCHM basically corresponded, with the superimposed interference of vibration modes and the van der Waals interaction between molecules concealing the expression of hydrogen bonds produced by bound water in the LCHM terahertz spectrum. In addition, the intensity of the characteristic peaks in the LCHM terahertz spectrum was higher because of its higher molecular polarity. In general, the analysis method combining THz-TDS and MPI provides a new theoretical reference for studying the relationship between biomolecules and water.

Terahertz waves are electromagnetic waves located between the microwave and infrared spectral regions, with a wavelength range of 0.03–3 mm (1). The timescales for the dynamics of water molecules affected by the presence of a molecule fall in the terahertz frequencies, specifically the sub-picosecond region (2). Therefore, in general, the changes of water of crystallization in biomolecules can be observed by the corresponding terahertz spectra (3–5). Pan and associates (6) found that the terahertz spectral characteristic peak position deviation of L-phenylalanine and its monohydrate originated from the hydrogen bonding produced by the water of crystallization in the monohydrate. Similarly, Yan and associates (7) studied the dehydration kinetics of methylene blue (MB) at different temperatures through the differences in the position of the characteristic absorption peaks of MB and its two hydrates in the terahertz spectra. These studies showed that most of the terahertz spectral differences between anhydrate and hydrates were mainly reflected in the position of the characteristic peaks. However, we found that the terahertz spectra characteristic peak positions of cysteine hydrochloride (LCH) and its monohydrate (LCHM) corresponded, and the latter’s intensity was higher, which is inconsistent with existing research.

Hydrochloric acid derivatives of cysteine generally include anhydrous form (LCH) and monohydrate form (LCHM). LCH is used as an antioxidant (8) and a scavenger of active oxygen in food chemistry (9). LCHM is an effective radiation protection drug (10) which can be used to prepare nonlinear optical crystals with excellent properties (11). Considering the importances of LCH and LCHM in medicine and optics, we have conducted a deeper study of their optical characteristics and other essential characteristics in the terahertz band. LCH and LCHM were selected as the experimental subjects, and their terahertz spectra were measured in the frequency range of 0.5 to 2.5 THz at 298 K. DFT and IGM were used to analyze the molecular interactions and the low-frequency vibration sources, respectively. MPI was calculated to study the influence of polar molecules on the terahertz spectra.

L-cysteine hydrochloride (LCH, assay ≥ 98%, CAS: 52-89-1) and L-cysteine hydrochloride monohydrate (LCHM, CAS: 7048-04-6) were the derived comparison sample; the high-density polyethylene (particles size 40–48 μm) powder was the substrate. All reagents were obtained from Sigma-Aldrich and without additional purification. The molecular structures of all reagents are shown in Figure 1 (12). All tablet densities must be uniform to obtain stable spectral data, so all reagents were passed through a 200-mesh sieve during preparation.

FIGURE 1: Molecular straw structures of LCH and LCHM. The atoms are designated as H (white), C (cyan cylinders), O (red), N (blue), S (yellow), and Cl (cyan spheres). Note the independent cyan spheres illustrated in LCH and LCHM are chlorine atoms.

Each reagent was mixed with polyethylene powder at a ratio of 1:5 (reagent: 40 mg; polyethylene: 200 mg). After pressing the mixtures at 30 MPa for 3 min, the tablets with 13 mm diameter and 1.9 ± 0.05 mm thickness were obtained. THz-TDS were performed by the CIP–TDS (Da Heng New Epoch Technology Inc.). The dynamic range for the whole setup was better than 70 dB, and the frequency resolution was less than 5 GHz. During the measurement, the THz path was fully sealed, and high purity nitrogen was continuously purged to keep the humidity ≤ 3.7% to avoid water absorption. To get a good signal-to-noise (S/N) ratio, 300 scans were recorded and averaged. The spectra collection was carried out in a constant temperature and humidity laboratory with a temperature of 25 °C and a humidity of 30%. An HDPE tablet (200 mg) was also measured, and the obtained spectra were used as the reference in calculating the sample absorbance spectra. The baseline spectra were stable since the particle size and tablet thickness were adjusted, so other spectral processes, such as normalization, were not required.

We needed to further process the obtained THz time-domain spectra to obtain optical parameters of the sample, such as absorbance, absorption coefficient, and refractive index. By applying the fast Fourier transform (FFT), we would obtain the corresponding frequency-domain signals, which contain the THz amplitude of sample (Esam [ω]) and reference (Eref [ω]). Then, according to the Fresnel formula (13), the THz absorbance A(ω) can be obtained as follows:

TGA was performed using PYRIS 1 (Perkin Elmer). The samples were put in a corundum crucible and heated from 25 to 400 °C with a heating rate of 10 °C/min under a nitrogen atmosphere. The thermogravimetric (TG) curves of each reagent are shown in Supplementary Materials Figure S1 (at the bottom of the article).

XRD patterns were performed using D8 Discover (Cu source, 50 kV voltage, 40 mA filament emission, Bruker Corporation). The samples were laid on the holder’s surface of the diffractometer and scanned through a 2θ range of 5–50°, with an angle of 0.2°.

Geometry Optimization and Vibration Analysis

The quantum chemistry calculation was based on solid-state density functional theory (DFT). All of the calculations were performed on the crystalline state within the generalized gradient approximation (GGA) at Perdew-Burke-Ernzerhof correlation (PBE) functional (14). Grimme’s dispersion-corrected method was used for DFT-D correction and norm-conserving pseudopotential, as implemented in CASTEP-20.11 (15) (Cambridge Sequential Total Energy Package). The crystalline structures of L-cysteine hydrochloride and its monohydrate were downloaded from Cambridge Crystallographic Data Centre (CCDC). The measured PXRD spectra of L-cysteine hydrochloride and its monohydrate were confirmed to be consistent with previously published unit cell parameters (16) (see Supplementary Materials Figure S2). There was no imaginary frequency in the calculated results for both LCH and LCHM, and no scaling factors were applied to the calculated frequency values. The quality of energy calculation was set to ultrafine. The plane-wave cut off energy was 830 eV.

Interaction Analysis and Polarity Analysis

To further explain the spectral differences, we combined Multiwfn (17) and VMD for independent gradient model (IGM) analysis, as proposed by Corentin and associates (18), and molecular polarity index (MPI) analysis (19). Unlike the traditional calculation methods, the interference between the atomic density gradients did not affect the phase when using IGM to calculate the initial molecular electron density gradient. The interaction area can be judged clearly through the graphical display of the three-dimensional real space function δg.

) represents the absolute value of each component of the 

 vector (still maintaining the vector form), and | | represents the modulus of the vector. 

 describes the interactions between all atoms in the system, and can be divided into δg

 (showing intrafragment interactions) and 

 (showing interactions between fragments). The 

 formula required in the research is as follows:

 is the fragment number of the molecule in the above formula. To visualize the calculated results when representing the interaction between molecules, the sign (

) rho function (20) was introduced, which is the symbol of the eigenvalue 

 of the electron density Hessian matrix, which can be abbreviated as the sign (

). The noncovalent interactions between molecules can be judged using the sign (

) rho function to project the isosurface of 

 in different colors. Among them, electrostatic attraction (hydrogen bond), van der Waals interaction, and electrostatic repulsion (nonbonded overlap) are expressed as blue, green, and red isosurfaces, respectively. It is worth noting that the IGM calculation must be based on the optimized molecular structure.

When discussing the relationship between samples’ polarity and their THz spectra, MPI was used (19) to quantify the understanding of polarity. This parameter can measure the system’s polarity based on the characteristics of the electrostatic potential distribution on the molecular surface. The integral form of the MPI formula is as follows:

In the above formula, the molecular surface 

 is the molecular electrostatic potential; and 

 is the molecular surface area. The uniformity of the charge distribution in the system reflects the polarity of the molecule. Uneven charge distribution will lead to extreme areas of the electrostatic potential on the molecular surface, which would increase the MPI. We calculated the MPI value of different substances based on Multiwfn and established a connection with the dielectric constant obtained from the spectrum.

Observation and Analysis of L-cysteine Hydrochloride and its Monohydrate The measured THz time-domain signals were converted into absorption spectra by FFT, and the absorption spectra of cysteine and its derivatives are plotted in Figure 2. Under the test conditions, there is a significant noise increase in the high-frequency part of the spectrum of these substances, which increased the difficulty of locating characteristic peaks. The THz absorption spectra of LCH and LCHM are characterized by absorption peaks at 1.15, 1.76, and 2.40 THz. The characteristic peak positions of the LCH and LCHM are almost the same, but the characteristic peak widths are slightly different at 1.76 THz. We found that the experimental spectral results of LCH and LCHM only have slight differences in the intensity of the characteristic peaks, which is quite different from previous research. To explain these phenomena, we conducted the following series of studies based on quantum chemical calculations.

FIGURE 2: The experimental THz absorption spectra of LCH and LCHM.

At first, geometric optimization and vibration analysis were performed on the unit cell so as to analyze the source of terahertz characteristic absorption through DFT. Before obtaining the vibrational properties of these samples, Grimme’s dispersion-corrected method was combined to optimize their structures. The THz experimental absorption spectra and calculated spectra of the samples are summarized in Figures 3a and 3b. It can be seen that the calculated spectrum is relatively consistent with the experimental spectrum, with the subtle differences mainly due to the differences in temperature and crystal structure between actual testing and calculation simulations. This phenomenon shows that the chosen basis set and method are suitable for the current experiment, and can be used for further analysis (6,21). Additionally, we combined the VMD to perform vibration analysis on the calculation results to distinguish different vibration modes.

FIGURE 3: The experimental and calculated THz absorption spectra of (a) LCH and (b) LCHM. Also showing terahertz spectroscopy-active vibration modes of (c) LCH and (d) LCHM, and relevant simulated peak positions.

The detailed vibration mode analysis results of the two substances are plotted in Table I and Figures 3c and 3d. For LCH, the experimental peaks at 1.15 THz and 1.76 THz correspond to the calculated peaks at 0.97 THz and 1.65 THz, respectively. The former was led by the rocking vibration of the molecule, and the latter was led by the translational vibration of the molecule along the 

-axis. For LCHM, the experimental peak at 1.15 THz corresponds to the calculated peak at 1.00 THz, which was led by the rocking vibration of the LCHM molecules, and the experimental peak at 1.76 THz corresponds to the calculated peaks at 1.51 THz and 1.69 THz, which were led by the translational vibration of LCHM molecules along the 

-axis and the torsional vibration along the 

-axis, respectively. It is worth mentioning that the vibration mode of the characteristic peak at 2.41 THz in the LCH and LCHM experimental spectra is uncertain, because the corresponding calculated peak may exceed the test range. Additionally, in the LCHM experimental spectrum, a broad experimental characteristic peak was formed by the superposition of two vibration modes, which may be led by various reasons. The temperature difference between the experiments and the calculations may be an influencing factor (22), and other factors still need further verification.

IGM Analysis of L-cysteine Hydrochloride and its Monohydrate 

In addition to sensing the low-frequency vibration of lattice and molecular skeleton, THz-TDS is also very sensitive to the change of weak interaction between molecules (23,24). To obtain the corresponding information on the interactions within and between the fragments of the substance molecules, the IGM was selected in this study because the interaction between molecules in the selected area can be plotted with the information provided by the IGM. Therefore, the type and strength of the interactions in the system are visualized. The unit cell structures of LCH and LCHM were optimized using the method in 3.1, and the optimized unit cell structures were analyzed by IGM. The IGM isosurface of the two substances is shown in Figure 4a and 4b. The electrostatic attraction (such as hydrogen bonds), van der Waals effect, and electrostatic repulsion (such as non-bonded overlap) were marked with blue, green, and red, respectively.

FIGURE 4: IGM isosurface maps for (a) LCH and (b) LCHM models. Intermolecular interaction type derived from the value of sign( 2) (blue: sign(2) < −0.005 au, sign(2) (green: −0.005 au < sign(2) < 0.005 au, red: sign(2) > 0.005 au). where represents the quantitative parameter of the critical point of weak interaction in the IGM theory. (c) Dielectric constant spectra (real part of the dielectric spectra) of LCH and LCHM. (d) The relative permittivity and the normalized MPI results of the LCH and LCHM.

There is an apparent green isosurface between the chlorine atom and the cysteine molecule in LCH, indicating that the van der Waals interactions between LCH molecules are dominant. The inside of LCHM is mixed with the combined effect of hydrogen bond and van der Waals, and the intermolecular interaction is still dominated by van der Waals. It should be noted that the van der Waals force between LCHM molecules is stronger than LCH. Although the unit cell conformation of the LCH and LCHM varies greatly, the type of intermolecular interaction is still dominated by van der Waals interaction.

Based on the IGM analysis, we found the reason for the high similarity between the LCH and LCHM terahertz experimental spectra. After the water molecule and cysteine hydrochloride are combined, the hydrogen bond does not constitute the dominant interaction between the molecules—on the contrary, hydrogen bonds formed a stronger mixed interaction dominated by van der Waals interaction. Therefore, it is difficult to find the difference in the position of the spectral characteristic peaks on the terahertz spectra of LCH and LCHM. It is worth noting that the characteristic peak intensities are different in the terahertz spectra of LCH and LCHM. Considering that water molecules are polar, the addition of water molecules may increase the polarity of the LCHM, thereby affecting the intensity of the corresponding terahertz spectrum. To verify this conjecture, firstly, we calculated the MPI of LCH and LCHM, and then obtained the relative permittivity of the two through terahertz spectra to obtain the polarities of the two and analyzed the influence of the polarity on their corresponding terahertz spectra.

Polarity Analysis of L-cysteine Hydrochloride and its Monohydrate

Because polar molecules like water molecules affect the absorption of THz (25), we further explored the influence of material polarity on the intensity of THz spectra. First, the dielectric constant spectra of the two substances are obtained through THz dielectric spectra. Previous studies have shown that the relative permittivity of a molecule is closely related to its polarity (19). Generally, the higher the polarity, the higher the relative permittivity of the molecule. We would take the average of the test data of three sets of parallel samples with uniform thickness to reduce the error caused by the tablet density or porosity. The experimental results are plotted in Figure 4c. It can be seen that the dielectric constants of the samples show that LCH is less than LCHM. To prevent the influence of high-frequency scattering effects, the relative permittivity data at 0.5 THz, 1.0 THz, and 1.5 THz are selected to measure the two substances’ polarity further; these values are summarized in Table II.

However, the relative permittivity is relatively unstable because of external influences. To verify the relative permittivity of the substances, we introduced MPI based on the calculation of the molecular electrostatic potential distribution (26). This parameter judges the polarity by using the characteristics of the electrostatic potential distribution on the molecular surface. We have calculated the MPI of the single molecule and the optimized unit cell, respectively. The MPI of the unit cell is the average value of the MPI calculated based on the unit cell in different vibration modes. The calculated results are also summarized in Table II. The order of polarity based on MPI is consistent with the result of relative permittivity.

To more intuitively observe the relationship between the polarity of the substance and the terahertz spectrum, the results obtained by the two methods are normalized and plotted in Figure 4d. Combining Table II and Figure 4d, we can see that the relative dielectric constants of the four materials at 1.0 THz and 1.5 THz are highly similar to MPI. However, the changing trend of the relative permittivity at 0.50 THz is quite different from the other two places because of the error caused by the background noise in the low-frequency spectrum (27). This result verified the relative consistency between the relative permittivity obtained by terahertz dielectric spectroscopy and the MPI results. The relative dielectric constants establish a connection between MPI and terahertz absorption spectra. Integrating the MPI and terahertz absorption spectra results, we can infer that the stronger the material’s polarity, the higher the intensity of the terahertz absorption characteristic peak, and the less its THz spectrum is affected by background noise and other factors. In addition, based on the experimental results of LCH and LCHM, this method shows high sensitivity.

The MPI results of single molecule and unit cell present consistent after normalization, which shows that MPI is hardly affected by molecular conformation and local polarity. Therefore, MPI can be used in terahertz spectroscopy because of its better anti-interference ability and stronger stability. Therefore, the polarity of the substance can be judged based on MPI before the spectrum test. The absorption intensity of the terahertz spectrum of the substance can then be predicted. Similarly, when analyzing terahertz spectra with similar structures, the difference in terahertz spectra can be predicted efficiently based on the MPI differences of substances.

In this study, we combined a series of quantum chemical calculations to conduct a systematic theoretical analysis of the difference between the terahertz spectra of LCH and LCHM. First, the characteristic vibration modes of the LCH and LCHM terahertz spectra were analyzed by DFT with dispersion correction. It was found that the characteristic peak at 1.76 THz in the experimental spectrum of LCHM is produced by the superposition of two vibration modes.

Second, the intermolecular interactions of LCH and LCHM were compared through IGM analysis, and it was found that both were dominated by the van der Waals effect. This finding indicates that the van der Waals interaction may conceal the expression of the hydrogen bonds generated by the water of crystallization in the LCHM terahertz spectrum.

Finally, we found that the more polar LCHM showed higher characteristic peak intensity in the terahertz absorption spectrum in further MPI analysis. Therefore, intermolecular interaction analysis and polarity analysis are also necessary for the study of terahertz spectroscopy. In short, this terahertz method provides a new theoretical reference for studying the relationship between biomolecules and water, which can have value in the fields of biology and medicine.

Figure S1 shows the PXRD information of L-cysteine hydrochloride and its monohydrate. Compared with the previously published unit cell parameters, we found that it is in good agreement with the experimental data. This shows that our experimental materials are reliable, and that the experimental materials’ lattice information is the same as the database.

FIGURE S1: PXRD data of (a) LCH and (b) LCHM.

Figure S2 shows the TGA information of L-cysteine hydrochloride and its monohydrate. It can be seen that the water molecules escaped when the LCHM was at 143 

C. This verifies the accuracy of the experimental materials.

FIGURE S2: TGA data of LCH (green) and LCHM (red).

(1) Liang, J.; Zhang, X.; Wang, N.; Chang, T.; Cui, H. L. Vibrational Spectra of Pyrazinamide and Isoniazid Studied by Terahertz Spectroscopy and Density Functional Theory. 

(2) Morales-Hernandez, J. A.; Singh, A. K.; Villanueva-Rodriguez, S. J.; Castro-Camus, E. Hydration Shells of Carbohydrate Polymers Studied by Calorimetry and Terahertz Spectroscopy. 

(3) Takahashi, M.; Ishikawa, Y. Terahertz Vibrations of Crystalline α-D-glucose and the Spectral Change in Mutual Transitions Between the Anhydride and Monohydrate. 

(4) Zhang, B.; Li, S.; Wang, C.; Zou, T.; Pan, T.; Zhang, J.; Xu, Z.; Ren, G.; Zhao, H. Terahertz Spectroscopic Investigation of Gallic Acid and its Monohydrate. 

(5) Zong, S. Q.; Ren, G. H.; Li, S. P.; Zhang, B.; Zhang, J. B.; Qi, W. P.; Han, J. G.; Zhao, H. W. Terahertz Time-Domain Spectroscopy of L-histidine Hydrochloride Monohydrate. 

(6) Pan, T.; Li, S.; Zou, T.; Yu, Z.; Zhang, B.; Wang, C.; Zhang, J.; He, M.; Zhao, H. Terahertz Spectra of l-phenylalanine and its Monohydrate. 

(7) Yan, S. H.; Zhang, H.; Yang, Z. B.; Tang, M. J.; Zhang, M. K.; Du, C. L.; Cui, H. L.; Wei, D. S. Transformation and Dehydration Kinetics of Methylene Blue Hydrates Detected by Terahertz Time-Domain Spectroscopy. 

(8) Li, T. T.; Wu, Q. X.; Zhou, Y. J.; Yun, Z.; Duan, X. W.; Jiang, Y. M. L-Cysteine Hydrochloride Delays Senescence of Harvested Longan Fruit in Relation to Modification of Redox Status. 

(9) Duan, X.; Zhang, H.; Zhang, D.; Sheng, J.; Lin, H.; Jiang, Y. Role of Hydroxyl Radical in Modification of Cell Wall Polysaccharides and Aril Breakdown During Senescence of Harvested Longan Fruit. 

(10) Ren, G.; Zong, S.; Zhu, Z.; Cheng, C.; Chen, L.; Zhou, L.; Zhang, J.; Liu, L.; Han, J.; Zhao, H. Far-Infrared Terahertz Properties of L-cysteine and its Hydrochloride Monohydrate. 

(11) Prasad, P. V.; Visweswara Rao, T. K.; Ramachandra Rao, K.; Satya Kamal, C.; Samuel, T. Studies on Influence of Cd(

) Ions in Unidirectional Growth and Characterization of l-Cysteine Hdrochloride Monohydrate Single Crystals. 

(12) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. 

(13) Dorney, T. D.; Baraniuk, R. G.; Mittleman, D. M. Material Parameter Estimation with Terahertz Time-Domain Spectroscopy. 

(14) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. 

(15) Clark, S. J.; Segall, M. D.; Pickard, C. J.; Hasnip, P. J.; Probert, M. I. J.; Refson, K.; Payne, M. C. First Principles Methods using CASTEP. 

Zeitschrift für Kristallographie - Crystalline Materials

(16) Chapman, R. P.; Bryce, D. L. A High-Field Solid-State 

Cl NMR and Quantum Chemical Investigation of the Chlorine Quadrupolar and Chemical Shift Tensors in Amino Acid Hydrochlorides. 

(17) Lu, T.; Chen, F. Multiwfn: A Multifunctional Wavefunction Analyzer. 

(18) Lefebvre, C.; Rubez, G.; Khartabil, H.; Boisson, J. C.; Contreras-Garcia, J.; Henon, E. Accurately Extracting the Signature of Intermolecular Interactions Present in the NCI Plot of the Reduced Density Gradient Versus Electron Density. 

(19) Liu, Z. Y.; Lu, T.; Chen, Q. X. Intermolecular Interaction Characteristics of the All-Carboatomic Rng, Cyclo[18]carbon: Focusing on Molecular Adsorption and Stacking. 

(20) Johnson, E. R.; Keinan, S.; Mori-Sanchez, P.; Contreras-Garcia, J.; Cohen, A. J.; Yang, W. Revealing Noncovalent Interactions. 

(21) Allen, J. L.; Sanders, T. J.; Horvat, J.; Lewis, R. A. Anharmonicity-Driven Redshift and Broadening of Sharp Terahertz Features of Alpha-Glycine Single Crystal from 20 K to 300 K: Theory and Experiment. 

(22) Chen, L.; Ren, G.; Liu, L.; Guo, P.; Wang, E.; Zhu, Z.; Yang, J.; Shen, J.; Zhang, Z.; Zhou, L.; et al. Probing NaCl Hydrate Formation from Aqueous Solutions by Terahertz Time-Domain Spectroscopy. 

(23) Ding, L.; Fan, W. H.; Song, C.; Chen, X.; Chen, Z. Y. Terahertz Spectroscopic Investigation of Salicylic Acid and Sodium Salicylate. 

(24) Wang, F.; Zhao, D.; Dong, H.; Jiang, L.; Liu, Y.; Li, S. Terahertz Spectra of DNA Nucleobase Crystals: A Joint Experimental and Computational Study. 

(25) Li, J.; Yao, Y. X.; Jiang, L. W.; Li, S.; Yi, Z. H.; Chen, X. Y.; Tian, Z.; Zhang, W. L. Time-Domain Terahertz Optoacoustics: Manipulable Water Sensing and Dampening. 

(26) Lu, T. Multiwfn Manual, version 3.7, Section 3.15.1, available at 

(27) Nickel, D. V.; Garza, A. J.; Scuseria, G. E.; Mittleman, D. M. The Isotropic Molecular Polarizabilities of Single Methyl-Branched Alkanes in the Terahertz Range. 

 are with the College of Textile Science and Engineering (International College of Silk) at Zhejiang Sci-Tech University, in Hangzhou, Zhejiang, PR China. Bin Yang is also with the Key Laboratory of Advanced Textile Materials and Manufacturing Technology of the Ministry of Education at Zhejiang Sci-Tech University. Direct correspondence to: 

Articles

This new terahertz method provides a theoretical reference for studying the relationship between biomolecules and water.

Terahertz Spectral Investigation of L-Cysteine Hydrochloride and its Monohydrate

Terahertz (THz) resonance absorption originates from intermolecular interactions, which are suitable for identifying amino acids with multiple isomers. L-threonine and L-allo-threonine are diastereomers with two characteristic peaks in the effective spectrum range of 1.0–2.3 THz, which are located at 1.42 and 2.14 THz for L-threonine (L-thr) and 1.63 and 2.16 THz for L-allo-threonine (L-allo-thr). Based on the density functional theory (DFT) of the crystal structures of L-thr and L-allo-thr, the vibration frequencies of 1.56, 1.87, 2.16 THz, and 2.22 THz were obtained, corresponding to the THz characteristic peaks. Through vibration model analysis, it was found that lattice and skeleton vibrations mediated by intermolecular hydrogen bonds play a crucial role in the THz response. Studying the experimental absorption spectra of different proportions co-crystallized mixtures and 1:1 physical mixture of L-thr and L-allo-thr, it was found that the characteristic peaks of the physical mixture include the characteristic peaks of the two diastereomers in the THz band, while amino-acid co-crystallized mixtures formed their own characteristic peaks depending on the proportion. The results show that the co-crystallized mixture composition of diastereomers can be quantitatively analyzed by THz time-domain spectroscopy (THz-TDS).

As one of the eight essential proteinogenic amino acids, L-threonine (L-thr) is required for body growth, uric acid metabolism, and the immune system. L-allo-threonine (L-allo-thr) is used as a starting material in the chemical synthesis of lysobactin, which displays strong antibiotic activity (1). L-thr and L-allo-thr are diastereomers and the molecular formula is C

, but the molecule configuration is different. The main difference between the two isomers is that the chirality of the center carbon atom C3 can be interchanged through the chiral conversion of the C

 atom (2,3). Another characteristic of the two isomers is that they can also be completely miscible to form a continuous set of co-crystallized mixtures (4,5). L-thr and L-allo-thr have the same functional groups (amino, carboxyl, and side-chain hydroxyl), and the same type of hydrogen bonds. Among them, proton donors, such as -OH and -NH

, and proton acceptors, such as oxygen and nitrogen atoms, can form not only intramolecular hydrogen bonds but also intermolecular hydrogen bonds. Therefore, the three structures of L-thr, L-allo-thr, and L-thr-L-allo-thr all have similar “framework” characteristics comprised of hydrogen bonds.

A terahertz (THz) wave refers to electromagnetic wave with frequencies in the 0.1–10 THz frequency band. The intermolecular and intramolecular interaction energies of many polar macromolecules are in tune with photons of frequencies in the THz band. For compounds with similar molecular structures and the same functional groups, their infrared (IR) spectroscopies are similar or the same, whereas THz can distinguish intermolecular hydrogen bonds, which is mainly because the IR resonance absorption originates from intramolecular vibrations. The molecular arrangement and intermolecular interactions of these compounds are different where most of the THz resonance absorption originates. Therefore, THz spectroscopy is more suitable for identifying such substances than IR spectroscopy (6,7). The 20 basic amino acids that comprise various proteins in organisms, except for glycine (the α-carbon atoms of other amino acids), are all chiral carbon atoms. Therefore, amino acids are isomers of various configurations. For various amino acids, the vibrational and torsional motions of groups of atoms, hydrogen bondings, van der Waals forces, and lattice vibrations have strong resonance absorption in the THz band, and they have their own distinct characteristic peaks, which can form a unique THz “fingerprint” spectra (8–10). For biomolecules with a large number of weak molecular interactions (such as hydrogen bonds), low-frequency vibration strongly affects their structure, conformation, and function because of the complexity of this type of vibration; the analysis of the THz spectrum lacks the corresponding molecular vibration model. Currently, solid density functional theory (DFT) simulation is considered to be an effective way to connect the predicted vibration mode and molecular structure. Simulated spectrum calculation refers to modeling with the aid of quantum mechanics calculation software to obtain a theoretical spectrum. Combined with experimental data and related visualization software, the vibration mode of absorption peaks can be identified, and the source of these absorption peaks can be explained.

In recent years, DFT has been well ap- plied to simulate the THZ spectrum of threonine molecules. In 2009, Wang and others used THz time-domain spectroscopy (THz-TDS) and Raman spectroscopy to test L-thr (11), compare THz absorption and Raman spectroscopy, and then combine them with quantum chemistry calculations to analyze threonine. In 2019, Li Wei and others studied the THz absorption spectra of three configurations of single molecule, dimer, and unit cell of L-thr samples (12) before comparing them with the calculated THz spectra of single molecule and dimer configurations. Their study determined that the calculated THz spectra of crystal cell configurations are in good agreement with the experimental spectra, indicating that the unit cell is the smallest structure that is more consistent with the sample structure and maintains the physical properties of threonine. The results showed that the use of THz spectroscopy could achieve rapid and accurate identification of L-thr and a better theoretical interpretation of its low-frequency vibration mode. Therefore, to understand how the configuration of the hydrogen bond network affects the low-frequency vibra- tions of the diastereoisomers L-thr and L-allo-thr, we further used THz spectroscopy and solid-state DFT calculations to compare their low-frequency vibration characteristics. For this study, it was helpful to comprehensively understand the diastereomer system of L-thr-L-allo-thr.

Currently, studies on binary solid solution systems formed by different amino acids roughly includes L-valine–L-iso-leucine (L-val-L-ile) (13,14), L-val-L-leucine (L-leu) (15), L-ile-L-leu and L-ala-L-ser (16), and L-thr-L-allo-thr (1,4), where the similarity of sizes and shapes of these amino acid molecules indicates that they can form solid solutions in their mixtures. However, most of them are nonequimolar discrete heterocompounds and limited solid solutions in the vicinities of pure components, whereas in the L-thr-L-allo-thr study system exists the phenomenon of the continuous isomorphic miscibility of the components, which is rarely found in chiral systems. In addition, in the synthesis or purification process of L-amino acid products, many impurities are often encountered, including pure compounds, physical mixtures, intermediate compounds, or co-crystallized mixtures. Therefore, a fast and accurate method is needed to determine whether a pure monomer or eutectic is formed (4,8,17). The THz absorption spectra is sensitive to the structural changes of diastereomers, taking L-thr and L-allo-thr as the research objects and analyzing the THz absorption of threonine and allo-threonine. The low-frequency vibration mode of them is calculated by using DFT, and the effects of intermolecular hydrogen bonds on different cell configurations and various low-frequency vibration modes in amino acid molecules are analyzed. To understand this process, it is helpful to further understand the molecular interaction of the co-crystallized mixtures L-thr-L-allo-thr.

Experimental Methods and Sample Preparation 

This experiment used the transmission TDS system produced by Daheng Optoelectronics. The THz emission source was a light guide antenna, and the THz detector was zinc telluride (ZnTe). The center wavelength was 800 nm, the laser pulse width was 100 fs, and the repetition frequency was set at a 84 MHz pulse. The effective THz spectral range was 0.2–3.0 THz. The sample was placed in a sample chamber filled with nitrogen, and the relative humidity was less than 1%, which was measured at room temperature.

A powder X-ray diffraction (PXRD) system from ThermoARL was used in this experiment. The working voltage was 40 kV. The working current was 40 mA. The scanning range was 5–50°. The scanning speed was 0.2°/step.

IR spectroscopy was performed using a Nicolet 5700 Fourier transform IR (FT- IR) spectrometer under the conditions of a pressure of 15 MPa and a temperature of room temperature. The scanning range was 400–4000 cm

, the number of scanning was 32 times, and the resolution was 4 cm

L-thr (CAS: 72-19-5, purity >99%) and L-allo-thr (CAS: 1509-34-08, purity >98%) were purchased from Aladdin. The samples used in the experiment were not purified further. The sample and polyethylene (PE) were mixed thoroughly at a ratio of 1:5. Then, the mixture was compressed into pellets with a diameter of 13 mm and a thickness of approximately 1 mm at room temperature under a pressure of 10 MPa for 3 min.

The theoretical calculations in this study were all performed using Castep-20.11 software (18). The generalized gradient approximation method (GGA) based on the plane wave was used for the calculation, and the exchange correlation functional used was the typical Perdew-Burke-Ernzerhof (PBE) functional method. The mass calculated by the energy was ultrafine, and the plane-wave cutoff energy was 750 eV. The total energy converges to 5.0 × 10

 eV/atom, and the maximum forces between atoms were less than 0.01 eV/Å. The crystal structure data of L-thr (identifier: 128538) and L-allo-thr (identifier: 1061093) come from the Cambridge Crystallographic Data Center (CCDC). The space group of L-thr and L-allo-thr is p212121 (

 = 4), and the corresponding lattice parameters are 

 = 5.110 Å, α, β, γ = 90°(19,20), and 

 = 5.14 Å, α, β, γ = 90° (4,21), respectively. The crystal structures of L-thr and L-allo- thr are shown in Figure 1.

FIGURE 1: The molecular structure of (a) L-thr (a) and (b) L-allo-thr.

To confirm the crystal structure of the amino acid samples in this experiment, PXRD technology was used, and the results are shown in Figure 2. The measured main peak positions are in good agreement with the spectra matched by the Jade 6.0 standard PDF library (L-thr (PDF #10-0751), L-allo-thr (PDF #05-0300)). The agreement shows that the quality of the samples used in this experiment was higher. The XRD spectrum of L-thr was different from that of L-allo-thr, which was caused by the different crystal structures of the two amino acids.

FIGURE 2: XRD patterns of (a) L-thr and (b) L-allo-thr.

THz Spectral Analysis and Cell Simulation Calculation of L-thr and L-allo-thr

PE has almost no absorption of the THz waves and can be used as the reference. The time domain spectra and frequency domain spectra of the reference and experimental samples are shown in the Figure 3a and 3b. Compared with the pure PE signal, the THz time-domain signal of samples containing amino acids have time delay and amplitude attenuation. The time delay is mainly because the refractive index of the THz wave passing through sample is larger than the reference, and the time delay increases with the increase of sample thickness. Meanwhile, the attenuation of signal intensity is caused by the absorption and random scattering of the THz wave by the samples. Using Matlab software, fast Fourier transform was performed on the THz time-domain spectra of the samples, and the THz frequency domain spectra can be obtained. Compared with PE, the frequency domain spectra of the L-thr and L-allo-thr samples not only has corresponding amplitude attenuation, but also has several distinct and unique characteristic absorption bands.

Figure 3 also shows the experimental and theoretical THz absorption spectra of L-thr (Figure 3c) and L-allo-thr (Figure 3d), respectively, with significant differences in the frequency range of 1.0–2.3 THz. The characteristic peaks of L-thr are located at 1.42 and 2.14 THz (11,12,22,23), and the characteristic peaks of L-allo-thr are located at 1.63 and 2.16 THz. It shows that L-thr and L-allo-thr can be distinguished by a THz fingerprint, which is consistent with the THz spectral results of L-thr in other literatures. Currently, there are no reports on the spectral studies of L-allo-thr. In addition, theoretical simulations can link the calculated vibration patterns to the molecular structure. The experimental and theoretical characteristic peaks of these two amino acids obtained by DFT crystal structure calculation based on GGA-PBE are listed in Table I. The theoretical characteristic peaks of L-thr are 1.56 and 2.16 THz, which corresponded to the experimental results of 1.42 and 2.14 THz. The theoretical characteristic peaks of L-allo-thr are 1.87 and 2.22 THz, corresponding to the experimental results of 1.63 and 2.16 THz. The two amino acids simulation peaks are blue shifted relative to the experimental characteristic peaks, and the largest error is that L-allo-thr is blue shifted 0.24 THz at 1.63 THz. According to the comparison of XRD data of L-thr and L-allo-thr and the optimized intermolecular hydrogen bonds in Table II, it can be seen that the intermolecular hydrogen bond lengths in their crystal cells are shortened. These shortened bond lengths further resulted in the increase of vibration frequency, the change of electron orbital density of amino acid molecules, and the internal structure reorganization. Thus, blue shift would be triggered (24,25). On the other hand, the DFT calculates the absorption frequency at the temperature of 0 K, whereas the experimental spectra are all measured at room temperature, and the difference of temperature will also cause blue shift (26–28).

FIGURE 3: (a) Time domain spectra and (b) frequency domain spectra of reference PE and amino acid samples. Experimental and theoretical THz absorption spectra of (c) L-thr and (d) L-allo-thr.

Hydrogen Bond Network Structure Analysis and Vibration Mode Attribution of L-thr and L-allo-thr 

Figure 4 shows the hydrogen bond network structure in unit cell of L-thr (Figure 4a) and L-allo-thr (Figure 4b), respectively. L-thr and L-allo-thr have the same atoms and number of atoms, forming three identical hydrogen bond types, namely O

, respectively. The threonine molecules in the crystal exist in the form of zwitterions, where two hydrogen bonds come from NH

 and one hydrogen bond from COO- (29). Both of them are all hydrogen bond networks formed between carbonyl, carboxyl, and amino groups through O···H-N and O···H-O, connecting with neighboring molecules (7,30–32). However, L-thr and L-allo-thr are diastereomers of each other, and there are certain differences in the hydrogen bond interaction between them. Between cells, L-thr requires three threonine molecules to form a three-dimensional (3D) hydrogen bond network, whereas L-allo-thr requires four allo-threonine molecules.

FIGURE 4: The crystal structure of unit cells of (a) L-thr and (b) L-allo-thr.

Figure 5 shows the calculated crystal structure vibration modes of L-thr at 1.56 (Figure 5a) and 2.16 (Figure 5b), and the structural modes of THz and L-allo-thr at 1.87 (Figure 5c) and 2.22 (Figure 5d) THz, respectively. A variety of vibration modes exist in the 1–2.3 THz region mainly because of the hydrogen bond interaction, leading to the specific structure in the crystal and resulting in the specific vibration behavior of different characteristic peaks in the THz frequency band (6). There are certain differences in the location of the THz characteristic peak and corresponding vibration modes between L-thr and L-allo-thr, which are mainly attributed to the deviation of their cell parameters, including the length of hydrogen bonds and the volume of the unit cell (33). It can be seen from Table I that the volume of the unit cell of L-allo-thr after optimization is much larger than that of L-thr. The reduction in the distance between the nitrogen and oxygen atoms in adjacent molecules in a unit cell leads to the enhancement of intermolecular interactions, especially the short-range force of hydrogen bonds. Small differences in the structure of molecules will cause changes in the hydrogen bond energy and bond length (7,33,34). Therefore, compared with L-allo-thr, the crystal structure of L-thr is denser, the noncovalent interaction is stronger, and the intermolecular motion is stronger. The detailed description of the calculated vibration modes of L-thr and L- allo-thr is shown in Table I. At 1.56 THz, L-thr is reflected in the swinging vibration of the amino and carboxyl groups connected with hydrogen bonds along the z-axis. The vibration mode at 2.16 THz was the carboxyl group in-plane swing vibration and the collective swing vibration of C

, whereas the vibration mode of L-allo-thr at 1.87 THz was mainly because of the collective swing vibration of all atoms along the 

-axis. The absorption peak at 2.22 THz was the rocking vibration in the base plane of R and the collective swing vibration of C

. It can also be seen from Figure 5 that the group connected with hydrogen bond in the crystal cell vibrates significantly. Therefore, within the effective THz frequency range, the vibration mode between L-thr and L-allo-thr molecules were mainly swing vibration modes centered on intermolecular hydrogen bonds.

FIGURE 5: The vibration modes of the crystal structure calculated by L-thr at (a) 1.56 THz and (b) 2.16 THz, and L-allo-thr at (c) 1.87 THz and (d) 2.22 THz, respectively.

Qualitative and Quantitative Analysis of the Co-crystallized Mixtures of L-thr and L-allo-thr

The L-thr and L-allo-thr powders were mixed in a fixed ratio and completely dissolved in deionized water, recrystallized in a vacuum oven at 50 °C, and formed a new phase. In addition, a 1:1 physical mixed tablet of L-thr and L-allo-thr was prepared. The two solid drugs were simply mixed together and pressed into tablets. To further determine the formation of co-crystallized mixtures, we also used PXRD technology to characterize the structure of these samples. The crystal structure changes after co-crystallized can be clearly seen from Figure 6a. The diffraction patterns of L-thr and L-allo-thr in Figure 6 have similar configurations; for the physical mixtures, the diffraction patterns contain most of the diffraction peaks of the two groups of pure enantiomers. For co-crystallized mixtures, their diffraction peaks tend to be wider and not as sharp as those measured by the pure components. Moreover, the total diffraction peaks are gradually biased to the pure L-allo-thr with the increase of the L-allo-thr component. The results show that L-thr and L-allo-thr form co-crystallized structures through noncovalent bonding successfully, especially driven by strong intermolecular hydrogen bonds, whereas the mixtures do not involve such hydrogen bonds or other noncovalent interactions (35).

FIGURE 6: (a) The XRD patterns of co-crystallized mixtures and physical mixture, (b) THz spectra of L-thr, L-allo-thr and their physical mixture, (c) THz spectra of cocrystallized mixtures (L-thr, 25%, 50%, 75%, L-allo-thr), (d) the correlation between the frequency offset of the co-crystallized mixtures’ characteristic peaks (Δf, relative to the characteristic peaks at 1.42 and 2.14 THz of L-thr at zero) and the composition ratio of L-allo-thr.

Figure 6b shows the THz absorption spectra of L-thr, L-allo-thr, and their 1:1 physical mixture. Figure 6c shows the THz absorption spectra of co-crystallized mixtures with different proportions. The THz spectra of the physical mixture has three characteristic peaks: one of the absorption peaks at 2.14 THz was relatively strong, and the other two absorption peaks at 1.42 THz and 1.63 THz were relatively weak because the spectra at 2.14 THz was affected by the common absorption of two amino acids (36). As a result, the absorbance was slightly improved. The results showed that the measured absorption peaks of the physical mixture corresponded to the experimental absorption peaks of the two pure diastereomeric components L-thr and L-allo-thr. However, the THz spectra of the eutectics contained two characteristic peaks. One of the weaker absorption peaks was located at 2.14 THz, and the other characteristic peaks of co-crystallized mixtures were respectively located at 1.44, 1.49, and 1.57 THz for 25%, 50% and 75% L-allo-thr. With the gradual increase of L-allo-thr components in co-crystal- lized mixtures, their characteristic peaks are closer to those of pure L-allo-thr (4). The significant difference here can be used for identifying co-crystallized mixtures in different proportions. The co-crystallized mixtures with different crystal structures have their own characteristic peaks, and the difference in the characteristic absorption peaks of the 1:1 physical mixture and co-crystallized mixtures was because of the THz experimental spectra caused by the collective low frequency vibration mode of the molecules (16).

In addition, Figure 6d is a linear fitting between the frequency offset of the co-crystallized mixtures characteristic peaks (Δf) and the composition ratio of L-allo-thr, which was relative to the characteristic peaks of L-thr at 1.42 and 2.14 THz. The results show a blue–red straight line, in which the fitting equation of Δf

 = 0.955 (1), and the fitting equation of Δf

 represents the relative frequency offset, x represents the composition ratio of L-allo-thr, 

 is the correlation coefficient of linear fitting, where the 

 of the red line, indicating that the fitted line between Δf

 and the composition ratio of L-allo-thr is in good agreement with the experimental data because with the increase of the absorption coefficient of the THz spectrum in the high-frequency region, the signal-to-noise (S/N) ratio decreases. This relationship enabled us to quickly determine the composition ratio of L-thr and L-allo-thr in co-crystallized mixtures according to the frequency offset of characteristic peaks (8,17,37,38).

The absorption spectra of the L-thr and L- allo-thr diastereomers were studied by THz-TDS, and the results calculated by DFT were in good agreement with the measured characteristic peaks. The cell configurations of L-thr, L-allo-thr, and the differences between various low-frequency vibration modes were analyzed from the perspective of hydrogen bond configurations. It was found that the hydrogen bond between the amino group, carboxyl group, and the side hydroxyl group between the cell molecules was the main cause of the molecular skeleton vibration. The THz spectra of the 1:1 physical mixture and co-crystallized mixtures of L-thr and L-allo-thr showed that the THz absorption characteristics of the physical mixture included the absorption characteristics of the two pure diastereomers, and the co-crystallized mixtures comprised of diastereomers in different proportions had their own distinct characteristic peaks. It was indicated that the THz-TDS technique provided a feasible method for the quantitative analysis of the co-crystallized mixtures formed by diastereomers L-thr and L-allo-thr.

(1) N. Taratin, H. Lorenz, D. Binev, et al, 

(2) H. Zhao, K. Hamase, A. Morikawa, et al, 

J. Chromatogr. B Analyt Technol. Biomed. Life Sci.

(4) D. Binev, N. Taratin, E. Kotelnikova, et al., 

(5) E.N. Kotelnikova, A.I. Isakov, L.Y. Kryuchkova, and A.A. Zolotarev, in 

Processes and Phenomena on the Boundary between Biogenic and Abiogenic Nature

 (Springer, Berlin, Germany, 2020), pp. 695–719.

(7) L.-j. Huang, X. Zhang, G. Wang, et al, 

(13) A.I. Isakov, E.N. Kotelnikova, S. Muenzberg, et al, 

(14) E.N. Kotelnikova, A.I. Isakov, and H. Lorenz, 

(15) A.I. Isakov, H. Lorenz, A.A. Zolotarev, Jr., et al, 

(16) E. Kotelnikova, R. Sadovnichii, L. Kryuchkova, et al, 

(17) M. Yamaguchi, F. Miyamaru, K. Yamamoto, et al, 

(18) S.J. Clark, M.D. Segall, C.J. Pickard et al, 

Acta Crystallogr. C Cryst. Struct. Commun. 

(22) W.-N. Wang, H.-Q. Li, Y. Zhang, et al, 

(23) Z.-f. Tang, H.-t. Lin, X.-w. Chen, et al, 

(27) J. Liang, X. Zhang, N. Wang, et al, 

Spectrochim. Acta A Mol. Biomol. Spectrosc.

(29) M.T. Ruggiero, J. Sibik, J.A. Zeitler, et al, 

Spectrochim. Acta A Mol. Biomol. Spectrosc. 

(31) T.B. Issa, F. Sayari, H. Ghalla, et al, 

(36) F. Huang, S. Zhuang, W. Liu, et al, 

Spectrochim. Acta. A Mol. Biomol. Spectrosc.

 are with the College of Materials and Textiles at Zhejiang Sci-Tech University, in Zhejiang, China. 

 is also with the Key Laboratory of Advanced Textile Materials and Manufacturing Technology of Ministry of Education at Zhejiang Sci–Tech University, in Zhejiang, China. Direct correspondence to: 

Studying the experimental absorption spectra of different proportions co-crystallized mixtures and 1:1 physical mixture of L-thr and L-allo-thr, it was found that the characteristic peaks of the physical mixture include the characteristic peaks of the two diastereomers in the THz band, while amino acid co-crystallized mixtures formed their own characteristic peaks depending on the proportion.

Terahertz Spectroscopic Analysis of Co-Crystallized Mixtures in an L-threonine Diastereomer System

The hydrogen bonding mode of 5-fluorouracil-proline co-crystal (5F-P) was investigated using terahertz time-domain spectroscopy (THz-TDS) and density functional theory (DFT). The THz spectrum of 5F-P is significantly different from the monomer of 5-fluorouracil (5- FU) and proline within the range of 0.25–2.6 THz. This result indicates that the two monomers react to form a new phase. Based on unit cell simulation, the theoretical spectrum of configuration 2 (N1–H1•••O4) is the most consistent with the experimental spectrum. Furthermore, configuration 2 carries the lowest amount of total energy and exhibits the most stable structure. Therefore, the hydrogen bonding mode of configuration 2 is closest to 5F-P. This work shows that the combination of THz-TDS and DFT can provide an effective method to investigate the hydrogen bonding mode of the new phase.

L-proline (L-Pro) has a unique structure among the 20 amino acids normally found in proteins because it is the only imino acid among them. The unique cyclic structure leads to the conformational rigidity and specific imino functionality compared to amino acids (1,2). The amino nitrogen binds with the side chain, with the α-carbon leading to a pyrrolidine ring. This structure restricts the conformations that proline can adopt within a peptide or protein, giving it a unique role in the secondary and tertiary structures of proline-containing proteins (3,4). Given this conformation, proline is often found in complex proteins. It plays a binding role in collagen, which has a high proline content. Proline-rich regions in proteins are also found in many kinase binding sites and are linked to many other cellular processes (5–8). Proline is an amphoteric substance with both carboxyl and amino groups. Its solubility in aqueous solutions is more soluble than other amino acids. Proline has also been shown to have strong water absorption in inelastic incoherent neutron scattering studies (9). Proline contains an imino (-NH-) group and a carboxyl (-COOH) group, which can form both intramolecular hydrogen bonds and intermolecular hydrogen bonds with other molecules (10).

In 1957, Heidelberger and others (11) first discovered 5-fluorouracil (5-FU). They replaced the hydrogen atom at the fifth position of uracil with a fluorine atom of similar size. The obtained fluoride is not only similar in volume to the original compound but also the carbon–fluorine bond formed is very stable, which is not easy to decompose during metabolism and can interfere with normal metabolism at the molecular level. 5-FU is an antimetabolic tumor drug, and its active metabolite inhibits thymidine synthase (TS), which is the primary mechanism for its antitumor effect (12). However, its actual efficacy is affected to a certain extent because of the low solubility of 5-FU. Because proline is nontoxic and soluble in water, it is an ideal co-crystal former (CCF). It can effectively improve the solubility of 5-FU by co-crystallization of 5-FU and proline (5F-P). However, its toxicity and anticancer need further research in medicine. The molecular structure of proline and 5-FU is shown in Figure 1.

FIGURE 1: Molecular structure of (a) proline and (b) 5-FU.

The drug co-crystal refers to a crystal formed by combining active pharmaceutical ingredient (API) and CCF in a particular proportion under the action of hydrogen bonds or other non-covalent bonds (13–15). 5-FU, as an active pharmaceutical ingredient, combines with proline to form a co-crystal through hydrogen bonding. The hydrogen bond is usually expressed as A—H•••B. Although hydrogen bonding is a weak interaction, it plays a vital role in affecting biomolecules’ physical and chemical properties. In recent years, researchers have done a lot of research on noncovalent interactions among molecules, including experimental methods such as infrared, Raman, and nuclear magnetic resonance spectroscopy, as well as calculation methods such as molecular dynamics simulation (16–18).

As a novel technique, terahertz (THz) spectroscopy is regarded as an advanced approach that studies new aspects of molecular structure and distinguishes polymorphic pharmaceutical, which is especially useful when studying biomolecules, which contain rich hydrogen-bond interactions (19–22). THz spectroscopy could directly collect helpful information on low-frequency vibrations shown within specific molecular structures and intermolecular hydrogen bonding interactions. Meanwhile, computational methods, such as solid-state density functional theory (DFT), have shown potential as reliable approaches to studying these low-energy lattice vibrational motions (23–26). Thus, the combination of DFT and experimental data is a very effective method to study the molecular geometry and intermolecular interactions. Applications of terahertz time-domain spectroscopy (THz-TDS) on weak interaction and co-crystal have been reported in previous studies. Sean P. Delaney and others demonstrated that the co-crystal of flufenamic acid and nicotinamide has better overall binding energy by THz-TDS and solid density (27). Chihiro Funaki and others proved that the existence of three different kinds of weak C-H•••O=C interactions in poly (ε-caprolactone) by THz-TDS (28). Qi Zhou and others characterized the co-crystal of carbamazepine with nicotinamide, saccharin, and fumaric acid by THz-TDS and DFT (29).

In this study, we successfully obtained the co-crystal of proline and 5-FU by evaporating the solvent. As a new phase, 5F-P co-crystal has no detailed report on its hydrogen bonding mode, which is significant for studying the properties of substances. Therefore, we investigated the hydrogen bonding mode of 5F-P by using THz-TDS and DFT. Furthermore, we also verified the accuracy of the hydrogen bonding model of 5F-P through the stability of the system.

L-Pro (purity > 98%), 5-FU (purity > 98%) and polyethylene (PE, particles size 40–48 μm) were purchased from Sigma-Aldrich. All chemicals were used without further purification.

Proline and 5-FU were placed in a 100 mL beaker at a 1:1 molar ratio (690 mg:780 mg) and added 70 mL of a methanol:ethyl and an acetate:acetone mixed solvent in a volume ratio of 1:1:1 to the mixture. It was then heated and stirred in a constant temperature water bath for about 2 h at 60 °C. Finally, the obtained solution was filtered under reduced pressure to remove insoluble matters to get clear solution. It was slowly evaporated at normal temperature for about three days to collect the solid phase. The obtained solid phase was dried overnight to obtain 5F-P in a constant temperature environment at 25 °C.

Samples and PE were weighed into agate mortar at 1:5 (weight of sample/weight of PE = 40 mg/200 mg) and ground gently into fine particles with a pestle to thoroughly mix samples with PE. The samples were then compressed into a pellet (about 1.9 mm in thickness and 13 mm in diameter) with a hydraulic press using 10 tons pressure for 3 min.

The THz-TDS used in the experiment was made by Daheng Optoelectronics. The center wavelength of the CIP-TDS spectrometer was 800 nm, and the repetition frequency was 84 MHz. The effective spectral band was 0.2–3 THz. The samples were placed into the THz system, purged with high-purity nitrogen and kept at the relative humidity <1%. The THz spectrum of each sample was measured at room temperature (20 °C), and PE was taken as a reference. Each spectrum was an average of three measurements. One measurement included 1024 scan numbers. The time-domain of the THz electric field was recorded for the reference and each sample. The fast Fourier-transform (FFT) operation obtained the time-domain spectrum, and the THz spectral absorption was obtained by dividing the sample frequency response by that of the reference.

The PXRD was purchased from Thermo ARL, Switzerland, for this experiment. The X-ray source was a copper shoe, and the working voltage was 40 kV. The working current was 40 mA. The scanning range was 5 to 50°, and the scanning speed was 0.2°/step.

The theoretical calculations in this study were all performed using the Cambridge Sequential Total Energy Package (CASTEP) (30). There was no virtual frequency in all simulation calculations, which indicates that the optimized structure achieved a stable structure corresponding to the lowest molecular energy. The exchange- related functional used generalized gradient approximation–Perdew-Burke-Ernzerhof (GGA–PBE) (31). The quality of the energy calculation was ultra-fine. The quality of the plane wave cut off energy was ultra-fine, and the energy cutoff was 830 eV. Norm conserving was selected as the pseudopotential. The quality of the 

-point set was fine. The convergence thresholds for energy, maximum force, stress, and maximum displacements was 5.0 × 10

 (eV/atom), 0.01 eV/Å, 0.02 GPa and 5.0 × 10

 Å, respectively. The Grimme correction was also added to the calculation task to evaluate weak interactions, such as hydrogen bonds among molecules (32).

Powder X-ray Diffraction (PXRD) of Proline, 5-FU and 5F-P 

Figures 2a, 2b, and 2c shows the PXRD patterns of proline, 5-FU and 5F-P. Figure 2d is the simulated PXRD based on the unit cell of 5F-P. The diffraction peaks of proline are located at 15.28°, 18.2°, 19.72°, 24.94°, 30.76°, 32.34°, and 39.96°, respectively. The diffraction peaks of 5-FU are located at 11.22°, 15.7°, 16.12°, 18.83°, 21.69°, 28.43°, 32.69°, and 36.95°. The results are consistent with previous reports (33,34). The PXRD of 5F-P is significantly different from the two monomers. The comparison of PXRD with the Materials Data Inc. (MDI) Jade 6.5 PDF standard card showed that the proline is orthorhombic with space group P212121, a = 11.626 Å, b = 9.024 Å, c = 5.246 Å, α = β = γ = 90° (PDF #21-1805), but the 5-FU is triclinic with space group P-1, a = 9.22Å, b=12.66Å, c=12.67Å, α=89.7°, β = 43.9°, γ = 98.6° (PDF #39-1860). Since 5F-P is a new phase, the crystal structure does not appear in the Cambridge Crystallographic Data Center (CCDC). We used the Reflex module to index the co-crystal based on the PXRD diffraction peak. The result after further refinement is that the 5F-P is triclinic, the space group is P-1 with a = 11.13 Å, b = 8.12 Å, c = 6.21 Å, α = 88.65°, β = 104.89°, γ = 107.43°. As shown in Figure 2d, the simulation PXRD pattern based on indexed data is consistent with the experimental pattern of 5F-P. Therefore, these results show that the simulated unit cell of 5F-P is reasonable.

FIGURE 2: Powder X-ray diffraction (PXRD) pattern of (a) proline, (b) 5-FU, (c) 5F-P and (d) the simulated PXRD pattern of 5F-P.

THz Absorption Spectra of Proline, 5-FU, and 5F-P

The THz experimental spectrum of monomers, co-crystal and mixture (mixed at 1: 1 molar ratio) are shown in Figure 3. The THz spectrum of proline exhibits a noticeable absorption peak at 1.99 THz, while 5-FU has a strong absorption at 2.46 THz, as shown in Figure 3a. There are few absorption peaks of proline and 5-FU. From Figure 3b, it can be seen that the THz absorption spectra of 5F-P display four visible absorption peaks, which are 0.98, 1.32, 1.55, and 1.99 THz, respectively. Furthermore, there is a weak absorption peak observed at 1.11 THz, which gives rise to a shoulder peak with the absorption peak at 0.98 THz. The THz spectrum of the mixture shows a weak absorption band at 1.99 THz. Besides, there is a suspicious absorption peak at 2.46 THz. It can be demonstrated that the absorption peak of the mixture is a simple superposition of two monomers. The THz absorption spectra of 5F-P is significantly different from the monomers and their mixture. Therefore, both the THz spectrum and PXRD show that the new solid phase is a co-crystal of the two reactants rather than a simple polymorphic transition of the single crystal. The crystal structure of the co-crystal has changed to some extent. Furthermore, a series of new absorption peaks appeared after the co-crystallization, which is attributed to the hydrogen bonding between the two molecules. Because the THz absorption peak of 5F-P is much more than two monomers, 5F-P produces more vibration modes. The formation of hydrogen bonds changes the arrangement of electrons, while the THz band is highly sensitive to hydrogen bonding interactions and collective vibrations of molecular lattices (35,36), thus resulting in corresponding changes in vibration intensity and frequency. The new absorption vibration peaks indicate that new hydrogen bonding interactions are produced in the co-crystal. Apart from that, the THz spectrum also can intuitively reflect the occurrence of intermolecular vibration.

FIGURE 3: The THz absorption spectra of (a) Pro, 5-FU. (b) The THz absorption spectra of co-crystal and mixture.

Natural Bond Orbital (NBO) Calculation of Charge Distribution

First of all, we have investigated potential intermolecular hydrogen bonding sites of 5F-P by the NBO calculation. NBO analysis is an effective tool for studying charge transfer, electronic transitions, and conjugate interactions. NBO atomic charge distribution was calculated at the B3LYP/6-311++G (d,p) level. Figure 4 shows the charge distribution of the stronger electronegative atoms in the proline (a) and 5-FU (b). There are two stronger polarity and negatively charged oxygen atoms O1 (-0.611e) and O2 (-0.694e) observed on the carboxyl group of proline. In addition, there is a stronger electronegative nitrogen atom N1 (-0.675e) spotted on the amino group. Also, there are strongly electronegative atoms on 5-FU, including two carbonyl oxygen atoms O1 (-0.564e), O2 (-0.612e) and two nitrogen atoms N1 (-0.650e), N2 (-0.605e). According to previous reports, both amino and carboxyl groups are more active groups. The highly electronegative atoms present on proline and 5-FU cause an uneven charge distribution. Therefore, there is the possibility of hydrogen bond formation between proline and 5-FU. In addition, proline is an effective co-crystal former. All these indicate that proline can easily form hydrogen bonds with 5-FU.

FIGURE 4: The atomic charge distribution of (a) Pro and (b) 5-FU. (The blue dashed line is hydrogen bonding). The various atoms are identified in the figure.

Table I shows the intermolecular hydrogen bond network of proline and 5-FU through previous reports (37,38). It can be seen that the imino and carboxyl groups on proline form intermolecular hydrogen bonds. The carbonyl oxygen atoms and imino groups on 5-FU participate in forming the intermolecular hydrogen bond network. This result is consistent with the NBO analysis in section 4.3. Therefore, the theoretical configuration is based on the formation of hydrogen bonds between strongly electronegative groups. In this study, 11 reasonable hydrogen bonding modes were obtained. The theoretical configuration of 5F-P was constructed based on these modes. We further optimized the unit cell of 11 theoretical configurations, two of which failed to achieve convergence. Therefore, we have obtained nine stable theoretical configurations of 5F-P and their hydrogen bond modes are shown in Figure 5. The hydrogen bonding modes of configurations 1–9 are N2–H2•••O4 (C1), N1–H1•••O4 (C2), N2–H2•••O3 (C3), N1–H1•••O3 (C4), N3–H12•••O2 (C5), N3–H12•••O1 (C6), N1–H1•••O3/ N3–H11•••O1 (C7), N1–H1•••O3/N3– H11•••O2 (C8), and N2–H2•••O3/N3– H11•••O2 (C9).

FIGURE 5: The hydrogen bonding mode of theoretical configuration (a) 1, (b) 2, (c) 3, (d) 4, (e) 5, (f) 6, (g) 7, (h) 8, (i) 9. The color code for atoms is given in the figure.

We have calculated the theoretical spectra based on nine different unit cells, as shown in Figure 6. There are four obvious absorption peaks in the theoretical spectrum of C1, which are 0.86, 1.16, 2.09, and 2.21 THz, respectively. There are six absorption peaks in the spectrum of C2, of which four absorption peaks are ap- parent, and two characteristic absorption peaks are weak. C3 has five obvious characteristic absorption peaks, located at 0.58, 1.15, 1.5, 1.99, and 2.35 THz. Compared with the other eight configurations, the theoretical absorption peak of C4 is the most, and there are seven absorption peaks. The absorption peaks of C5 are located at 0.55, 1, 1.33, 1.6, and 2.08 THz. There are four weak ab- sorption peaks in the theoretical spectrum of C6, which include one relatively strong absorption peak located at 2.28 THz. The number of theoretical absorption peaks of C7, C8, and C9 is 6, 5, and 5. There are significant differences in the theoretical spectra of different hydrogen bonding modes. This is because different hydrogen bond modes lead to different strengths of hydrogen bonds in the unit cell. The hydrogen bond itself is a noncovalent interaction. The stronger the noncovalent interaction, the more intense the intermolecular motion-induced (39).

FIGURE 6: The THz theoretical spectra of (a) to (i) as C1 to C9, respectively.

Moreover, the matching degree between theoretical and experimental spectra calculated by different hydrogen bonding modes are significantly different. As shown in Table II, the corresponding absorption peak numbers of the nine theoretical spectra are 3, 5, 3, 3, 4, 3, 4, 3, and 4, respectively. Based on the number of matching absorption peaks and the shifting degree of absorption peaks, it can be concluded that the hydrogen bonding modes of C2, C5, C7, and C9 are closer to the actual situation. It is worth noting that the N1–H1•••O4 hydrogen bond mode in C2 is theoretically the most consistent with the hydrogen bond of 5F-P.

The four systems of C2, C5, C7, and C9 have atoms. The covalent bonds in the unit cell are also the same. The difference among them lies in the mode of hydrogen bonding. Because the hydrogen bond is only a non-covalent interaction and not a chemical bond, its strength cannot affect the overall stability. Besides, a complete system always tends to the lowest energy because the lower the total energy of a substance, the weaker its ability to react. This substance can continue to exist in the most stable state. Therefore, the more stable a substance is, the lower the total energy it carries (40). In this study, the total energy of the four better configurations (C2, C5, C7, and C9) was calculated by PBE+Grimme to compare their stability and further verify the accuracy of the conclusion in section 4.5. Table III shows that the total energy is 2 < 7 < 5 < 9 before correction, and the total energy after dispersion correction is the same as 2 < 7 < 5 < 9. Therefore, the system with C2 as the unit cell is the most stable. The intermolecular hydrogen bonding mode of 5F-P is most inclined to C2. This result is consistent with the results of THz spectroscopy.

The comparison of the theoretical and experimental spectrum of 5F-P, as shown in Figure 7. The absorption peaks of 0.98, 1.11, 1.32, 1.55, and 1.99 THz in the experimental spectrum correspond to 1.01, 1.13, 1.35, 1.55, and 2.2 THz of the theoretical spectrum. It shows that the theoretical spectrum calculated by C2 simulates all the experimental absorption peaks of 5F-P. However, the theoretical absorption peak at 0.63 THz does not appear in the experimental spectrum. The reason is that the effect of water molecules cannot be eliminated entirely from the actual process and the low signal-to-noise ratio (S/N) of the equipment. As shown by the black arrow, the absorption peak at too low or too high frequency is annihilated by noise, making it difficult to detect accurately. Except for the absorption peak at 1.55 THz, other theoretical absorption peaks are blue-shifted to varying degrees relative to the experimental absorption peak. The reasons can be primarily attributed to the difference between the experimental (293 K) and the simulation temperature (0 K) (41,42). The change in temperature affects the size of the unit cell (43). The peak position and peak width of the theoretical spectrum will change accordingly. Because of the thermal expansion effect, the intermolecular distance increases with the temperature increasing, as does the bond length. The thermal expansion of the crystal lattice weakens the strength of the intermolecular interactions, resulting in the red-shift of the absorption peak (44,45). On the contrary, the absorption peak is blue-shifted. Therefore, the temperature difference affects unit cell simulation. It makes the frequency position of each absorption peak fluctuate to a certain extent, which is directly manifested as the theoretical absorption peak of C2 moving to a higher frequency. Some vibration modes exhibit an obvious red-shift with the increase of temperature, while others are less affected by temperature (46). The theoretical absorption peak at 1.55 THz does not have a blue-shift. This may be because the absorption peak originates from intramolecular and intermolecular vibrations in the THz band. The weak interactions have little effect on the vibration mode at this frequency. Therefore, the change in temperature cannot significantly affect the position of this absorption peak. In addition, another reason is the difference between actual and theoretical crystal structure. The simulation is based on the perfect crystal structure. However, it is difficult to reach the ideal degree in the experiment (47–49). The temperature effect and sample difference make the crystal parameters of experiment and theory not match exactly. It is manifested directly in the difference between experimental and theoretical spectra.

FIGURE 7: The THz experimental spectra of (green) 5F-P, and (red) the theoretical spectra of C2.

Figure 8 shows the vibration modes of the theoretical absorption peaks in the unit cell, corresponding to 0.98, 1.11, 1.32, 1.55, and 1.99 THz of the experimental absorption peaks, respectively. The vibration of C2 at 1.01 THz is mainly caused by the librational vibration of 5-FU along the b-axis. The vibration at 1.13 THz is attributed to the strong reverse translational vibration of Pro and 5-FU. The vibration at 1.35 THz is similar to 1.01 THz, mainly caused by the librational vibration of 5-FU along the b-axis, but the one-sided translational vibration of the five-membered ring on proline at 1.35 THz is more obvious. The absorption peak at 1.55 THz, assigned in rotational vibration of the five-membered ring on proline. In addition, it can be seen that the hydrogen bond vibration of N1–H1•••O4 between proline and 5-FU is close to zero. This result proves that the theoretical absorption peak at 1.55 THz does not blue-shift because the weak interaction has little effect at this position. The vibration at 2.2 THz is caused by the translational vibration of proline and unilateral 5-FU in the same direction.

FIGURE 8: The vibration mode of C2 at (a) 1.01, (b) 1.13, (c) 1.35, (d) 1.55, and (e) 2.2. THz (Black arrows represent translational vibration, red arrows represent librational vibration and blue arrows represent rotational vibration). The color code for atoms is given in the figure.

In this work, an investigation was conducted into the hydrogen bonding modes of 5F-P by combining THz-TDS and DFT. We found that the THz spectrum of 5F-P in the range from 0.25 to 2.6 THz is significantly different from the two monomers. Based on NBO analysis and hydrogen bond network structure, we constructed nine stable theoretical configurations of 5F-P. The THz theoretical spectra of differ- ent hydrogen bond modes are obtained by unit cell simulation, and the theoretical spectra of configuration 2 is the most consistent with the experimental spectra. The total energy of configuration 2 is the lowest and the structure is the most stable, which was determined by calculating the total energy of the system. Therefore, it can be seen that the proline and 5-FU in the 5F-P are connected by N1–H1•••O4 hydrogen bonding based on THz spectroscopy and system stability. These results demonstrate that THz-TDS can well characterize the formation of new phases. The combination of THz-TDS and DFT can effectively study the connection mode of hydrogen bonds.

This work shows that the combination of THz-TDS and DFT can provide an effective
method to investigate the connection mode of hydrogen bonds and thus aid in characterizing compounds such as this potential drug candidate.

Xun Zhang

Articles

Terahertz Spectral Investigation of L-Cysteine Hydrochloride and its Monohydrate

Terahertz Spectroscopic Analysis of Co-Crystallized Mixtures in an L-threonine Diastereomer System

The Hydrogen Bonding Mode of 5-Fluorouracil-Proline Co-Crystal (5F-P) Studied by Terahertz (THz) Spectroscopy and Density Functional Theory (DFT)