Molecular dynamics simulations of a coarse-grained model are used to study the formation mechanism of periodic mesoporous silica over a wide range of cationic surfactant concentrations. This follows up on an earlier study of systems with low surfactant concentrations. We started by studying the phase diagram of the surfactant–water system and found that our model shows good qualitative agreement with experiments with respect to the surfactant concentrations where various phases appear. We then considered the impact of silicate species upon the morphologies formed. We have found that even in concentrated surfactant systems—in the concentration range where pure surfactant solutions yield a liquid crystal phase—the liquid-crystal templating mechanism is not viable because the preformed liquid crystal collapses as silica monomers are added into the solution. Upon the addition of silica dimers, a new phase-separated hexagonal array is formed. The preformed liquid crystals were found to be unstable in the presence of monomeric silicates. In addition, the silica dimer is found to be essential for mesoscale ordering at both low and high surfactant concentrations. Our results support the view that a cooperative interaction of anionic silica oligomers and cationic surfactants determines the mesostructure formation in the M41S family of materials.

Platinum-based materials play an important role as electrocatalysts in energy conversion technologies. Graphene-supported Pt nanoclusters were recently found to be promising electrocatalysts for fuel-cell applications due to their enhanced activity and tolerance to CO poisoning as well as their long-term stability toward sintering. However, structure–function relationships that underpin the improved performance of these catalysts are still not well understood. Here, we employ a combination of empirical potential simulations and density functional theory (DFT) calculations to investigate structure–function relationships of small PtN (N = 2–80) clusters on model carbon (graphene) supports. A bond-order empirical potential is employed within a genetic algorithm to go beyond local optimizations in obtaining minimum-energy structures of PtN clusters on pristine as well as defective graphene supports. Point defects in graphene strongly anchor Pt clusters and also appreciably affect the morphologies of small clusters, which are characterized via various structural metrics such as the radius of gyration, average bond length, and average coordination number. A key finding from the structural analysis is that the fraction of potentially active surface sites in supported clusters is maximized for stable Pt clusters in the size range of 20–30 atoms, which provides a useful design criterion for optimal utilization of the precious metal. Through selected ab initio studies, we find a consistent trend for charge transfer from small Pt clusters to defective graphene supports resulting in the lowering of the cluster d-band center, which has implications for the overall activity and poisoning of the catalyst. The combination of a robust empirical potential-based genetic algorithm for structural optimization with ab initio calculations opens up avenues for systematic studies of supported catalyst clusters at much larger system sizes than are accessible to purely ab initio approaches.

We have performed replica-exchange reaction ensemble Monte Carlo simulations to study the low-energy crystalline structures of a reactive model of silica. We have utilized a model of silica polymerization based on the reactive assembly of semiflexible tetrahedral units developed by us previously to reproduce silica bulk moduli as well as self-assembly of amorphous silica gels and nanoparticles. Our implementation of replica-exchange Monte Carlo involves simulating several system copies, each with its own value of the equilibrium constant controlling silica condensation/hydrolysis reactions, which are essential for building higher-order network structures and eventually crystals. These replica-exchange simulations were found to traverse energy landscapes from amorphous to crystalline phases, yielding the dense silica polymorphs α-cristobalite, β-cristobalite, and keatite, as well as the nanoporous silica materials SOD and EDI and nanoporous phosphates with DFT and ATT structures. Simulated crystal structures were confirmed by computing X-ray patterns for comparison with known XRD data. The behavior of this model opens the door to future simulation studies of the free energy barriers controlling these crystallization processes.

We use parallel tempering Monte Carlo simulations to search for crystalline states of a lattice model of silica polymerization in the presence of structure directing agents (SDAs). Following previous work where we have discretized continuous space into a body-centered cubic (bcc) lattice, we have modeled tetrahedral molecules (T(OH)4) as corner-sharing tetrahedra on a bcc unit cell. The SDAs were represented as quasi-spherical species with diameters of 6.4 and 10.4 Å to study the effect of SDA size on the resulting crystal structures. Our parallel tempering Monte Carlo simulations produce fully connected crystalline structures finding the emergence of 3D microporous materials with SDAs occupying the pore spaces and 2D layered materials with SDAs occupying the gallery space in between layers. We have found that the strength of SDA–oxygen attraction plays a significant role in directing final micropore structures. For relatively strong attractions (>1.2 kcal/mol SDA–oxygen contacts) we have found only 2D layered materials; for attractions below this cutoff we observed 3D microporous crystals; and for no attraction—modeling the SDA as a quasi-hard sphere—we again found only 2D layered materials. In the space of 3D microporous crystals, we have also found that using larger SDAs or a lower concentration of a given SDA generate crystals with larger rings.

The reaction ensemble Monte Carlo method was used to model the self-assembly and structure of silica nanoparticles found in the initial stages of the clear-solution synthesis of the silicalite-1 zeolite. Such nanoparticles, which comprise both silica and organic structure-directing agents (OSDAs), are believed to play a crucial role in the formation of silica nanoporous materials, yet very limited atomic-level structural information is available for these nanoparticles. We have modeled silica monomers as flexible tetrahedra with spring constants fitted in previous work to silica bulk moduli and OSDAs as spheres attracted to anionic silica monomers. We have studied one-step and two-step formation mechanisms, the latter involving the initial association of silica species and OSDAs driven by physical solution forces, followed by silica condensation/hydrolysis reactions simulated with reaction ensemble Monte Carlo. The two-step process with preassociation was found to be crucial for generating nearly spherical nanoparticles; otherwise, without preassociation they exhibited jagged, ramified structures. The two-step nanoparticles were found to exhibit a core–shell structure with mostly silica in the core surrounded by a diffuse shell of OSDAs, in agreement with SANS and SAXS data. The Qn distribution, quantifying silicon atoms bound to n bridging oxygens, found in the simulated nanoparticles is in broad agreement with 29Si solid-state NMR data on smaller, 2 nm nanoparticle populations. Ring-size distributions from the simulated nanoparticles show that five-membered rings are prevalent when considering OSDA/silica mole fractions (∼0.2) that lead to silicalite-1, in agreement with a previous IR and modeling study. Nanoparticles simulated with higher OSDA concentrations show ring-size distributions shifted to smaller rings, with three-membered silica rings dominating at an OSDA/silica mole fraction of 0.8. Our simulations show no evidence of long-range silicalite-1 order in these nanoparticles.

Ab initio and classical simulations were used to study the equilibrium and fluctuating ring diameters for all-silica zeolites SOD, FER, and MFI over the temperature range 300–900 K. Such simulations are important for understanding and predicting zeolite/guest fit, especially for relatively bulky guest species, e.g., those derived from biomass. We simulated equilibrium zeolite structures, IR spectra, thermal expansion coefficients, and ring breathing vibrations to investigate the competition between negative thermal expansion and enhanced vibrational amplitudes with increasing temperature. We find that although negative thermal expansion tends to shrink equilibrium ring sizes with increasing temperature, this trend is nullified by considering ring breathing vibrations, giving effective pore sizes that are roughly constant with temperature, and larger than those extracted from X-ray data. Several force fields were tested and a modified BKS force field was found to give the best agreement with the simulated properties listed above, especially for MFI. Our results are consistent with previous work suggesting that effective zeolite ring sizes are underestimated by using oxygen ionic radii for estimating atomic excluded volume.

We simulated structural and dynamical properties of imidazoles tethered to aliphatic backbones to determine how chain length influences the competition between extended hydrogen-bond networks and imidazole reorientation dynamics. We performed molecular dynamics simulations on hypothetical solids using the GAFF Amber force field over the temperature range 300–800 K, for chain lengths varying from monomers to pentamers. We investigated the effect of heterogeneity by simulating monodisperse and polydisperse solids with the same average chain length. We computed hydrogen-bond cluster sizes and percolation ratios; orientational order parameters associated with imidazole rings, tethering linkers, and backbones; and orientational correlation functions for imidazole rings. We found the surprising result that chain-length heterogeneity negligibly affects system density at standard pressure, the fraction of percolating hydrogen-bonded clusters, and the order parameters for backbone, linker, and imidazole ring. Decreasing oligomer chain length from pentamers to shorter chains decreases the tendency to form percolating hydrogen-bond networks while dramatically decreasing imidazole ring reorientation times from a broad range of 100–700 ps for pentamers down to 20 ps for monomers, hence quantifying the competition between hydrogen-bond cluster size and reorientation rate. The computed orientational order parameters suggest the following hierarchy of structural excitations: imidazole ring reorientation in the range 400–500 K, linker motion around 500–600 K, and backbone relaxation at 600–700 K in this model. The question remains for this class of systems which of these motions is crucial for facile proton transport.

We present Monte Carlo simulations of a lattice model describing silica polymerization with an emphasis on the transition between gel states and nanoparticle states as the pH and silica concentration are varied. The pH in the system is controlled by the addition of a structure-directing agent (SDA) of the type SDA+(OH–). The silica units are represented by corner-sharing tetrahedra on a body-centered cubic lattice and the SDA+ species by single sites with near-neighbor repulsions. We focus on two systems: one with a low silica concentration with composition comparable to that of the clear solution silicalite-1 zeolite synthesis and a high silica concentration system that leads to gel states. In the dilute system, clusters have a core–shell structure, with the core predominantly comprised of silica with some SDA+ cations, surrounded by a shell of only SDA+ cations. Moreover, the average cluster size gradually decreases from 2 to 1.6 nm with increasing pH. The concentrated system forms a gel that remains stable to increasing pH up to about 9.2. At pH values in the range of 9.2–10, the gel transforms to nanoparticles of size around 1.0 nm, surprisingly smaller than those in the dilute system. We also study the evolution of the Qn distribution (a measure of the silica network structure) for both systems and obtain good agreement with 29Si NMR data available for the concentrated system.

We applied density functional theory to investigate the mixed aldol condensation of acetone and formaldehyde in acid zeolites HZSM-5 and HY, as a prototypical bond-forming reaction in biofuel production. We modeled the acid-catalyzed reaction in HZSM-5 and HY in two steps: keto–enol tautomerization of acetone and bimolecular condensation between formaldehyde and the acetone enol. For both acid zeolites, the keto–enol tautomerization of acetone was found to be the rate-determining step, consistent with the accepted mechanism in homogeneous acid-catalysis. Convergence studies of the activation energy for keto–enol tautomerization, with respect to cluster sizes of HZSM-5 and HY, exhibit rather different convergence properties for the two zeolites. The keto–enol activation energy was found to converge in HY to ∼20 kcal/mol for a cluster with 11 tetrahedral atoms (11T cluster), which does not complete the HY supercage. In contrast, the activation energy for HZSM-5 reaches an initial plateau at a value of ∼28 kcal/mol for clusters smaller than 20T and then converges to ∼20 kcal/mol for clusters of size 26T or greater, well beyond the completion of the HZSM-5 pore. As such, completing a zeolite pore surrounding a Brønsted acid site may be insufficient to converge activation energies; instead, we recommend an approach based on converging active-site charge.

We performed kinetics experiments and quantum calculations to investigate the reaction of furan to benzofuran catalyzed by the acidic zeolite HZSM-5, which is a key step in the conversion of biomass to biofuels through catalytic fast pyrolysis. The reaction was studied experimentally by placing the zeolite in contact with solution-phase furan and detecting the benzofuran product over the temperature range 270–300 °C, yielding an apparent activation energy of 72 ± 3 kJ/mol. The reaction was modeled in gas and zeolite phases to determine the energetics of the following two competing pathways: a Diels–Alder mechanism often assumed in interpretations of experimental data and a ring-opening pathway predicted by the chemoinformatic software RING. Quantum calculations on the zeolite/guest system were performed using the ONIOM embedded cluster approach. We computed the energetics of reactants, products, and all intermediate steps. Locating relevant transition states fell beyond our computational resources because of system size and the ruggedness of the energy landscape. The Diels–Alder mechanism in the gas phase was found to pass through a high-energy intermediate roughly 380 kJ/mol above the reactant energy, which reduces to approximately 200 kJ/mol in HZSM-5. In contrast, the ring-opening mechanism passes through a gas-phase intermediate roughly 500 kJ/mol above the reactant energy, which falls to approximately 50 kJ/mol in HZSM-5. The energy of the ring-opening mechanism over HZSM-5 fits into the experimentally determined energy “budget” of 72 ± 3 kJ/mol. These experimental and computational results highlight the importance of the ring-opening mechanism for this key step in making biofuels. Our results strongly indicate that, in the cavities of HZSM-5, the condensation of two furan molecules to form benzofuran and water does not proceed by a Diels–Alder reaction between the reactants.Keywords: biofuels; catalytic pyrolysis; Diels−Alder; furan; HZSM-5; ONIOM; QM/MM; zeolite catalysis

We have constructed and applied an algorithm to simulate the behavior of zeolite frameworks during liquid adsorption. We applied this approach to compute the adsorption isotherms of furfural–water and hydroxymethyl furfural (HMF)–water mixtures adsorbing in silicalite zeolite at 300 K for comparison with experimental data. We modeled these adsorption processes under two different statistical mechanical ensembles: the grand canonical (V–Nz–μg–T or GC) ensemble keeping volume fixed, and the P–Nz–μg–T (osmotic) ensemble allowing volume to fluctuate. To optimize accuracy and efficiency, we compared pure Monte Carlo (MC) sampling to hybrid MC–molecular dynamics (MD) simulations. For the external furfural–water and HMF–water phases, we assumed the ideal solution approximation and employed a combination of tabulated data and extended ensemble simulations for computing solvation free energies. We found that MC sampling in the V–Nz–μg–T ensemble (i.e., standard GCMC) does a poor job of reproducing both the Henry’s law regime and the saturation loadings of these systems. Hybrid MC–MD sampling of the V–Nz–μg–T ensemble, which includes framework vibrations at fixed total volume, provides better results in the Henry’s law region, but this approach still does not reproduce experimental saturation loadings. Pure MC sampling of the osmotic ensemble was found to approach experimental saturation loadings more closely, whereas hybrid MC–MD sampling of the osmotic ensemble quantitatively reproduces such loadings because the MC–MD approach naturally allows for locally anisotropic volume changes wherein some pores expand whereas others contract.

We have used Monte Carlo simulations to study the formation of the MCM-41 mesoporous silica material, with a new lattice model featuring explicit representations of both silicic acid condensation and surfactant self-assembly. Inspired by experimental syntheses, we have adopted the following two-step “synthesis” during our simulations: (i) high pH and low temperature allowing the initial onset of mesostructures with long-range order; (ii) lower pH and higher temperature promoting irreversible silica condensation. During step (i), the precursor solution was found to spontaneously separate into a surfactant–silicate-rich phase in equilibrium with a solvent-rich phase. Lamellar and hexagonal ordering emerged for the surfactant–silicate-rich mesosphases under different synthesis conditions, consistent with experimental observations. Under conditions where silica polymerization can be neglected, our simulations were found to transform reversibly between hexagonal and lamellar phases by changing temperature. During step (ii), silica polymerization was simulated at lower pH using reaction ensemble Monte Carlo to treat the pH dependence of silica deprotonation equilibria. Monte Carlo simulations produced silica–surfactant mesostructures with hexagonal arrays of pores and amorphous silica walls, exhibiting Qn distributions in reasonable agreement with 29Si NMR experiments on MCM-41. Compared with bulk amorphous silica, the wall domains of these simulated MCM-41 materials are characterized by even less order, larger fractions of 3- and 4-membered rings, and wider ring-size distributions.

We modeled nascent decomposition processes in cellulose pyrolysis at 327 and 600 °C using Car–Parrinello molecular dynamics (CPMD) simulations with rare events accelerated with the metadynamics method. We used a simulation cell comprised of two unit cells of cellulose Iβ periodically repeated in three dimensions to mimic the solid cellulose. To obtain initial conditions at reasonable densities, we extracted coordinates from larger classical NPT simulations at the target temperatures. CPMD-metadynamics implemented with various sets of collective variables, such as coordination numbers of the glycosidic oxygen, yielded a variety of chemical reactions such as depolymerization, fragmentation, ring opening, and ring contraction. These reactions yielded precursors to levoglucosan (LGA)—the major product of pyrolysis—and also to minor products such as 5-hydroxy-methylfurfural (HMF) and formic acid. At 327 °C, we found that depolymerization via ring contraction of the glucopyranose ring to the glucofuranose ring occurs with the lowest free-energy barrier (20 kcal/mol). We suggest that this process is key for formation of liquid intermediate cellulose, observed experimentally above 260 °C. At 600 °C, we found that a precursor to LGA (pre-LGA) forms with a free-energy barrier of 36 kcal/mol via an intermediate/transition state stabilized by anchimeric assistance and hydrogen bonding. Conformational freedom provided by expansion of the cellulose matrix at 600 °C was found to be crucial for formation of pre-LGA. We performed several comparison calculations to gauge the accuracy of CPMD-metadynamics barriers with respect to basis set and level of theory. We found that free-energy barriers at 600 °C are in the order pre-LGA < pre-HMF < formic acid, explaining why LGA is the kinetically favored product of fast cellulose pyrolysis.

The potential of tailored nanopores to transform technologies such as drug delivery, biofuel production, and optical-electronic devices depends on fundamental knowledge of the self-assembly of ordered nanoporous solids. Atomic-level geometries of critical nuclei that lead to such solids have remained hidden in the nanoscale blind spot between local (<0.5 nm) and collective (>5 nm) probes of structure. Heroic efforts at molecular simulation of nanopore formation have provided massive libraries of hypothetical structures;(1-5) however, to date no statistical simulation has generated a crystallization pathway from random initial condition to ordered nanoporous solid, until now. In this work, we show that a recently developed atomic lattice model of silica and related materials can form ordered nanoporous solids with a rich variety of structures including known chalcogenides, zeolite analogs, and layered materials. We find that whereas canonical Monte Carlo simulations of the model consistently produce the amorphous solids studied in our previous work, parallel tempering Monte Carlo gives rise to ordered nanoporous solids. The utility of parallel tempering highlights the existence of barriers between amorphous and crystalline phases of our model. Moreover, the self-assembly or nanoporous crystalline phases in the model open the door to detailed understanding of nanopore nucleation.Keywords: crystalline; nanoporous materials; parallel tempering Monte Carlo; self-assembly;

We have studied base strengths of nitrogen-substituted (nitrided) zeolites with faujasite (FAU) structure by calculating sorption energies of probe molecules (BF3 and BH3) using density functional theory with mixed basis sets applied to embedded clusters. BH3 was found to be a better probe of base strength because it does not introduce competing metal−fluorine interactions that obfuscate trends. In all cases, the base strengths of nitrided zeolites (denoted M−N−Y) were found to exceed those of the corresponding standard M−Y zeolites, where M = Li, Na, K, Rb, or Cs charge-compensating cations. We have found that for a particular Si:Al ratio, BH3 sorption energies vary in the order Li < Na < K ∼ Rb ∼ Cs. Sorption energy and hence base strength was found to decrease with increasing Si:Al ratio from 1 to 3 beyond which the base strength was found to increase again. The initial regime (1 < Si:Al < 3) is consistent with the prevailing understanding that the base strength increases with Al content, while the latter regime (Si:Al > 3) involves the surprising prediction that the base strength can be relatively high for the more stable, high-silica zeolites. In particular, we found the sorption energy in Na−N−Y (Si:Al = 11) to be nearly equal to that in (Si:Al = 1). Taken together, these results suggest that K−N−Y (Si:Al = 11) optimizes the balance of activity, stability, and cost.

We have modeled structures and energetics of anhydrous proton-conducting wires: tethered hydrogen-bonded chains of the form ···HX···HX···HX···, with functional groups HX = imidazole, triazole, and formamidine; formic, sulfonic, and phosphonic acids. We have applied density functional theory (DFT) to model proton wires up to 19 units long, where each proton carrier is linked to an effective backbone to mimic polymer tethering. This approach allows the direct calculation of hydrogen bond strengths. The proton wires were found to be stabilized by strong hydrogen bonds (up to 50 kJ/mol) whose strength correlates with the proton affinity of HX [related to pKb(HX)] and not to pKa(HX) as is often assumed. Geometry optimizations and ab initio molecular dynamics near 400 K on imidazole-based proton wires both predict that adding a proton to the end of such wires causes the excess charge to embed into the interior segments of these wires. Proton translocation energy landscapes for imidazole-based wires are sensitive to the imidazole attachment point (head or feet) and to wire architecture (linear or interdigitated). Linear imidazole wires with head-attachment exhibit low barriers for intrawire proton motion, rivaling proton diffusion in liquid imidazole. Excess charge relaxation from the edge of wires is found to be dominated by long-range Grotthuss shuttling for distances as long as 42 Å, especially for interdigitated wires. For imidazole, we predict that proton translocation is controlled by the energetics of desorption from the proton wire, even for relatively long wires (600 imidazole units). Proton desorption energies show no correlation with functional group properties, suggesting that proton desorption is a collective process in proton wires.

We present a new model and method for the Monte Carlo simulation of silica polymerization in aqueous solution. We focus on the idea that silica structures are built from corner sharing tetrahedra and these tetrahedra are the basic units of the model. Rather than use a reactive force field, the assembly of tetrahedral units is accomplished via Monte Carlo simulation in the reaction ensemble. The simplicity of the model and the use of the reaction ensemble make it possible to study silica polymerization for quite large system sizes, reaching a high degree of condensation under ambient conditions. We find that the reaction ensemble Monte Carlo simulation protocol can provide a description of the overall polymerization kinetics, after making some key assumptions. Very good agreement is obtained when comparing simulated and experimental evolutions of the Qn distribution as a function of both time and degree of condensation, indicating an approximately linear relationship between physical time and number of Monte Carlo steps up to about 5 h. Analyses of cluster-size and ring-size distributions reveal that polymerization proceeds in the following main stages: oligomerization forming small units (0–1 h), ring formation (1–2.6 h), cluster aggregation (2.6–5.6 h), and finally cross-linking of the aging gel at later times.

We examine the effect of synthesis conditions on the degree of substitution and retention of crystallinity and microporosity during ammonia treatment of Y zeolite. Our objective is to find a reproducibly optimal synthesis protocol to make nitrogen substituted (nitrided) zeolites. We find that a temperature of approximately 750 °C at a mean ammonia flow rate of 600 cm3/min produces substitution without loss of crystallinity or microporosity. We have investigated the effects of several synthesis parameters using fifteen different synthesis protocols. The most important parameter in the synthesis is the ammonia flow rate; we recommend keeping the flow rate as high as possible. Calculations of the NMR and vibrational spectra of substituted zeolites are used to investigate the overall utility of these techniques for confirming the presence of nitrogen in the zeolite. Infrared and Raman spectroscopy were not found to be diagnostic, but they can corroborate the presence of nitrogen in the framework. We recommend a combination of 29Si MAS NMR spectroscopy, X-ray diffraction, and high-resolution adsorption to test any new reactor design; these techniques can establish nitrogen substitution in the framework, as well as crystallinity and microporosity in the product.

We present an extension of the simple-cubic lattice model developed by Jorge et al. [ J. Am. Chem. Soc. 2005, 127, 14388] of nanoparticle growth in the clear solution synthesis of silicalite-1 (MFI). We have implemented the model on a body-center cubic (bcc) lattice with second-neighbor repulsions, to generate a four-coordinate network that mimics the tetrahedral structure of silica. With this low-coordination lattice model we observe that the nanoparticles are metastable, possessing a core−shell structure with mostly silica in the core and templates forming a shell. Nanoparticle size is found to increase with temperature and decrease with solution pH, in qualitative agreement with results from experiment and the previous lattice model study. The low-coordination model makes it possible to model porosity in the silica core of nanoparticles. We use this feature to investigate the extent of template penetration into the silica core, a level of nuance missing in experimental data on the core−shell model. We find that template penetration is rare for bulky templates. We discuss the implications of this result for understanding the role of these nanoparticles in the growth of MFI, especially in light of recent experiments on the long-time behavior of nanoparticle suspensions.

We have developed a model for silica polymerization at ambient temperatures and low densities and have studied this using reactive Monte Carlo simulations. The model focuses on SiO4 coordination with the energetics of hydrolysis and condensation reactions treated via the reaction ensemble. The simplicity of the model makes large system sizes accessible on a modest computation budget, although it is necessary to make additional assumptions in order to use the reactive Monte Carlo method as a simulation of the system dynamics. Excellent agreement for the evolution of the Qn distribution is obtained upon comparing the simulation results to experimental observations. The analysis of simulation trajectories provides mechanistic insight into the polymerization process, showing the following three regimes: oligomerization (0−1 h), ring formation (1−2.6 h), and cluster aggregation (2.6−5.6 h).Keywords (keywords): Monte Carlo simulations; polymerization kinetics; silica gel; silica polymerization; sol−gel process;

We have developed a new grand canonical molecular dynamics (GCMD) algorithm to study microwave (MW) heating effects on competitive mixture sorption and have applied the method to methanol and benzene in silicalite zeolite. The new algorithm combines MW-driven molecular dynamics with grand canonical Monte Carlo (GCMC), the latter modeling adsorption/desorption processes. We established the validity of the new algorithm by benchmarking single-component isotherms for methanol and benzene in silicalite against those obtained from standard GCMC, as well as against experimental data. We simulated single-component and mixture adsorption isobars for conventional and MW-heated systems. In the case of the single-component isobars, we found that for dipolar methanol, both the MW and conventional heated isobars show similar desorption behavior, displaying comparable loadings as a function of molecular temperature. In contrast, nonpolar benzene showed no desorption upon exposure to MWs, even for relatively high field strengths. In the case of methanol/benzene mixtures, the fact that benzene is transparent to the MW field allows the selective desorption of methanol, giving rise to loading ratios not reachable through conventional heating.

We have modeled permeation through anisotropic zeolite membranes with nanoscopic defects that create shortcuts perpendicular to the transmembrane direction (x). We have found that the dimensionless ratio Dy/(kdΔy) can be used to estimate whether the shortcuts contribute significantly to the overall flux. Here Dy is the diffusion coefficient for motion in the plane of the membrane, kd is the rate of desorbing into defect voids, and Δy is the spacing between adjacent defects. For values of Dy/(kdΔy)⪢1, we find that shortcuts increase the flux by significant amounts. The magnitude of the flux is increased as the imperfection spacing Δy is decreased. For small values of Δy, permeation through shortcuts becomes sorption-limited so that decreasing Δy further does not increase the flux through a single shortcut. However, as Δy is decreased, the concentration of shortcuts increases, thereby increasing the total contribution of the shortcuts to the flux. We have found regimes where increasing Δy or decreasing Dy decreases the overall flux, showing that permeation can be diffusion-limited by motion perpendicular to the transmembrane direction.

We have performed embedded-cluster calculations using density functional theory to investigate mechanisms of nitrogen substitution (nitridation) in HY and silicalite zeolites. We consider nitridation as replacing Si–O–Si and Si–OH–Al linkages with Si–NH–Si and Si–NH2–Al, respectively. We predict that nitridation is much less endothermic in HY (29 kJ/mol) than in silicalite (132 kJ/mol), indicating the possibility of higher nitridation yields in HY. To reveal mechanistic details, we have combined for the first time the nudged elastic band method of finding elusive transition states, with the ONIOM method of treating embedded quantum clusters. We predict that nitridation of silicalite proceeds via a planar intermediate involving a ring with pentavalent Si, whereas nitridation of HY is found to proceed via an intermediate similar to physisorbed ammonia. B3LYP/6-311G(d,p) calculations give an overall barrier for silicalite nitridation of 343 kJ/mol, while that in HY is 359 kJ/mol. Although the overall nitridation barriers are relatively high, requiring high temperatures for substitution, the overall barriers for the reverse processes are also high. As such, we predict that once these catalysts are made, they remain relatively stable.The mechanism of nitridation in HY and silicalite is revealed using denstiy functional theory. The barriers for forward and backward processes are large, indicating that nitrided zeolites are stable once formed.Download high-res image (57KB)Download full-size image

We have modeled the formation kinetics of nitrogen-substituted (nitrided) zeolites HY and silicalite; we have also modeled the stability of nitrided sites to heat and humidity. These kinetic calculations are based on mechanisms computed from DFT-computed pathways reported in our previous work. Reactant ammonia and product water concentrations were fixed at various levels to mimic continuous nitridation reactors. We have found that zeolite nitridation — replacing Si–O–Si and Si–OH–Al linkages with Si–NH–Si and Si–NH2–Al, respectively — proceeds only at high temperatures (>600°C for silicalite and >650°C for HY) due to the presence of large overall barriers. These threshold temperatures are in good agreement with experiments. Nitridation yields were found to be sensitive to water concentration, especially for silicalite where nitridation is more strongly endothermic. As a result, overall nitridation yields in silicalite are predicted to be much lower than those in HY. The stability of nitrided sites was investigated by modeling the kinetics of nitridation in reverse, going back to untreated zeolite plus ammonia. Using 10 h as a benchmark catalyst lifetime, nitrided silicalite and HY half-lives exceeded 10 h for temperatures below 275 and 500 °C, respectively, even at saturation water loadings. As such, our calculations suggest that nitrided silicalite and HY zeolites require high temperatures to form, but once formed, they remain relatively stable, auguring well for their use as shape-selective base catalysts.Nitrided HY and silicalite are predicted to remain stable below threshold temperatures (500 °C and 275 °C, respectively) even at saturation water loadings, warranting their use as shape-selective basic catalysts.Download high-res image (93KB)Download full-size image