In silico Evaluation of Cucurbit[6]uril as a Potential Detector for Cocaine and Its Adulterants Lidocaine, Caffeine, and Procaine

Abstract

Illicit drugs and their trafficking require worldwide efforts in investigation, detection, and control. Colorimetric tests are often applied to identify drugs. Cocaine has some well-known adulterants that can provide a false positive response. Cucurbit[6]uril (CB[6]) has been suggested as a potential detector for cocaine and other illicit drugs. This work uses in silico methods to evaluate the use of CB[6] to detect cocaine and these interfering substances. More specifically, this work analyzes different possibilities of CB[6] complexation with cocaine, lidocaine, caffeine, and procaine and compares the results achieved for cocaine and its adulterants. Different methodologies were employed: quantum chemistry was investigated through DFT B3LYP/TZVP (density functional theory-Becke, three-parameter, Lee-Yang-Parr with triple zeta valence plus polarization basis set) and the semi-empirical methods Austin model 1 (AM1), parametric methods 3, 6, and 7 (PM3, PM6, PM7), and Recife model 1 (RM1). We used these methodologies intending to compare the reasonability and reproducibility of the results in the gas phase condition. Solvent influence was studied by molecular dynamics (MD) simulations. Results showed that CB[6] does not bind to these substances, as judged from the positive values of binding free energy obtained with all methods. DFT and MD were the most reliable methods whereas semiempirical ones were not reproductible in describing these systems. Results also showed that interactions are not specific, so CB[6] does not provide a good response for cocaine detection.

Keywords:
cucurbit[6]uril; cocaine; lidocaine; caffeine; procaine; in silico methods

Introduction

Drug trafficking and consumption are worldwide concerns, and integrated efforts are needed to combat them. Among illicit drugs, cocaine has a potential supply-driven market in expansion, while coca bush cultivation has increased concomitantly with the global drug manufacture.¹

Correct drug identification is essential to ensure human rights: it can prevent the arrest of innocents and avoids the miscarriage of justice.^2,3 It is also important in harm reduction projects because there are advantages in testing drugs before they are consumed. Many harm reduction agencies are exploring techniques for illicit drug testing to identify and, when possible, to quantify their constituents, thereby allowing their users to make informed decisions.^3,4

Drug characterization comprises two steps: presumptive (in situ) and confirmatory (in lab) tests. Presumptive tests are less time-consuming and less expensive. However, most of these tests consist of colorimetric tests and can give doubtable results.⁵ The Scott test is the most widely employed for cocaine identification in law enforcement actions. Nevertheless, it can return false positives because cocaine is not sold in its pure form. On the other hand, confirmatory tests are more trustful, despite being more expensive as they require specialized equipment, trained examiners, and complex sample preparation, besides they can be destructive.^2,6 Therefore, new materials and methods are constantly being searched in an attempt to obtain trustful identification at a low cost.^7-11

Cucurbiturils (CB) can form host-guest complexes with various drugs consisting of small organic or inorganic molecules.¹² Hydrophobic effects facilitate drug encapsulation within the cucurbituril cavity, and further hydrogen bonding or ion-dipole interactions with the cucurbituril portal stabilize the resulting complex.¹³ The defining structural features of the CB[n] family are their highly symmetric structure bearing negatively charged carbonyl rims and a hydrophobic cavity.^14,15 Cucurbituril (cucurbit[6]uril, or CB[6]) is a hexameric macrocyclic compound originating from self-assembly during an acid-catalyzed condensation reaction between glycoluril and formaldehyde. Curcubituril is commonly abbreviated as CB[6] or even CB6, Q[6], Q6, or Cuc6 in some cases; ‘6’ represents the number of glycoluril units in the macrocycle.^16,17 The CB[n] cavity can recognize amines and azaheterocyclic compounds via hydrophobic effects, size effects, and ion-dipole interactions at the CB[n] portals.^18,19 The rigid structure and the ability to form stable complexes with molecules and ions make CB[6] an attractive building block to construct supramolecular architectures.²⁰ CB[6] can also facilitate cycloaddition. CBs have been used in drug identification studies for opioids, marijuana and cocaine, 3,4-methylenedioxymethamphetamine (MDMA) and amphetamine-type stimulants.^14,21-23

Previous in silico methods based on the semi-empirical methods PM7 (parametric method 7)²² and RM1 (Recife model 1)²³ have been applied to study cucurbituril complexes. GAFF (generalized Amber force field) have also been used to examine how acyclic cucurbituril (aCBs) and opioids interact.²¹ The DFT-B3LYP (density functional theory-Becke method, three-parameter, Lee-Yang-Parr) hybrid functional and the 6-31G(d) basis set have been applied to study reactive blue 19 dye adsorption onto CB[6] and CB[8].²⁴ Microsecond time scale molecular dynamics simulations have aided the investigation into the binding enthalpies of CB[7] with different guests in aqueous solution.²⁵ Adsorption of hydrogen sulfide by CB[7] has been studied by ab initio van der Waals density functional (vdW-DF) calculations.²⁶

This work aimed to employ in silico methods to evaluate the CB[6] ability to form a complex with cocaine and some of its adulterants: lidocaine, caffeine, and procaine. Quantum chemistry and molecular dynamics approaches were used. Our intention was to verify the differential energy involved with CB[6] and the molecules of interest as an indicator of the interaction between these compounds.

Methodology

In silico methods include several methodological approaches. A combination of theoretical and experimental work is usually used to find agreement among them. Despite the importance of experimental evaluation, computational methods can provide valuable and detailed information that is not possible to obtain directly from empirical observation.²⁷ However, the theoretical tool must provide an accurate result to reproduce experimental findings. That is why we decided to carry out these analyses by combining different methods: semiempirical, DFT and MD (molecular dynamics) simulations.

Quantum chemistry study

Crystallographic structures of the molecules studied herein were obtained from The Cambridge Structural Database (CCDC).^28-31 Figure 1 presents the studied system. The CCDC structures were previously optimized with UFF (universal force field) in the Avogadro v.1.1.1^32,33 software. This force field can optimize the geometry of all elements and can properly handle inorganic materials and organometallic materials.³⁴

Figure 1
Study system.

After optimization with UFF, the molecules were submitted to semi-empirical optimization with the MOPAC software by applying five different methods: Austin model 1 (AM1), parametric methods 3, 6, and 7 (PM3, PM6, PM7), and Recife model 1 (RM1).^35-40 All the molecules were re-optimized with the hybrid DFT-B3LYP⁴¹ with a single zeta Ahlrichs basis set, and then single point refinement was carried out with a triple zeta Ahlrichs basis set (TZVP).⁴² All calculations were performed by ORCA v4.0.1.2 software.⁴³

To characterize the stationary state, we have used the crystallographic structures for all molecules as a starting point for further optimization. No imaginary frequencies were found, indicating the achievement of a minimum energy structure for each molecule.^44-48 The thermodynamic functions and the final energies of the complexes were obtained.²⁴

To evaluate the tendency of the molecules to form a complex with CB[6], the hypotheses listed in Table 1 were investigated from a theoretical viewpoint. The table also defines the respective abbreviations. In all cases, the energy tendency was evaluated by the binding energy (ΔE_bind) equation 1:^49,50

Table 1
Drug positions in relation to CB[6]

Position	Abbreviation
Cocaine and CB[6]
Cocaine in the center	COC_P01
Amino group oriented toward the center	COC_P02
Aromatic ring oriented toward the center	COC_P03
Ester carbonyl group ester oriented toward the center	COC_P04
Lidocaine and CB[6]
Lidocaine in the center	LID_P01
Amino group oriented toward the center	LID_P02
Aromatic ring oriented toward the center	LID_P03
Caffeine and CB[6]
Caffeine in the center	CAF_P01
Five-member ring oriented toward the center	CAF_P02
Aromatic ring oriented toward the center	CAF_P03
Procaine and CB[6]
Procaine in the center	PROC_P01
Amino group oriented toward the center	PROC_P02
Aromatic ring oriented toward the center	PROC_P03

• CB[6]: cucurbit[6]uril.

(1) $\Delta E_{\mathrm{bind}}=E\left(\mathrm{complex}\right)-\left(E\left(\mathrm{CB}\left[6\right]\right)+E\left(\mathrm{Mol}\right)\right)$

where Mol = cocaine, lidocaine, caffeine, or procaine.

MD simulations

Molecule parametrization

CB[6] and the studied drugs were parametrized according to a protocol described in more detail elsewhere.⁵¹ Briefly, the molecules were subjected to structural optimization followed by single-point energy calculations in Gaussian 09; the HF/6-31(d) method and the Merz-Kollman scheme were employed.⁵² Point charges were derived through the restricted electrostatic potential (RESP) fit available in the antechamber program from Amber 16 suite.⁵³ The molecules were then fully parametrized with the aid of the GAFF.^53,54

MD simulations and calculation of binding enthalpies

For each CB[6]-drug system, different manually-generated conformations were subjected to MD simulations to determine their time stabilities and their respective binding enthalpies (ΔH_bind). Each studied system was embedded into a cubic box and solvated in methanol to match the experimental conditions.^23,55-60 The box size was adjusted so that the minimum distance between the box walls and the solute surface was 10.5 Å and the number of methanol molecules was exactly 630 in all cases. The methanol explicit solvation model available in Amber 16 was used to conduct the MD simulations.^53,61 No counterions were added because the solute molecules were all neutral. Then, the systems were subjected to a sequence of energy minimization (EMs), equilibration, and production runs according to a modified version of a protocol reported by Fenley et al.,²⁵ which allowed the ΔH_bind values to be calculated through a multibox strategy. In addition, some other adaptations were made during the equilibration stages. Differently from the tutorial, isothermal-isochoric (NVT) equilibration was accomplished for 500 ps, while the solute heavy atoms were kept restrained with a spring constant (k) of 10 kcal Å^-2 mol^-1. Heating was carried out by linearly increasing the temperature from 10 to 298.15 K; this final value was kept constant in all subsequent simulations. A subsequent 500-ps isothermal-isobaric (NPT) equilibration was performed by applying the same previous restraints on the solute heavy atoms, which were then decreased down to k = 2 kcal Å^-2mol^-1 in a stepwise fashion during three consecutive 500-ps NPT equilibrations (k = 8, 6, and 2 kcal Å^-2 mol^-1, respectively). Next, in accordance with the tutorial, 40-ns unrestrained NPT equilibrations were conducted for all the systems to determine the respective average box sizes. Finally, 1 µs NVT production runs were accomplished for the stable systems after their box sizes were modified to match the respective average volumes calculated from the precedent 40-ns MD simulations. Three replicate MD simulations were run for each CB[6]-drug conformation by following the same sequence of steps previously described, and different random velocities were assigned during the heating process. All the EMs and heatings were carried out with Amber 16 pemd.MPI, while NPT equilibrations and production runs were conducted with Amber 16 pmemd.cuda.⁵³ Treatment of nonbonded interactions and temperature and pressure control, among others, were identical to those employed by Fenleyet al.²⁵ For the free drugs and CB[6], the average blocking approach was applied to estimate the standard error of the mean (SEMs) of their average potential energies (E_ptot) required for ΔH_bind calculations, as only one productive run was carried out for these systems. The gmx analyze program available in GROMACS v.5.1.3 was employed for the previous analyses.⁶² The E_ptot mean values for the CB[6]-drug poses were obtained by averaging the results for the replicate MD simulations and using only those frames in which the drug was encapsulated. SEMs were estimated using the following formula: SEM = SD/√3, where SD is the standard deviation of the mean values in the three replicate MD simulations. Finally, the errors were propagated in order to calculate the SEMs of the reported ΔH_bind values.

Free energy calculations

The attach-pull-release (APR) protocol reported by Henriksen et al.⁶³ was used to calculate the binding free energies (ΔG_bind) of the studied CB[6]-drug complexes that remained stable during, at least, one of their respective 1 µs MD simulation replicates. Last frames generated during such stable trajectories were employed as the starting structures to perform the APR simulations, which were automatically carried out by means of an updated set of python scripts that enable the use of various non-aqueous solvents, including methanol, during the system setup. Except for the solvent, the MD simulations were conducted by the default inputs of the APR protocol in most cases. The solvent input was set to methanol, and 900 molecules of the latter solvent were added to each solvation box. The production runs of each window were simulated for up to 40 ns, the first 10 ns being discarded from subsequent free energy calculations. Hydrogen mass repartitioning (HMR) was employed to increase the simulation time step to 4 fs. Radial conformational restraints, i.e., jacks, were also set at the entrance of the CB[6] cavity to allow smooth exit of the encapsulated drug molecule.⁶³ Three pairs of oxygen atoms (O,O) and six pairs of nitrogen atoms (N,N), all lying in radially opposite positions at the entrance ring where the drug was placed, were chosen to define the jacks. In all cases, the jack distance and the force constant were set to 10 Å and the default value (25 kcal Å^-2 mol^-1), respectively. The complexes with a six-member ring inserted in CB[6] were the exception: the jack distance was set to 11.5 Å due to the wider molecular shape of the enclosed moiety. The free energies and their uncertainties were determined by the thermodynamic integration (TI) approach and blocking analysis, respectively. The equilibration and production runs were performed with Amber 16 pmemd.MPI and pmemd.cuda, respectively.⁵³

Trajectory analysis and visualization

Root mean square deviations (RMSDs) for the drug heavy atoms were calculated throughout the concatenated heating, NPT equilibration and production runs for every trajectory, by fitting all frames to the CB[6] heavy atoms and taking the minimized starting structure as reference. As a complementary way to assess complex stability, the solvent accessible surface area (SASA) of CB[6] was monitored for every frame during all MD simulations of CB[6]-drug poses. The bondi set of van der Waals atomic radii and the MOLSURF program available in Amber 16 were employed in such calculations, which, together with those of RMSD values, are integrated into the cpptraj module.⁵³ All figures depicting the molecules were created with PyMOL v2.1.0.⁶⁴

Results

Quantum chemistry methods

Table 2 summarizes all the results for the calculations performed with DFT methods for those conformations that achieved convergence. No relevant information was obtained from semi-empirical methods because they did not reach convergence for most conformations. Semiempirical results are shown in the Supplementary Information (SI) section. We observe that the DFT calculation returned positive values for the Gibbs free energy, evidencing that encapsulation was not spontaneous even when exothermic values were obtained for enthalpic contribution.

Table 2
Calculated values with triple zeta Ahlrichs basis set (TZVP)

	Electronic energy ΔEEbind / (kcal mol-1)	Enthalpy ΔHbind / (kcal mol-1)	Gibbs energies ΔGbind / (kcal mol-1)	Entropic contribution TΔSbind / (kcal mol-1)
COC_P01	15.1165	15.5390	31.1182	17.852
COC_P03	14.9764	15.9978	29.7191	15.994
LID_P01	29.2124	31.0687	49.7492	18.680
CAF_P01	14.2829	15.5320	32.5367	17.004
CAF_P02	13.5701	14.7878	31.7852	16.997
PROC_P01	31.7086	33.5801	53.7916	20.211
PROC_P03	1.6990	3.4646	20.1780	16.713

• T: temperature.

Figure S1 (available in SI section) illustrates the final conformation obtained by B3LYP/TZVP. The final conformations obtained by semi-empirical methods are provided in SI section too.

Molecular dynamics

We defined each starting conformations of the analyzed systems according to Table 1. In this section we will assess the stability of the encapsulated drugs by using MD simulations and thermodynamic calculations based on the generated ensembles.

Figure 2 shows the outcomes of the MD simulations conducted for the different poses of cocaine in complex with CB[6]. Note that poses COC_P01 and COC_P03 converged to equivalent conformations during their respective MD simulation replicates that remained stable up to 1 µs. DFT calculations evidenced the same behavior (Table 2 and Figure S1, SI section). In fact, for the CB[6]-cocaine complexes, the only binding mode that remained stable in at least one long MD simulation replicate was the mode bearing the phenyl group inside the host cavity. In the three remaining poses, cocaine either dissociated from CB[6] or converged to the same previous stable structures.

Figure 2
Time evolution of various starting CB[6]-cocaine structures during MD simulations. (a) and (b) show for each pose the evolution of SASA of the CB[6] molecule and the RMSD values with respect to the drug heavy atoms and taking as reference its position in the minimized starting structure, which is depicted in (c). The insets in certain graphs revealed steep variation in RMSD values after the restraints applied during equilibration were removed. (d) Structural representation of the last frames collected from each replicate MD simulation of the analyzed poses.

MD simulations carried out for different initial poses of the CB[6]-lidocaine complex identified a stable binding mode up to the completion of the simulation time in some replicates, while others showed the same binding mode dissociating from the host. Figure 3 depicts structural representations of the initial (t = 0) and final frames of the MD simulations in each case, together with the SASA of the CB[6] molecule and RMSD time profiles calculated with respect to the initial position of the drug heavy atoms. The insets in certain graphs demonstrated that RMSD values varied widely after restraints applied during equilibration were removed. Out of the three poses analyzed in this work, two of them (LID_P01 and LID_P02) converged to the same final structures in certain replicates, in which the diethylamine moiety of lidocaine was inserted in CB[6]. DFT calculations revealed similar final convergence to LID_P01. Conversely, the lidocaine pose in which the aryl group was placed into CB[6] readily dissociated during a shorter 42-ns unrestrained MD simulation (LID_P03).

Figure 3
Time evolution of various CB[6]-lidocaine starting structures during MD simulations. (a) and (b) show for each pose the evolution of SASA of the CB[6] molecule and the RMSD values with respect to the drug heavy atoms and taking as reference its position in the minimized starting structure, which is depicted in (c). The insets in certain graphs revealed steep variation in RMSD values after the restraints applied during equilibration were removed. (d) Structural representation of the last frames collected from each replicate MD simulation of the analyzed poses.

Figure 4 displays the results for caffeine. In this case, we identified a stable conformation of the CB[6]-caffeine complex in at least one MD simulation replicate, in which the drug’s six-membered ring was inserted in CB[6], whereas the other two poses were unstable during their respective MD simulations (Figure 4).

Figure 4
Time evolution of various CB[6]-caffeine starting structures during MD simulations. (a) and (b) show for each pose the evolution of SASA of the CB[6] molecule and the RMSD values with respect to the drug heavy atoms and taking as reference its position in the minimized starting structure, which is depicted in (c). The insets in certain graphs revealed steep variation in RMSD values after the restraints applied during equilibration were removed. (d) Structural representation of the last frames collected from each replicate MD simulation of the analyzed poses.

On the other hand, procaine presented two stable conformations according to Figure 5, each of them bearing a different group lying within the CB[6] cavity. Note that poses PROC_P01 and PROC_P02 converged to equivalent conformations during all their respective MD simulation replicates. The diethylamine moiety inserted in CB[6] was obtained from two different starting poses during independent MD simulations.

Figure 5
Time evolution of various CB[6]-procaine starting structures during MD simulations. (a) and (b) show for each pose the evolution of SASA of the CB[6] molecule and the RMSD values with respect to the drug heavy atoms and taking as reference its position in the minimized starting structure, which is depicted in (c). The insets in certain graphs revealed steep variation in RMSD values after the restraints applied during equilibration were removed. (d) Structural representation of the last frames collected from each replicate MD simulation of the analyzed poses.

After analyzing stability, we selected the stable conformations of each host-guest system to conduct thermodynamic calculations. We subjected five complexes to the APR protocol to calculate their respective ΔG_bind values (see Figures S3 and S4, SI section). In addition, we accomplished complete thermodynamic characterization of such complexes by estimating their ΔH_bind values through the independent multibox approach described in Methodology section, which in turn allowed us to calculate the entropic contributions (TΔS_bind). Table 3 summarizes the results of the previous calculations for all the analyzed systems.

Table 3
Thermodynamic parameters associated with the formation of CB[6]-drug complexes

Energy componentsa	Complexb
	LID_P01 and LID_P02	CAF_P03	PROC_P01 and PROC_P02	PROC_P03	COC_P01 and COC_P02
Wattach / (kcal mol-1)	118.0 ± 0.5c	474.7 ± 0.2	109.8 ± 0.2	460.8 ± 0.1	459.9 ± 0.4
Wpull / (kcal mol-1)	3.7 ± 0.4	-23.9 ± 0.3	5.5 ± 0.5	4.1 ± 0.2	0.8 ± 0.2
Wrelease / (kcal mol-1)	-119.1 ± 0.1	-466.0 ± 0.1	-119.1 ± 0.1	-466.0 ± 0.1	-465.8 ± 0.2
Wrelease-std / (kcal mol-1)	-7.1	-7.1	-7.1	-7.1	-7.1
ΔGbind / (kcal mol-1)	4.5 ± 0.6	22.4 ± 0.4	11.0 ± 0.5	8.1 ± 0.2	12.3 ± 0.4
ΔHbind / (kcal mol-1)	3.5 ± 0.2	47.3 ± 0.2	-0.3 ± 0.2	3.6 ± 0.2	5.8 ± 0.2
-TΔSbind / (kcal mol-1)	1.0 ± 0.6	-24.9 ± 0.4	11.3 ± 0.5	4.5 ± 0.3	6.5 ± 0.4

• ^a
  The energy components associated with different APR protocol stages are represented as W_attach, W_pull, W_release, and W_release-std, which stand for the work of attaching the restraints to the bound complex, pulling the guest away from the host, releasing the restraints when the host-guest complex is unbound, and removing the guest rotational restraints and jacks at standard concentration (1 M), respectively. The binding free energy (ΔG_bind) was calculated by adding the previous components with opposite sign.⁶³ Binding enthalpy (ΔH_bind) was predicted as explained in Methodology section, and -TΔS_bind is simply the difference between ΔG_bind and ΔH_bind;
• ^b
  the structures of the complexes chosen for thermodynamic analyses correspond to the last frames of a stable 1 µs MD simulation replicate, and are labelled according to the nomenclature presented in Figures 2-5. Poses that evolved to equivalent final binding modes were analyzed together;
• ^c
  mean values ± SEMs are shown in each case.

As noted in Table 3, all the ΔG_bind values were positive, indicating that encapsulation of the studied drugs was not spontaneous. The same behavior was evident from DFT and semiempirical results. Furthermore, both MD and DFT calculations predicted that all the binding processes were endothermic, except for the slightly exothermic insertion of the diethylamine moiety of procaine in CB[6]. Therefore, we did not expect any interactions with a net stabilizing effect to arise during host-guest binding. On the other hand, note that, except for the CB[6]-caffeine complex in MD calculation, the formation of all the other complexes was entropically unfavorable according to both DFT and MD calculations. Finally, having unfavorable complexes that displayed time stability during at least one of their respective long MD simulation replicates was not contradictory. In fact, the analyzed poses were manually inserted in CB[6], and the eventual exit process had to overcome a large energetic barrier associated with host ring opening. Hence, the systems remained trapped in unstable states during 1 µs MD simulations, which seems to depend on the initial velocities assigned to the atoms, as for almost all systems, except for procaine, we were able to sample dissociation events in that time scale. We expect that, in principle, longer simulation times will lead to the eventual dissociation of all complexes and no encapsulation will be observed afterwards. These simulations, however, are not required because the estimated free energies already predict that the analyzed drugs are not encapsulated by CB[6].

Discussion

The Scientific Working Group for the Analysis of Seized Drugs (SWGDRUG),^4,65 provides the most widely accepted and followed recommendation for forensic illicit drug evaluation. According to them, there are categories of techniques for forensic characterizations (Table 4). The main recommendation is that, when a validated technique belonging to Category A is used, at least one other technique from the other categories should also be employed. Moreover, when a technique belonging to Category A cannot be used, at least three other different techniques must be applied, being two of them of Category B (and not correlated). Color tests are in Category C and they cannot be used uniquely to detect drugs, since they can provide doubtable results (false positives or false negatives). Electrochemical techniques, otherwise, are no longer considered in this recommendation.

Table 4
Techniques recommended by Scientific Working Group for the Analysis of Seized Drugs (SWGDRUG) for drug characterization

Category A	Category B	Category C
Infrared spectroscopy	capillary electrophoresis	color tests
Mass spectrometry	gas chromatography	fluorescence spectroscopy
Nuclear magnetic resonance spectroscopy	ion mobility spectrometry	immunoassay
Raman spectroscopy	liquid chromatography	melting point
X-ray diffractometry	microcrystalline tests	ultraviolet spectroscopy
	pharmaceutical identifiers
	thin layer chromatography
	cannabis only: macroscopic examination and microscopic examination

Cucurbiturils have been used as a material able to encapsulate many substances, among them drugs of abuse. The potential of fluorescent acyclic cucurbituril (aCB) to bind opioids has been investigated.²¹ Marijuana and cocaine detection with piezoelectric devices chemically modified with CB[6] has been studied.²² CB[6]-modified electrodes have been shown to detect MDMA by voltammetry.²³ CB[7] has been applied as a sensor to detect amphetamine-type stimulants.¹⁴ Some of these studies use cucurbituril compounds as an alternative to colorimetric tests. However, there are no studies that make a differential approach to compare drugs with their common adulterants. Therefore, this work intended to study cocaine in a comparative fashion with lidocaine, caffeine, and procaine.

Results for theoretical methods showed that the semiempirical approach did not produce reliable and reproducible results for this specific system. That is why we have decided to improve the information with DFT and MD simulations. Despite there are papers published^23,38 with semi-empirical approach to experimental observation, we have demonstrated that for this specific system semi-empirical is not accurate enough to support laboratory observation.

The results concerning final conformations in DFT and MD calculations provided similar results for cocaine. Calculations converged to the same structures, showing that encapsulation preferentially occurred with the aromatic ring inserted in the cavity. Semi-empirical results were not effective since most of the positions did not reach into convergence (see SI section). For lidocaine, DFT and MD calculations also led to similar results. In this case, the amino group inserted in the CB[6] cavity was the most probable conformation.

DFT and MD calculations gave opposite results for caffeine. Convergence to CAF_P01 and CAF_P02 in the gas phase had similar energies. CAF_P03 was the only conformation possible in the presence of the solvent. For procaine we have similar results for both calculations. For DFT, the final conformation when the molecule was in the center (PROC_P01) resembled the conformation obtained for both PROC_P01 and PROC_P02 by MD.

In terms of electronic energy, the energy was similar for the cocaine final DFT conformations. The results for enthalpies showed that cocaine had similar behavior for both DFT and MD calculations. Lidocaine had positive values for DFT and MD calculations. However, the DFT enthalpic contribution was higher when compared with cocaine. For MD, the opposite trend emerged. For caffeine, the energetic results were not comparable because DFT and MD calculations provided different convergences. For DFT and MD calculations, procaine enthalpy values had opposite trends regarding dimension, but they gave essentially an endothermic behavior (slightly negative for PROC_P01 and PROC_P02 during MD).

As for free energy, both methods, as well as all the converged structures returned positive values, indicating an unfavorable binding process. DFT showed that PROC_P03 had the lowest energy values. All the molecules in the center of the CB[6] cavity displayed the highest energies. MD calculations provided the lowest value for LID_P01 and LID_P02, which gave the same final conformation when lidocaine left the center of the cavity. Except for caffeine results obtained by MD, none of the entropic contributions were favorable.

Conclusions

Investigation is only effective when scientific methods are applied to the real case. No doubt must remain about the nature of the suspicious substance. In this paper, the main objective was to investigate how CB[6] is able to differentiate cocaine and its common adulterants. It is important since there are many studies in literature that consider these substances to detect drugs as an alternative to colorimetric tests. We used in silico methods to study the energetic interaction between CB[6] and cocaine or its adulterants lidocaine, caffeine, and procaine. There was partial agreement between the MD and DFT calculations for cocaine and lidocaine. Nevertheless, the same final conformations were achieved, in which the process was found to be endothermic, not spontaneous, and entropically unfavorable. For caffeine and procaine, the opposite behavior emerged, showing that solvent influence made the difference in these cases. MD calculations provided slightly exothermic value for procaine and a favorable entropic contribution for caffeine. Calculations revealed the aromatic ring inside the CB[6] being the most probable interaction. No calculation showed that cucurbituril was specific for cocaine. Therefore, despite published works on cocaine detection with CB[6], this material should not be considered as a suitable detector for this drug. Overall, these results suggest the inability of CB[6] to encapsulate the drugs studied herein.

Supplementary Information

Supplementary data are available free of charge at https://jbcs.sbq.org.br as PDF file.

Acknowledgments

The authors thank the Brazilian Agencies: Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq, grant 465450/2014-8), Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP, grants 2016/19766-1 and 2016/24587-9) and Coordenação de Aperfeiçoamento de Pessoal de Nível Superior-Brazil (CAPES), finance code 001. The DFT simulations were performed at CENAPAD-SP (Centro Nacional de Processamento de Alto Desempenho em São Paulo). We also thank Cynthia Maria de Campos Prado Manso for language editing and proofreading.

References

• ¹
  United Nations Office on Drugs and Crime (UNODC); Global Overview of Drug Demand and Supply World Drug Report 2019; United Nations: Vienna, 2019, available at https://wdr.unodc.org/wdr2019/prelaunch/WDR19_Booklet_2_DRUG_DEMAND.pdf, accessed in November 2020.
  » https://wdr.unodc.org/wdr2019/prelaunch/WDR19_Booklet_2_DRUG_DEMAND.pdf
• ²
  Penido, C. A. F. O.; Pacheco, M. T. T.; Lednev, I. K.; Silveira, L.; J. Raman Spectrosc. 2016, 47, 28.
• ³
  Kloosterman, A.; Mapes, A.; Geradts, Z.; van Eijk, E.; Koper, C.; van den Berg, J.; Verheij, S.; van der Steen, M.; van Asten, A.; Philos. Trans. R. Soc., B 2015, 370, 20140264.
• ⁴
  Harper, L.; Powell, J.; Pijl, E. M.; Harm Reduct. J. 2017, 14, 52.
• ⁵
  Tsujikawa, K.; Iwata, Y. T.; Segawa, H.; Yamamuro, T.; Kuwayama, K.; Kanamori, T.; Inoue, H.; Forensic Sci. Int. 2017, 270, 267.
• ⁶
  Hondebrink, L.; Nugteren-van Lonkhuyzen, J. J.; van der Gouwe, D.; Brunt, T. M.; Drug Alcohol Depend. 2015, 147, 109.
• ⁷
  Nelson, H. C.; Gardner, E. A.; Matteo, D.; J. Forensic Sci. 2011, 56, 736.
• ⁸
  Marcelo, M. C. A.; Mariotti, K. C.; Ortiz, R. S.; Ferrão, M. F.; Anzanello, M. J.; Microchem. J. 2016, 127, 87.
• ⁹
  Oguri, K.; Wada, S.; Eto, S.; Yamada, H.; Jpn. J. Toxicol. Environ. Health 1995, 41, 274.
• ¹⁰
  Tsumura, Y.; Mitome, T.; Kimoto, S.; Forensic Sci. Int. 2005, 155, 158.
• ¹¹
  Nelson, D. L.; Cox, M. M.; Princípios de Bioquímica de Lehninger, vol. 1, 6^th ed.; Artmed: São Paulo, 2014.
• ¹²
  Walker, S.; Kaur, R.; McInnes, F. J.; Wheate, N. J.; Mol. Pharm. 2010, 7, 2166.
• ¹³
  Walker, S.; Oun, R.; McInnes, F. J.; Wheate, N. J.; Isr. J. Chem. 2011, 51, 616.
• ¹⁴
  Jang, Y.; Jang, M.; Kim, H.; Lee, S. J.; Jin, E.; Koo, J. Y.; Hwang, I.-C.; Kim, Y.; Ko, Y. H.; Hwang, I.; Oh, J. H.; Kim, K.; Chem 2017, 3, 641.
• ¹⁵
  Assaf, K. I.; Nau, W. M.; Chem. Soc. Rev. 2015, 44, 394.
• ¹⁶
  Masson, E.; Ling, X.; Joseph, R.; Kyeremeh-Mensah, L.; Lu, X.; RSC Adv. 2012, 2, 1213.
• ¹⁷
  Pandey, S.; Soni, V. K.; Choudhary, G.; Sharma, P. R.; Sharma, R. K.; Supramol. Chem. 2017, 29, 387.
• ¹⁸
  Cong, H.; Li, C.-R.; Xue, S.-F.; Tao, Z.; Zhu, Q.-J.; Wei, G.; Org. Biomol. Chem. 2011, 9, 1041.
• ¹⁹
  Isaacs, L.; Acc. Chem. Res. 2014, 47, 2052.
• ²⁰
  Lee, J. W.; Samal, S.; Selvapalam, N.; Kim, H.-J.; Kim, K.; Acc. Chem. Res. 2003, 36, 621.
• ²¹
  Shcherbakova, E. G.; Zhang, B.; Gozem, S.; Minami, T.; Zavalij, P. Y.; Pushina, M.; Isaacs, L. D.; Anzenbacher, P.; J. Am. Chem. Soc. 2017, 139, 14954.
• ²²
  de Menezes, M. M. M. T.; Balbino, M. A.; Castro, A. S.; Izabel, C.; Demets, G. J. F.; Oliveira, O. V.; Ipólito, A. J.; McCord, B. R.; Oliveira, M. F.; SM J. Forensic Res. Criminol. 2017, 1, 1003.
• ²³
  Tadini, M. C.; Balbino, M. A.; Eleoterio, I. C.; de Oliveira, L. S.; Dias, L. G.; Demetes, G. J. F.; Oliveira, M. F.; Electrochim. Acta 2014, 121, 188.
• ²⁴
  Xie, X.; Li, X.; Luo, H.; Lu, H.; Chen, F.; Li, W.; J. Phys. Chem. B 2016, 120, 4131.
• ²⁵
  Fenley, A. T.; Henriksen, N. M.; Muddana, H. S.; Gilson, M. K.; J. Chem. Theory Comput. 2014, 10, 4069.
• ²⁶
  Ganji, M. D.; Danesh, N.; RSC Adv. 2013, 3, 22031.
• ²⁷
  Arslanian, A. J.; Dearden, D. V. In Cucurbiturils and Related Macrocycles; Kim, K., ed.; Royal Society of Chemistry: Cambridge, UK, 2019.
• ²⁸
  Hrynchuk, R. D.; Barton, R. J.; Robertson, B. E.; Acta Crystallogr., Sect. A: Found. Crystallogr. 1981, 37, C72.
• ²⁹
  Bambagiotti-Alberti, M.; Bruni, B.; di Vaira, M.; Giannellini, V.; Guerri, A.; Acta Crystallogr., Sect. E: Struct. Rep. Online 2007, 63, o768.
• ³⁰
  Enright, G. D.; Terskikh, V. V.; Brouwer, D. H.; Ripmeester, J. A.; Cryst. Growth Des. 2007, 7, 1406.
• ³¹
  Groom, C. R.; Bruno, I. J.; Lightfoot, M. P.; Ward, S. C.; Acta Crystallogr., Sect. B: Struct. Sci., Cryst. Eng. Mater. 2016, 72, 171.
• ³²
  Hanwell, M. D.; Curtis, D. E.; Lonie, D. C.; Vandermeersch, T.; Zurek, E.; Hutchison, G. R.; J. Cheminf. 2012, 4, 17.
• ³³
  Avogadro: An Open-Source Molecular Builder and Visualization Tool, version 1.XX; available at https://avogadro.cc/, accessed in November 2020.
  » https://avogadro.cc/
• ³⁴
  Rappe, A. K.; Casewit, C. J.; Colwell, K. S. Goddard III, W. A.; Skiff, W. M.; J. Am. Chem. Soc. 1992, 114, 10024.
• ³⁵
  Stewart, J. J. P.; MOPAC 2016; Stewart Computational Chemistry, Colorado Springs, CO, USA, 2016.
• ³⁶
  Rocha, G. B.; Freire, R. O.; Simas, A. M.; Stewart, J. J. P.; J. Comput. Chem. 2006, 27, 1101.
• ³⁷
  Dewar, M. J. S.; Zoebisch, E. G.; Healy, E. F.; Stewart, J. J. P.; J. Am. Chem. Soc. 1985, 107, 3902.
• ³⁸
  Stewart, J. J. P.; J. Mol. Model. 2013, 19, 1.
• ³⁹
  Stewart, J. J. P.; J. Mol. Model. 2007, 13, 1173.
• ⁴⁰
  Stewart, J. J. P.; J. Comput. Chem. 1989, 10, 209.
• ⁴¹
  Becke, A. D.; J. Chem. Phys. 1993, 98, 5648.
• ⁴²
  Weigend, F.; Ahlrichs, R.; Phys. Chem. Chem. Phys. 2005, 7, 3297.
• ⁴³
  Neese, F.; Comput. Mol. Sci. 2012, 2, 73.
• ⁴⁴
  Ramachandran, K. I.; Deepa, G.; Namboori, K.; Computational Chemistry and Molecular Modeling: Principles and Applications; Springer Verlag: Berlin, Heidelberg, 2008.
• ⁴⁵
  Bortoli, M.; Dalla Tiezza, M.; Muraro, C.; Pavan, C.; Ribaudo, G.; Rodighiero, A.; Tubaro, C.; Zagotto, G.; Orian, L.; Comput. Struct. Biotechnol. J. 2019, 17, 311.
• ⁴⁶
  Sinn, S.; Spuling, E.; Bräse, S.; Biedermann, F.; Chem. Sci. 2019, 10, 6584.
• ⁴⁷
  Gentile, D.; Floresta, G.; Patamia, V.; Nicosia, A.; Mineo, P. G.; Rescifina, A.; Org. Biomol. Chem. 2020, 18, 1194.
• ⁴⁸
  de Oliveira, O. V.; Costa, G. C.; Costa, L. T.; J. Phys. Chem. B 2018, 122, 12107.
• ⁴⁹
  Bujak, R.; Gadzała-Kopciuch, R.; Nowaczyk, A.; Raczak-Gutknecht, J.; Kordalewska, M.; Struck-Lewicka, W.; Markuszewski, M. J.; Buszewski, B.; Talanta 2016, 146, 401.
• ⁵⁰
  Li, X.; Zhu, Y.; Liu, C.; Lin, X.; Zhang, W.; Tang, M.; Struct. Chem. 2017, 28, 1631.
• ⁵¹
  Hernández González, J. E.; Hernández Alvarez, L.; Pascutti, P. G.; Valiente, P. A.; Proteins: Struct., Funct., Bioinf. 2017, 85, 1666.
• ⁵²
  Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery Jr., J. A.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J.; Gaussian 09, Revision A.02; Gaussian Inc., Wallingford, CT, 2009.
• ⁵³
  Case, D. A.; Walker, R. C.; Cheatham III, T. E.; Simmerling, C.; Roitberg, A.; Merz, K. M.; Luo, R.; Darden, T.; Wang, J.; Duke, R. E.; Roe, D. R.; LeGrand, S.; Swails, J.; Cerutti, D.; Monard, G.; Sagui, C.; Kaus, J.; Betz, R.; Madej, B.; Lin, C.; Mermelstein, D.; Li, P.; Onufriev, A.; Izadi, S.; Wolf, R. M.; Wu, X.; Götz, A. W.; Gohlke, H.; Homeyer, N.; Botello-Smith, W. M.; Xiao, L.; Luchko, T.; Giese, T.; Lee, T.; Nguyen, H. T.; Nguyen, H.; Janowski, P.; Omelyan, I.; Kovalenko, A.; Kollman, P. A.; Amber 2016; University of California, San Francisco, 2016.
• ⁵⁴
  Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A.; J. Comput. Chem. 2004, 25, 1157.
• ⁵⁵
  https://pubchem.ncbi.nlm.nih.gov/compound/3676, accessed in November 2020.
  » https://pubchem.ncbi.nlm.nih.gov/compound/3676
• ⁵⁶
  Remington’s Pharmaceutical Sciences, 14^th ed.; Osol, A; Hoover, J. E., eds.; Mack Publishing Co.: Easton, 1970.
• ⁵⁷
  Collinsworth, K. A.; Kalman, S. M.; Harrison, D. C.; Circulation 1974, 50, 1217.
• ⁵⁸
  Lu, H.; Chen, X.; Zhan, C.-G.; J. Phys. Chem. B 2007, 111, 10599.
• ⁵⁹
  Li, P.; Zhao, K.; Deng, S.; Landry, D. W.; Helv. Chim. Acta 1999, 82, 85.
• ⁶⁰
  O’Neil, J.; Heckelman, P. E.; Koch, C. B.; Roman, K. J.; J. Am. Chem. Soc. 2007, 129, 2197.
• ⁶¹
  Caldwell, J. W.; Kollman, P. A.; J. Phys. Chem. 1995, 99, 6208.
• ⁶²
  Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.; Hess, B.; Lindahl, E.; SoftwareX 2015, 1-2, 19.
• ⁶³
  Henriksen, N. M.; Fenley, A. T.; Gilson, M. K.; J. Chem. Theory Comput. 2015, 11, 4377.
• ⁶⁴
  DeLano, W. L.; PyMOL; DeLano Scientific, San Carlos, 2002.
• ⁶⁵
  Scientific Working Group for the Analysis of Seized Drugs (SWGDRUG); SWGDRUG Recommendations, version 7.1; United States Department of Justice, Washington DC, USA, 2016.

Publication Dates

• Publication in this collection
  29 Mar 2021
• Date of issue
  Apr 2021

History

• Received
  17 July 2020
• Accepted
  20 Nov 2020