Prediction of Solubilities and Crystal Shapes with Molecular Dynamics Simulations

Crystallization is a paramount purification technique in the pharmaceutical industries. Typically, pharmaceutical products are crystallized from solution as crystalline powders. The growth of these crystals in solution is influenced by the choice of their growth environment such as solvent, antisolvent, and supersaturation, and is acutely complicated to model. Despite the importance in industry, crystallization phenomena are still not exhaustively understood. Purely experimental approaches for the design of the drug substance’s crystallization process are often tedious trial-and-error endeavors and by no means the most resource efficient approaches. Molecular dynamics simulations have emerged as a practical tool to obtain insights into the complex crystallization phenomena on the atomistic level. Molecular dynamics can go, where experiments reach their limits. The simulations allow the study of the interplay between the crystal surface, solute, solvent, and antisolvent on the nano-scale and enable to unveil important mechanistic insights and predict thermodynamic as well as kinetic system properties. In this work, molecular dynamics simulation setups were developed and used for the prediction and analysis of solubility and crystal shapes, which are both important properties in the design of crystallization processes. The introduced simulation methodologies lead to improvements of the predictive capacities of molecular dynamics simulations used to study solubilities and crystal shapes. A simulation methodology for the prediction of solubilities of organic molecules is introduced. The simulation methodology is based upon the growth and dissolution of molecules at kink sites. The energy difference between the molecule crystallized at the kink site and dissolved in solution, allows to identify the solubility as the solute mole fraction, where this energy difference is zero. To efficiently sample the growth and dissolution of solute molecules at the kink site, well-tempered Metadynamics were used through a collective variable. The collective variable captures the slow degrees of freedom of the growth process, namely the solute's diffusion and adsorption at the kink site as well as the desolvation of the kink site and vice versa for the dissolution. The methodology is showcased for urea and naphthalene in a variety of solvents and reproduces well the experimental trends. The kink growth simulation setup is extended to predict the solubility of organic salts in complex growth environments. The simulation setup is demonstrated for anhydrous sodium acetate grown in a variety of solvent-antisolvent mixtures. The simulation results are compared to experimental results obtained within the research group, and the agreement between the two approaches is excellent. From the direct comparison between simulations and experiments, the advantages of combining molecular dynamics simulations with experiments for the study of complicated crystallization processes of organic salts is delineated in detail. Furthermore, a molecular dynamics simulation setup is presented for the prediction of crystal shapes of organic crystals. The simulation setup thereby focuses on the growth of crystal layers and allows to estimate the growth rates of crystal faces, which in turn enable the prediction of the steady state crystal shape. To bring the crystal growth within the timescale of molecular dynamics simulations, well-tempered Metadynamics are applied with a collective variable that is spatially constrained to the growing surface layer. This collective variable is designed such that it quantifies the most important slow degree of freedom, i.e. the crystallinity of the growing surface layer, and also allows efficient enhanced sampling. The simulation setup is demonstrated for naphthalene grown from ethanol solution. The simulations are in line with experiments and accurately predict the crystal shapes in dependence of supersaturation. The mechanistic details of crystal growth are delineated as well for all the relevant crystal faces of naphthalene. In summary, the presented results obtained with the reported molecular dynamics simulation techniques are capable of predicting important crystallization properties and can be utilized to gain a deeper understanding of the crystallization process on the nano-scale. The simulation setups, reported in this work, are general and can be implemented and further advanced to rationalize, guide, and speed up experimental crystallization process designs and screening campaigns