Molecular Structure, Electronic Properties, Reactivity (ELF, LOL, and Fukui), and NCI-RDG Studies of the Binary Mixture of Water and Essential Oil of Phlomis bruguieri

Essential oils are volatile oil-like liquids with a characteristic strong smell and taste. They are formed in plants and are then extracted. Essential oils have extremely strong physiological and pharmacological properties, which are used in the medicine, cosmetics, and food industries. In this study, the molecules caryophyllene oxide, β-pinene, 1,8-cineol, α-cubebene, and β-caryophyllene, which are the molecules with the highest contents in the essential oil of the plant mentioned in the title, were selected and theoretical calculations describing their interactions with water were performed. Because oil–water mixtures are very important in biology and industry and are ubiquitous in nature, quantum chemical calculations for binary mixtures of water with caryophyllene oxide, β-pinene, 1,8-cineol, α-cubebene, and β-caryophyllene were performed using the density functional theory (DFT)/B3LYP method with a basis of 6–31 G (d, p). Molecular structures, HOMO–LUMO energies, electronic properties, reactivity (ELF, LOL, and Fukui), and NCI-RDG and molecular electrostatic potential (MEP) on surfaces of the main components of Phlomis bruguieri Desf. essential oil were calculated and described.


Introduction
Phlomis is a genus belonging to the family Lamiaceae. It is known by the local names "çalba" and "shalba". The genus Phlomis includes more than 100 species of plants and shrubs distributed in Europe, Asia, and North Africa, and 52 taxa, including 6 species, 12 natural hybrids, and 34 endemic taxa, are known to grow in Turkey [1,2]. In traditional Turkish folk medicine, the flowers and/or leaves of Phlomis species are widely used as herbal teas, as a tonic, carminative, and appetite stimulant, and for stomach ache. Besides these traditional uses, Phlomis plants have been shown to have antidiabetic, ulcerogenic, antimicrobial, anti-inflammatory, antinociceptive, antimutagenic, immunosuppressive, and free radical scavenging properties. Furthermore, phytochemical investigations of members of the Phlomis genus have shown that these plants contain iridoids, flavonoids,

Experimental Study
The plant material used in the study was collected from the natural environment during the flowering period of the plant and dried in the shade. The oil obtained from the plants was extracted using the hydrodistillation method. The oil obtained by hydrodistillation of the air-dried plant materials was injected into a GC/GC-MS device after the distillation process was completed. The retention indices of the obtained components were used as reference points to calculate the retention indices and the chemical compounds, retention times, and mass spectra of the essential oil were determined.

Theoretical Study 2.2.1. HOMO-LUMO Analysis
Noncovalent interactions in a binary mixture play an important role in determining the mixture's thermodynamic properties and reactivity [18,19]. Mixtures of simple substances with water are subjects of active study both experimentally and theoretically [20,21].
Studies of the components of essential oils using experimental and theoretical methods play an important role in understanding their impacts and potential uses [13,30,31]. First, after optimizing the structures, we performed a HOMO-LUMO analysis.
The highest occupied molecular orbital and the lowest unoccupied molecular orbital are very important parameters in quantum chemistry [32]. This information is also used along with the boundary electron density to predict the most reactive positions in various systems [33]. The highest occupied molecular orbital and the lowest unoccupied molecular orbital are the main orbitals involved in chemical stability [34]. The energy difference between the HOMO and LUMO orbitals is called the energy gap, which is an important stability characteristic of structures. In Figure 1, green and red colors represent the positive and negative phases, respectively. The energy values of the HOMO (donor) and LUMO (acceptor), and their energy gap, reflect the chemical activity of a molecule [35,36]. In general, an atom with a higher HOMO density should have a greater ability to donate an electron, while an atom with a higher occupied LUMO density should have a greater tendency to gain an electron [37]. Three-dimensional plots of the HOMO and LUMO orbitals for the binary mixtures of the essential oil components with water are indicated in Figure 1.

HOMO-LUMO Analysis
Noncovalent interactions in a binary mixture play an important role in determining the mixture's thermodynamic properties and reactivity [18,19]. Mixtures of simple substances with water are subjects of active study both experimentally and theoretically [20,21].
Studies of the components of essential oils using experimental and theoretical methods play an important role in understanding their impacts and potential uses [13,30,31]. First, after optimizing the structures, we performed a HOMO-LUMO analysis.
The highest occupied molecular orbital and the lowest unoccupied molecular orbital are very important parameters in quantum chemistry [32]. This information is also used along with the boundary electron density to predict the most reactive positions in various systems [33]. The highest occupied molecular orbital and the lowest unoccupied molecular orbital are the main orbitals involved in chemical stability [34]. The energy difference between the HOMO and LUMO orbitals is called the energy gap, which is an important stability characteristic of structures. In Figure 1, green and red colors represent the positive and negative phases, respectively. The energy values of the HOMO (donor) and LUMO (acceptor), and their energy gap, reflect the chemical activity of a molecule [35,36]. In general, an atom with a higher HOMO density should have a greater ability to donate an electron, while an atom with a higher occupied LUMO density should have a greater tendency to gain an electron [37]. Three-dimensional plots of the HOMO and LUMO orbitals for the binary mixtures of the essential oil components with water are indicated in Figure 1.
According to the data given in Table 1, the introduction of a water molecule led to changes in the values of the energy gap and other physical parameters. Thus, for αcubebene, a decrease in the size of the energy gap was observed when a water molecule was introduced into the system. A similar change was also observed for β-pinene. For β-caryophyllene, caryophyllene oxide, and 1,8-cineole, the observations were different. For these substances, the introduction of a water molecule increased the size of the energy gap.
The data obtained partially agree with the literature data. Thus, in [39], for pure α-cubebene, the observed values of E HOMO , E LUMO , and Egap corresponded to −5.69, 0.90, and 6.59, respectively; for a mixture with water, these values were −5.96, 0.51, and 6.47 eV, respectively. When comparing the data given in [39] and Table 1, it is apparent that the introduction of one water molecule into a system led to changes in most physical parameters. For example, the value of ω for pure α-cubebene was 0.44 [39], and for its binary mixture in water was 1.1481.
The introduction of a water molecule also affected other physical parameters (Table 1). In general, there was a correlation with the previously described phenomenon, in that those parameters that increased for α-cubebene and β-pinene when a water molecule was introduced into the system decreased for β-caryophyllene, caryophyllene oxide, and 1,8-cineole, and vice versa.
Thus, it can be assumed that the introduction of water into these systems slightly changed the values of the energy gap for the substances under study. The data obtained were consistent with the results given in [20,21].

MEP Analysis
The net electrostatic effect produced by the overall charge distribution (electrons + nuclei) of a molecule at a point in space around it is known as the molecular electrostatic potential. It relates the total charge distribution to dipole moments, partial charges, electronegativity, and the site of chemical activity of the molecule [40]. In Figure 2, different electrostatic potential values on the surface are represented by different colors: red represents the areas with the most negative electrostatic potential, blue represents the areas with the most positive electrostatic potential, and green represents areas with zero potential.
The most positive and negative regions, represented by blue and red colors on the MEP map, are the preferred locations for nucleophilic and electrophilic attacks, respectively [41]. Such maps are very useful when studying molecular structure and physicochemical characteristics because they use color grading to represent areas of positive, negative, and neutral electrostatic potential [42].
In all the molecules we studied, the electrostatic potential was primarily limited to oxygen, which is a reactive region for electrophilic attack ( Figure 2). It should be noted that the red region was located above the water oxygen in the binary mixtures under study, and the blue region was located mainly above the hydrogen atom in the water ( Figure 2).

Energy of Population Density of States
Since nearby orbitals in a boundary region may have quasi-degenerate energy levels, using only the HOMO and LUMO to describe the boundary orbitals may not be complete enough [43]. The DOS spectrum explains the contribution of electrons to the conduction and valence bands. The spectra of the investigated binary mixtures, which illustrate how many states are available at different energy levels, are presented in Figure 3.
A virtual orbital is also known as an acceptor orbital because it is not occupied. On the other hand, a filled orbital is known as a donor orbital. DOS is positive for binding interactions, negative for antibinding interactions, and zero for no binding interaction [44].

Energy of Population Density of States
Since nearby orbitals in a boundary region may have quasi-degenerate energy levels, using only the HOMO and LUMO to describe the boundary orbitals may not be complete enough [43]. The DOS spectrum explains the contribution of electrons to the conduction and valence bands. The spectra of the investigated binary mixtures, which illustrate how many states are available at different energy levels, are presented in Figure 3. A virtual orbital is also known as an acceptor orbital because it is not occupied. On the other hand, a filled orbital is known as a donor orbital. DOS is positive for binding interactions, negative for antibinding interactions, and zero for no binding interaction [44]. According to the data shown in Figure 3, the TDOS (black) and PDOS frag. 1 (red) curves had the same shape and size for all the studied substances from −21.77 to −15.57 eV. In the range of −15.56 to 5.44 eV, these curves had similar shapes, but different sizes, while the TDOS curve, as a rule, was above the PDOS frag. 1 curve. It should be noted that PDOS frag. 2 (blue) had very different shapes and sizes from PDOS frag. 1. Thus, the introduction of water molecules strongly affected the TDOS, PDOS, and OPDOS data of the substances under study.

Molecules
Center of TDOS (eV) According to the data shown in Table 2, the introduction of a water molecule led to a shift in the center of the PDOS towards a lower value for all the substances under study. Thus, for α-cubebene, a shift from −6.718485 to −5.640946 eV was observed when a water molecule was added.

NCI-RDG Analysis
RDG is an efficient method used to study noncovalent interactions in real space based on electron density and its derivatives [45]. This approach has been shown to be able to distinguish between hydrogen bonds, van der Waals interactions, and steric repulsion in small molecules and molecular complexes. To calculate the RDG values, we used the following equation: By analyzing the low-density gradient, regions of low electron content, which are responsible for weak interactions, can be confirmed; similarly, the high-density gradient value is used to identify strong interactions. Plots of RDG versus ρ (low-gradient peakdensity values) seem to be an indicator of the strength of an interaction. The product of the electron density ρ and the sign of the second eigenvalues λ 2 has been proposed as a tool to distinguish between a wide range of interaction types. The sign λ 2 is a critical parameter used to distinguish between bound (λ 2 < 0) and noncoupling (λ 2 > 0) interactions. On the basis of (λ 2 )ρ, interactions can be divided into three types. The first is an attraction (hydrogen bond or dipole-dipole) and is characterized by negative sign values (λ 2 )ρ. The second assumes a nonbinding interaction of a strong repulsion type, which is confirmed by large and positive values of the sign (λ 2 )ρ. The third type, namely weak interactions (van der Waals interactions), is associated with values close to zero. Figure 4 shows the predicted RDG plots for the test substances. The interactions between water and α-cubebene, β-caryophyllene, caryophyllene oxide, 1,8-cineole, and βpinene were described by different types of interactions. The first-red color-was typical for all the studied substances, indicating a strong repulsion within the system. The secondgreen-showed the van der Waals interactions. The third-blue color-corresponded to strong hydrogen bond interactions and was clearly observed between the water molecules and caryophyllene oxide and 1,8-cineole.
an attraction (hydrogen bond or dipole-dipole) and is characterized by negative sign values (λ2). The second assumes a nonbinding interaction of a strong repulsion type, which is confirmed by large and positive values of the sign (λ2). The third type, namely weak interactions (van der Waals interactions), is associated with values close to zero. Figure 4 shows the predicted RDG plots for the test substances. The interactions between water and α-cubebene, β-caryophyllene, caryophyllene oxide, 1,8-cineole, and βpinene were described by different types of interactions. The first-red color-was typical for all the studied substances, indicating a strong repulsion within the system. The second-green-showed the van der Waals interactions. The third-blue colorcorresponded to strong hydrogen bond interactions and was clearly observed between the water molecules and caryophyllene oxide and 1,8-cineole.

ELF and LOL Analyses
Surface analysis based on covalent bonding was performed using an electron localization function (ELF) with a bump map and LOL maps. A bump map with a large or narrowed peak area describes the electronic environment of an atom. They offer a

ELF and LOL Analyses
Surface analysis based on covalent bonding was performed using an electron localization function (ELF) with a bump map and LOL maps. A bump map with a large or narrowed peak area describes the electronic environment of an atom. They offer a valuable opportunity to identify pairs of electrons on the molecular surface [46,47]. ELF and LOL have similar chemical mapping because they depend on the kinetic energy density. However, ELF arises from the electron pair density; LOL simply exhibits a gradient of localized orbitals and is used when localized orbitals overlap [46]. Images of the studied substances in ELF and LOL formats are shown on the color tint map and contour maps of the hydrogen bond region in Figure 5. The ELF map was built in the range of 0.0-1.0; however, delocalized electrons are projected in a reduced range (<0.5) [48,49]. If an electron position dominates in the electron density, then the LOL phase reaches the upper region (>0.5) [50,51]. Due to the presence of a covalent bond or nuclear layer, a high position of electrons is created in this region, as demonstrated by its high value in this region. The electron cloud was delocalized around some carbon atoms in the studied substances, which is represented by blue areas. Critical points and their trajectories, chemical bonding, and chemically significant regions (red and orange) are clearly shown on the ELF map, localized mainly around hydrogen atoms and present in all substances under study. It should be noted that the central regions of some hydrogen atoms in all the studied substances were white, which indicates that the electron density exceeded the upper limit of the color scale (0.80) ( Figure 5).

Fukui Analysis
The condensed Fukui functions are determined to characterize the regioselectivity of an atom 'r' in a molecule. For this purpose, the Fukui functions (f + (r), f − (r), f 0 (r), and ∆f(r)) were calculated using Multiwfn [52]. where q (N) (r) is the atomic charge on site 'r' in a neutral system; q (N−1) (r) is the atomic charge on site 'r' in a cationic system; and q (N+1) (r) is the atomic charge on site 'r' in an anionic system. Additionally, the dual identifier ∆f(r) is the difference between nucleophilic and electrophilic attacks at a specific site and can be calculated according to the following equation [54][55][56]: If this value is negative, the particular area can be considered an electrophilic attack site. If this value is positive, then it may be a nucleophilic attack site.

Collection of the Plant Material
The plant material used in the study was collected in Turkey. It was collected from its natural habitat in June, which is the flowering period of the plant, and dried in the shade. The plant specimens are preserved at Bingol University.

Isolation of Essential Oil and GC-MS Analysis
The plants used in this study were air-dried. The oil was extracted from the plants using the hydrodistillation method. The air-dried aerial parts (200 g) of the plants were hydrodistilled for three hours in a Clevenger-type apparatus. The organic layer in the collection vial was injected into the GC/GC-MS instrument after the distillation was completed.
GC-MS was used to examine the essential oil. The MS instrument was a Shimadzu GCMS-QP2020 model. The mass range was 40-330 m/z and the ionization energy was 70 eV. An RXI-5MS capillary column (30 m × 0.25 mm × 0.25 m) was used. The carrier gas was helium, with a constant column flow rate of 1 mL/min. The column oven temperature program was set to a temperature of 40 • C with a hold time of two minutes and a rate of change of 3 • C/min. The final temperature was 240 • C The injection volume was set to 1 l and the mode was set to split (split ratio 1:10 or 1:100). A 3.5 min solvent delay was used for hexane samples. The mass spectrometric parameters were as follows: full scan mode, 20,000 amu/s scan speed, and 50 spectrum sampling frequency. The temperatures at the interface and the ion source were 250 • C and 200 • C, respectively.
Alkanes were used as reference points to calculate the retention indices (RI). The chemical compounds of the essential oils were identified by comparing their retention times (RT) and mass spectra to those reported in MS libraries (Wiley MS library, New York, NY, USA) and the literature (NIST 20 and Wiley Libraries) [57]. A traditional library search ignores retention parameters, and the process consists solely of comparing spectra. When searching libraries in this study, a combination of storage indexes were used, making compound identification easier and more reliable. The device's retention index spectral libraries were also used in this study. The same analysis method as described for the column library was used for better results.
The identified constituents of the essential oil are listed in Table 5.

Theoretical Study
DFT calculations were performed using Gaussian 09W software [76] and the GausView 0.5 program [77]. The optimized molecular structures of the binary mixtures of α-cubebene, β-caryophyllene, caryophyllene oxide, 1,8-cineole, and β-pinene with water were calculated with the 6-31 G(d,p) basis set using Becke's three-parameter hybrid functional (B3) [78] and the Lee-Yang-Parr (LYP) correlation function [79]. Molecular electrostatic potential surface maps and HOMO-LUMO plots were plotted using the same method. Using the Multiwfn software [80], total density of state (TDOS), partial density of state (PDOS), and overlap population density of state (OPDOS) analyses were carried out for the binary mixtures. Furthermore, scatter plots of the reduced density gradients (RDGs) and noncovalent interactions (NCI) were constructed using Multiwfn software and visualized using VMD software [81]. Finally, ELF, LOL, and Fukui functions were calculated using Multiwfn 3.8 software.

Conclusions
Essential oils have a wide range of beneficial properties. In this work, we isolated the essential oil from Phlomis bruguieri. It was found that caryophyllene oxide, β-pinene, 1,8-cineol, α-cubebene, and β-caryophyllene showed the maximum yields among the substances found in the essential oil of this plant. For these substances, calculations describing their mixtures with water were carried out, since this information also has industrial significance. It was shown that the introduction of water into these systems affected the electronic properties of the molecules, including changing the size of the energy gaps. Molecular structure, electronic properties, reactivity (ELF, LOL, and Fukui), and NCI-RDG data were obtained for the tested substances. The energy gaps observed in the DOS and HOMO-LUMO spectra for the studied substances were comparable and consistent with each other. The results obtained in this work can be applied in various fields from aromatherapy [82] and medicine [83] to encapsulation and food chemistry [84,85].
In addition, the theoretical data obtained in this work may help others to understand the noncovalent interactions of some components of essential oils with water.