THC sensing from breath is challenging because it is present in a very small concentration in vapor phase. Successful breath assays of THC have required that picogram amounts of THC be collected from aerosol droplets in the breath that are expelled when the subject breathes out.Other complicating factors of THC detection include its fast metabolism in the body and the other chemically similar compounds that are present in cannabis, such as cannabinodiol and tetrahydrocannabinolic acid , a metabolyte of THC. Efforts to detect THC from breath have largely relied on mass spectrometry and, more recently, a fluorescence-based assay.In this work we propose using metal-organic frameworks to adsorb THC prior to detection, which can be carried out using existing optical,electrochemical, or fluorescence-based sensing methods. MOFs are porous, crystalline materials composed of metal nodes linked by organic ligands in a three-dimensional network.They have a high surface area for adsorption and a high degree of chemical tunability that makes possible the synthesis of a framework with ideal properties for a given application. The huge design space of MOFs is reflected in the large and growing number of experimentally synthesized MOFs that have been reported to date, nearing 100 000 structures.With this chemical diversity of structures, MOFs have exhibited unmatched performance in applications such as adsorption, catalysis, and sensing.MOFs have been used as a detection medium for biologically relevant molecules such as chemical warfare agents and toxic industrial chemicals.Sensing in biological systems poses particular challenges for MOFs,multi tier vertical growing including interference of water with the adsorbate of interest, and the thermal and chemical stability of MOFs in the presence of bodily fluids such as blood or saliva.
One solution is to use hydrophobic materials that would minimize adsorption of water in order to detect the analyte or to adsorb water from the sample before it comes in contact with the MOF. More generally, detection of larger biologically-relevant molecules in MOFs is interesting because the adsorption energy of larger molecules may be more dependent on pore shape than adsorption of small gas molecules such as N2 and CO2 . This allows for potentially higher selectivity when the pore shape can be optimized for an adsorbate molecule. Simulation of larger molecules can also pose computational challenges. Properties such as diffusion and adsorption become more expensive to compute, and the difference in computational time can be significant, especially when screening a library of structures for a particular application. Strategies to reduce computational cost while still capturing the decisive properties for an application can be useful. Toward this end, we tested two models of THC of different complexity. With these adsorbate models, we investigated 645 MOF structures that have been filtered from a large set 5171 materials, that combine the CoRE MOF database and the MOF-74 analog database. Geometric analysis of the the porous frameworks was carried out with the existing Zeo++ package and a novel method using 2D projections of the adsorbate and pores to analyze pore accessibility toward asymmetric probe molecules. We used NVT and Widom insertion simulations to compute the adsorption energy and Henry coefficient of THC in these structures, and to determine whether using a simplified model would be an appropriate screening strategy for THC.
We identified motifs of adsorption in the most promising structures. Finally we computed the THC/H2O selectivity in these materials and determined which might be suitable for THC detection in the presence of water. Short-range dispersion interactions had a cutoff distance of 12 Å without applying tail corrections. The unit cells of all frameworks were replicated so that perpendicular widths between opposite faces were larger than twice the cutoff distance. The positions of the framework atoms were kept fixed for all simulations. In this study, two models of THC are used: the THC-full model and the THC-head model. The THC-head model includes only the bulky head group of the THC molecule, and is fully rigid. The THC-full model treats the aromatic head of the THC molecule as rigid, but also includes the hydrocarbon tail, which is flexible. We used Monte Carlo methods to screen adsorption of THC in our materials library using both THC-head and THC-full models. The tail of the THC-full model was simulated using configurational-bias Monte Carlo .We ran Widom insertion simulations in the desolvated framework and NVT simulations with a loading of 1 molecule per simulation cell. In NVT simulations we chose moves with equal probability of translation, rotation, reinsertion, and reinsertion-in-place, where the last is a CBMC move that applies to the THC-full model. For THC, we used blocking spheres generated in Zeo++ with a probe radius of 2.5 Å. This is smaller than the size of the THC molecule, and was chosen conservatively in order not to exclude pores which might be accessible to THC considering a minimal vibration of the crystal. A more detailed method to determine pore accessibility is detailed in the following section. MC simulations were run in the RASPA package,considering a temperature of 309 K, the human body temperature. Our models and simulations were first tested in depth on a single framework, the well studied Mg-MOF-74.
This preliminary analysis led to the conclusions that Coulombic interactions can be neglected for the investigation of THC adsorption, and the need of 10 times more MC Widom insertions than MC NVT steps to achieve a comparable convergence. We report the full analysis on Mg-MOF-74 in the Supporting Information. To check for selectivity of THC adsorption in the presence of water vapor, we computed the Henry coefficient of water using MC Widom insertions following the protocol ofMoghadam et al.Point charges for these materials were derived using the EQeq method,that we recently benchmarked.To model water, we used the TIP4P model with a unitedatom Lennard-Jones interaction site. We computed KH for water in these materials using 1 000 000 Widom insertion steps and blocking spheres computed for a probe radius of 1.5 Å. Since we are including MOFs which have a coordinatively unsaturated metal site one has to consider that modelling interactions using only Coulombic and dispersion interactions often leads to an underestimation of the binding of polar gas molecules at these sites. Our group spent a considerable effort to address this issue, by proposing tailor made corrections to tune the force field for some specific gas-OMS interactions.However, a reliable and systematic correction for all the combinations of polar molecules and OMSs has not yet been proposed, and therefore it is of common practice to take the results from standard interactions with a grain of salt and analyze overall trends over a large set of structures.We are interested in restricting our screening only to those materials through which THC can percolate, as well as determining whether THC can pass through a particular channel in a material. Most geometric methods to evaluate pore accessibility typically use idealized spherical probe molecules for their analyses. For example, the pore limiting diameter indicates the size of the largest sphere which can percolate through the material,rolling vertical grow racks giving an idea of the barriers to diffusion. This is useful for adsorbates which are approximately spherical or at least have a circular cross-section, however, THC is asymmetric from all view points and is therefore characterized by multiple characteristic lengths: As seen in Figure 3.2, when the THC molecule is viewed down its longest axis, it has a cross section with dimensions of roughly 10 Å and 7 Å along its longer and shorter axes, respectively. Therefore, THC should be able to traverse channels with a PLD greater than its longest dimension, while it would not be able to fit through a channel whose PLD is smaller than its shortest dimension. For intermediate pore sizes, accessibility depends on the precise pore geometry as compared with the geometry of the probe molecule. We developed a method which can be used to determine whether an asymmetric molecule can percolate through a porous material, by comparing the shape of the probe molecule to the cross-sectional geometry at the narrowest part of the pore. The algorithm, which is described in detail in the Methods section, computes the minimum area of overlap between the probe and the pore. We present our results applying this algorithm to our screening database in this section.The mmen-Mg2 MOF and its expanded analogs. In particular, characterizing the nature of the amine dynamics in the absence of CO2 can provide insight into design rules for materials for optimal CO2 capture performance under different conditions.
A central feature of the amine-appended MOFs is the lability of the coordinate bond between the amine functional groups and the metal nodes. NMR spectroscopy has indicated that for mmen-Mg2 at timescales relevant for CO2 capture, there is only one 15N peak corresponding to the functional group N,N’-dimethylethylenediamine . Since only one amine group can be coordinated to the metal while the other is free, the evidence that there is only one chemical environment for 15N suggests that there is an exchange of amine group coordination occuring faster than the timescale of NMR.27 Most classical force field models assume that bond networks remain fixed throughout the course of a simulation, which would not be appropriate for this system. On the other end of the spectrum, quantum mechanical methods do not assume a bonding network, but the computational cost of such methods is prohibitive for modeling a system on the length and timescales necessary to study amine functional group dynamics. We therefore turn to reactive force field methods, which incorporate a bond order term in the total energy function of the system, where the bond order is calculated at each time step based on the pairwise distances between atoms. Many reactive force field methods exist, some which have been expanded for simulation of a variety of systems, such as ReaxFF, and others which are designed for a class of materials, such as REBO and AIREBO for carbon/hydrogen systems.We chose the ReaxFF force field method because it has been used widely for both inorganic and organic systems.The ReaxFF parameters are optimized using GARF field program developed by Andres Jaramillo-Botero.GARFfield takes a training set of geometries with energies calculated by DFT and conducts a genetic algorithm search to match the ReaxFF energies to those calculated using DFT. The training set used to optimize the Gen 2 force field consisted of a set of unit cell expansions and contractions of the bare Mg2 MOF in order to improve the stability of the framework on its own. The Gen 3 force field, still under development, uses a training set that comes from a ReaxFF simulation of the entire mmen-Mg2 MOF using the Gen 2 force field, which allows the amines to adopt different configurations so that the interaction energies between the amines and the framework atoms can be tuned. The activation enthalpy and entropy for mmen inter conversion are kJ mol−1 and J mol−1 K−1 , respectively, where the error represents a 95% confidence interval. The enthalpy and entropy measured via NMR by experimental collaborators are3.6 kJ mol−1 and −82 J mol−1 K−1 . The simulated and measured activation enthalpies are close to within the uncertainty range. However, the activation entropies are very different, indicating that while the rates scale similarly with temperature, the simulated reaction rate is orders of magnitude larger than the experimental reaction rate according to the Eyring-Polanyi equation . The difference in activation entropy for mmen is attributed to the excess mobility of the simulated framework as compared to the experimental system. As noted, the 1st generation force field allows for a notable deformation of the pore volume and shape while the system is equilibrating, and in the equilibrated system there is dynamic buckling of the linkers, changing the coordination environment around the metal nodes. This altered coordination environment appears to make the activated state of the “hopping” inter conversion more entropically favored in simulation compared to experiment. In both simulation and experiment, however, the transition state has lower entropy than either the reactant and product state. This supports the proposed reaction mechanasm as involving a transition state where both ends of the diamine are coordinated to the metal, and therefore have restricted mobility compared to when one amine end is “free.” The closeness in values of the activation enthalpy also suggests that the activation complex observed in simulation is the same as seen in experiment.