Quantum chemical calculations for predicting the partitioning of drug molecules in the environment

Abstract

Regional and temporal trends in legal and illicit drug use can be tracked through monitoring of municipal wastewater, ambient air, indoor air, and house dust. To assess the analytical result for the selected environmental matrix, reliable information on the partitioning of the target substance between the different compartments is required. The logarithmic partition coefficients octanol/water (log K_OW), octanol/air (log K_OA) and air/water (log K_AW) are usually applied for this purpose. Most drug molecules are semi-volatile compounds with complex molecular structures, the handling of which is subject to legal regulations. Chemically, they are often acids, bases, or zwitterions. Consequently, the physical and chemical properties are in most cases not determined experimentally but derived from quantitative structure–activity relationships (QSARs). However, the lack of experimental reference data raises questions about the accuracy of computed values. It therefore seemed appropriate and necessary to calculate partition coefficients using alternative methods and compare them with QSAR results. We selected 23 substances that were particularly prominent in European and US drug reports. Different quantum mechanical methods were used to calculate log K_OW, log K_OA, and log K_AW for the undissociated molecule as a function of temperature. Additionally, the logarithmic hexadecane/air partition coefficient log K_HdA ≡ L and the logarithmic vapor pressure of the subcooled liquid log P_L were determined in the temperature range 223 < T/K < 333. Despite the sometimes high variability of the parameters, it is possible to estimate how an investigated substance distributes between air, water and organic material.

---

Environmental significance

According to the United Nations, drug use is growing faster than the global population. Both the increasing demand and the diversity of supply are exacerbating a social problem that is certainly cause for serious concern. A proven method in forensic drug analysis is detection in wastewater, soil, house dust, etc. However, a realistic assessment of the distribution of drugs in environmental matrices requires reliable values for partition coefficients and vapor pressure. These are often unavailable because the measurements are complex or because of legal hurdles. Therefore, 23 representative substances were selected from the vast amount of available drugs, and the physical properties were calculated using quantum mechanical methods. The results were subjected to a critical analysis for their plausibility.

1 Introduction

Psychoactive substances are closely linked to human cultural history. As early as 3000 BC, hallucinogenic mushrooms of the genus Psilocybe spp., whose main ingredient is psilocybin, were consumed in Mesoamerica.¹ Also in South America, chewing coca leaves has been a popular practice for centuries to combat fatigue and hunger. The sleep-inducing and pain-relieving effects of opiates were already recognized in prehistoric times and became a frequent component of medicine from the 16th century onwards. The development of modern chemistry in the 19th century made it possible to produce amphetamines, morphine and heroin, which were often used in the military sector.² The first prohibitions on drug trafficking in the 19th century were motivated by trade policy. Health-related measures were introduced starting in the 1920s. Since the 1980s, the illicit drug market has been flooded with synthetic products known as so-called designer drugs. For example, Gerona³ lists 21 structurally different synthetic cannabinoids. Fentanyl was first synthesized in 1960, and a wide range of fentanyl analogues are now also available.⁴ For years, it has been evident that synthetic chemistry is progressing faster than the drug regulations by authorities. Nitrous oxide (N₂O), known as laughing gas, is not only used as an anaesthetic in medicine or in the kitchen to whip cream, but more recently also abused in some countries without restrictions as a so-called lifestyle drug.⁵

One method for tracking drug usage is the analysis of environmental matrices. Wastewater not only allows to identify which drugs are being consumed locally, but also to identify delivery routes and emerging or declining trends.^6–8 In the indoor environment, house dust is a suitable matrix for screening purposes.⁹ Of equal importance is the question of how a chemical substance distributes between different environmental matrices. Outdoors, the simplest model includes the compartments air, water, and organic phase.¹⁰ In reality, however, the dynamics are considerably more complex. For example, semi volatile airborne substances require consideration of the partitioning between the gas and the particle phase^11,12 or the aqueous and gas phase.⁸ Analogously, the same applies to indoor areas with gas phase, particle phase, house dust, clothing and other surfaces.^13–17

To assess the environmental equilibrium partitioning of a substance, parameters are required that describe molecular interactions between different phases. These are essentially the octanol/water partition coefficient (K_OW), the octanol/air partition coefficient (K_OA), the air/water partition coefficient (K_AW) and the saturation vapor pressure of the subcooled liquid (P_L).¹⁰ A frequently used quantity in linear-free-energy-relationships (LFER) is the hexadecane/air partition coefficient K_HdA and its logarithmic form log K_HdA ≡ L, respectively.¹⁰ However, these parameters are often difficult to measure or require significant experimental effort. Therefore, early attempts were made to predict physical and chemical properties from the molecular structure.¹⁸ With the increase in computing power, these models became more complex.^19,20 Advanced quantitative-structure–activity-relationships (QSAR) models use molecular descriptors, often in conjunction with machine learning.^21–23 Quantum mechanical (QM) methods provide a more fundamental approach to obtain these quantities of interest. These methods can, among other things, predict the solvation energy ΔG_solv in the respective solvents of interest.

Compared to other substances, the properties of drug molecules have not been investigated as thoroughly experimentally – often due to their complex structure and respective legal requirements. The use of theoretical methods is therefore common, which naturally raises questions about their accuracy. It is known that popular prediction tools such as EpiSuite and SPARC provide unreliable values for large molecules.²⁴ On the other hand, the aforementioned quantum chemical methods are often not suffering from such problems, but require advanced know how and computational effort.

For this study, we selected 23 molecules representing a broad spectrum of currently legally and illegally consumed drugs. This includes both the long-known and newly emerging substances. The partition coefficients log K_OW, log K_OA, log K_AW, log K_HdA, and P_L were calculated using quantum chemical methods. The obtained results were compared with data from popular QSAR/QSPR tools and databases. The comparison data had partly experimental and partly predictive background. Since some molecules are weak bases or acids, the pK_a values were also taken into account. Our results are not only of interest regarding the possibilities and limitations of using QM and QSAR computed parameters for estimating the partitioning of drug molecules in the environment. We additionally discuss the temperature dependence of these parameters in the interval 283–308 K based on the calculated temperature-dependent free energy of solvation ΔG_solv.

2 Compounds and methods

2.1 Drug molecules

For the selection of drug molecules, three main references were used: (i) the report of the German Bundeskriminalamt (BKA) (Federal Criminal Police Office) on drug-related crime in Germany 2023,²⁵ (ii) the World Factbook: Illicit Drugs of the U.S. Central Intelligence Service (CIA)²⁶ and (iii) the Commonly Used Drugs Charts of the U.S. Institute on Drug Abuse (NIDA).²⁷ From the multitude of possibilities, we selected 23 molecules. The criteria included their relevance for production and consumption, their future significance, their diversity of structural differences, their detectability in air, water, and bio-monitoring, as well as the required computational effort. All molecules investigated are listed in Table 1 with their abbreviation and CAS number; their respective structure is shown in Fig. 1. Boiling points and pK_a values were taken from the OPERA database.²⁸

Table 1 Investigated drugs, their abbreviations and CAS numbers. Except for N₂O, the boiling points and pK_a were taken from OPERA.²⁸ The boiling point of N₂O was taken from Rumble et al.,²⁹ n.a. = not available

Drug	Abbr.	CAS no.	M (g mol^−1)	BP (°C)	pKa
1-Benzylpiperazine	BZP	2759-28-6	176.26	289	9.80
4-Hydroxybutanoic acid	HBA	591-81-1	104.11	219	5.23
Amphetamine	AMP	300-62-9	135.21	203	9.87
Cathinone	KHAT	71031-15-7	149.19	261	6.39
Cocaine	COC	50-36-2	303.35	332	7.78
Midomafetamine (formerly ecstasy)	MDMA	42542-10-9	193.25	285	9.96
Fentanyl	FEN	437-38-7	336.47	353	7.15
Heroin	HER	561-27-3	369.42	398	8.28
Ketamine	KET	6740-88-1	237.74	312	8.08
Lysergic acid diethylamide	LSD	50-37-3	323.42	361	5.66
Mephedrone	MPD	1189805-46-6	177.24	274	6.63
Mescaline	MES	54-04-6	211.26	294	9.75
Methadone	MTD	76-99-3	309.45	353	6.76
Methamphetamine	METH	537-46-2	149.23	308	9.97
Methaqualone	MEQ	72-44-6	250.30	332	2.23
Methylphenidate	MEP	113-45-1	233.31	307	10.54
Morphine	MOR	57-27-2	285.34	385	8.06
Nitrous oxide	N2O	10024-97-2	44.01	−89	n.a.
Pethidine	PET	57-42-1	247.33	295	8.70
Phencyclidine (angel dust)	PCP	77-10-1	243.39	288	10.49
Psilocybin	PSY	520-52-5	284.25	286	3.63
Δ^9-Tetrahydrocannabinol	THC	1972-08-3	314.47	328	7.80
Trifluoromethylphenylpiperazine	TFPP	15532-75-9	230.23	289	6.91

---

Fig. 1 Structures of the investigated drugs with their abbreviations, see Table 1 for the chemical names.

2.2 Theory

We consider the partitioning of a substance between two adjacent phases α and β. The ratio of the equilibrium concentrations is known as the partition coefficient K_αβ = c_α/c_β. In the case of missing solute–solute correlations, K_αβ is directly related to the free energy of transfer ΔG_αβ of solute between the phases α and β via^30,31where R is the molar gas constant and T the temperature. According to eqn (1), calculation of ΔG_αβ allows for a direct determination of the partition coefficient K_αβ.³² If on the other hand, α denotes the pure ideal gas phase with concentration c_α = p_α/(RT) and β its pure condensed phase with concentration c_β = ρ_β/M, eqn (1) can directly be used to calculate the vapor pressure P_L ≡ p_α (analogous to ref. 33) viawhere M is the molar mass and ρ_β the density of phase β. If T < T_fus of the substance, P_L in eqn (2) refers to the vapor pressure of the subcooled liquid. For arbitrary phases β, the density of the phase ρ_β is not always known. For this reason, we set it to the density of liquid water of 997 kg m⁻³. This of course introduces an error, which, however, is lower than 0.1 kcal mol⁻¹ in terms of ΔG_αβ and thus negligible in terms of the accuracy of our computational workflow.^33,34 In addition, we use the geometries obtained in octanol for the ones in the subcooled liquid. Most drug-like compounds in their condensed (subcooled-liquid) state have a relative dielectric permittivity ε_r ≈ 8–15, close to that of 1-octanol.³⁵ We therefore optimize geometries in an octanol-like dielectric as a proxy for each compound's own liquid. We do not optimize each drug in its individual ε_r because experimental compound-specific values are generally unavailable; although ε_r can be computed in principle, there is no straightforward, robust protocol for diverse drug-like molecules, the required calculations are computationally demanding, and the setup is often substance-specific.^36–38

2.2.1 Ensemble-averaged free energies. To obtain the quantities of interest herein, the free energy of a substance in all phases (gas, water, octanol, and hexadecane) has to be known. As non-rigid molecules often have multiple relevant conformers contributing to the total free energy, we need to consider the ensemble of conformers that contribute to the total free energy of that substance. The following theory and equations are based on previous work and are explained in more detail there (e.g. ref. 34, 39 and 40). The free energy of a single conformer i is given by the gas-phase electronic energy E_el, the thermostatistical correction, and the possibly needed solvation free energy contribution ΔG_solv(T).

G_i(T) = E_el,i + G_trv,i(T) + ΔG_solv,i(T) (3)

The electronic energy can be obtained with any electronic structure method – and is obtained in our case with Density Functional Theory (DFT). The thermostatistical correction at temperature G_trv,i(T) accounts for translation, rotation, and vibration degrees of freedom and including the zero-point vibrational and volume work terms. To obtain the free energy of the substance from its individual conformers, the individual conformer free energies G_i need to be Boltzmann-weighted via eqn (4).

The free energy of a conformer ensemble should additionally include a conformational free energy part G_conf that stems from conformational entropy −TS_conf as a result of mixing multiple populated conformers.³⁴ This contribution is, however, most often neglected due to the huge amount of computation cost needed to accurately determine S_conf.⁴¹

Generally, the enthalpic contributions are very similar for gas- and solution-phases, similarly to the vibrational entropy contribution S^gas_v ≈ S^solv_v.⁴² The terms left are thus just the change in rotational and translational entropy. Because this change is generally already implicitly included in the solvation free energy of quantum chemical solvation models due to the parameterization on experimental data, we neglect the additional thermostatistical correction G_trv,i(T).

2.2.2 Accounting for solvation effects. To obtain free energy of a conformer (eqn (3)) in solution, the solvation contribution ΔG_solv (i.e., the solvation free energy) is needed. The solvation free energy describes the change in free energy when transferring a substance from the gas phase to a liquid.⁴³ In computational chemistry, it is often decomposed into several contributionswhere ΔG_ES are the electrostatic solvent–solute interactions (i.e., polarization), ΔG_non-ES denote non-electrostatic contributions (e.g., cavity formation), ΔG_N is the nuclear relaxation term describing the change in geometry when transferring a solute from the gas to liquid phase, and is the so-called standard-state correction.^44–46 This term is given bywhich depends on the ideal gas constant R, the standard concentration c^⊖ = 1 mol L⁻¹ and standard pressure p^⊖ = 1 atm and amounts to 1.89 kcal mol⁻¹ at 298.15 K.⁴⁷

The nuclear relaxation contribution, ΔG_N, captures the change in energy associated with the solute's structural reorganization upon transferring from the gas to the solution phase.⁴⁸ Concretely, it is computed as the difference between the gas phase electronic energy of the molecule at its solution phase optimized geometry, and the gas-phase electronic energy at its gas-phase optimized geometry. Because the gas phase geometry is a minimum on the gas phase potential energy surface, any other geometry – including the one favored in solvation – can only have equal or higher gas-phase energy. Consequently, ΔG_N will always be equal or larger than zero. Note that by performing full geometry optimizations both in the gas phase and in solution and using those respective electronic energies, the nuclear relaxation contribution is automatically included. Generally, it is very desirable to use the optimized geometries in all respective phases. In cases like, e.g., strong zwitter ions, the solute in solution will be zwitterionic, whereas in gas phase it will be neutral. Not using the respective optimized geometries will yield large errors for those systems. However, many implicit solvation models are parameterized on gas phase structures only, partially absorbing the nuclear relaxation into the parameterization and the model itself.⁴⁴ It is, however, not possible to correctly account for the nuclear relaxation of very complex solutes in this way, which can lead to significant errors.⁴⁸ In the following, we will shortly introduce the used implicit quantum-chemical solvation models and their workings.

2.2.2.1 SMD. The universal solvation model based on electron density⁴⁴ is an improvement over conductor-like polarizable continuum model (CPCM)^49,50 by adding a so-called cavity-dispersion solvent-structure term, that accounts for interactions that regular CPCM does not account for like the cavity creation and short-range dispersion interactions.
2.2.2.2 COSMO-RS and openCOSMO-RS. The conductor-like screening model for real solvents^51,52 is also based on CPCM, however, it uses the obtained surface charges of solute and solvent. These surface charges are then used to compute pairwise interaction energies and combinatorial contributions within an empirical statistical thermodynamic framework, from which it directly predicts the solvation free energy. openCOSMO-RS is the fully open-source implementation of COSMO-RS and both differ mainly in their parameterization.

Besides the three mentioned models, the aforementioned conductor-like polarizable continuum model^49,50 is tested and the results are shown in the SI. As CPCM only accounts for electrostatic interactions, it is not suitable to compute accurate solvation free energies or partition coefficients.

2.3 Computational workflow

To assemble the solvent-specific conformer–tautomer ensembles described in Section 2.2.1, we use the conformer–rotamer ensemble sampling tool, CREST, which uses meta-dynamics to gently bias molecules out of their local energy minima and into unexplored regions of the prototropic-tautomeric and conformational space.^53,54 Sampling with CREST using the GFN2-xTB tight-binding method, allowing an extensive exploration at a low cost. Because GFN2-xTB is a semi-empirical tight-binding method, however, a higher level optimization of the resulting ensemble is needed afterwards. For this, we use the command-line energetic sorting tool, CENSO, which carries out full DFT geometry optimizations and re-ranking.^34,55 This multi-level approach – broad, meta-dynamics-driven exploration followed by high-level DFT refinement – yields accurate ensembles at a reasonable cost. This workflow is carried out for each substance in each phase separately, yielding a refined ensemble for each phase, optimized at DFT (r²SCAN-3c) level. These structures are then used to carry out the calculation of solvation free energies. The workflow is shown in Fig. 2.

Fig. 2 Sketch of the employed CREST + CENSO workflow. The workflow is run separately for each drug and solvent combination.

2.4 Computational details

Quantum chemical calculations were performed with xTB 6.7.0 ⁵⁶ and ORCA 6.0.1.⁵⁷ The workflows used Open Babel 3.1.0,⁵⁸ CREST 3.0,^54,59 CENSO 2.0,³⁴ and MolBar 1.1.0.⁶⁰ The sampling and subsequent optimization was done using GFN2-xTB⁶¹ and the r²SCAN-3c composite method.⁶² Range-separated hybrid calculations are on the ωr²SCAN-D4/def2-TZVPPD level.^63–67 The matching def and def2 effective small core potentials (ECPs)^68,69 for heavy elements with Z > 36 were generally employed for all calculations. Matching general-purpose auxiliary basis sets are constructed on the fly using Stoychev et al.'s automatic generation of auxiliary basis sets in ORCA.⁷⁰ The RIJCOSX^71–73 approximation was used for all hybrid calculations. COSMO-RS^49,51,74 is calculated using TurboMole 7.9.0 ⁷⁵ with COSMOtherm C30-1601 and uses per default BP86/def-TZVP level of theory.^76,77 openCOSMO-RS⁵² results are obtained using ORCA and utilize BP86/def2-TZVPD. Solvation contributions of the CPCM^50,78 and SMD⁴⁴ models are obtained using the r²SCAN-3c composite method. The ALPB solvation model⁴⁸ was employed for all semi-empirical solution phase optimizations and the CPCM solvation model for all r²SCAN-3c-based optimizations.

2.5 Databases and literature data

The Open (Quantitative) Structure–activity/property Relationship App (OPEn(q)saRApp = OPERA)²⁸ provides QSAR/QSPR models for chemical properties and was primarily used for comparison with the quantum chemically calculated data. The workflows behind OPERA consist of complex algorithms,^79,80 whose structure is explained by accompanying documents directly in the database. The basis is an implemented dataset from which a training set and a test set are randomly generated. The desired property of the respective target substance is computed using the k-nearest neighbor (kNN) method with k = 5. This algorithm identifies the five substances from the dataset whose molecular descriptors are closest to the target substance and uses the Euclidean distance as metric.²⁸ If the target substance itself is among the nearest neighbors, the value is marked as “experimental”. This can be checked using the CompTox database⁸¹ (note that we used Version 2.5.3 – April 8, 2025), which also provides statistical parameters for the reliability of the value calculated by OPERA. Experimental data were particularly available for log K_OW, so we differentiated accordingly in this case (see SI). We also noticed that the log K_OW values often match the values of Hansch et al.⁸² In several cases, Hansch et al.⁸² refer to inaccessible sources, so the quality of the log K_OW data could not always be verified. For 13 substances we found experimentally determined Linear Solvation Energy Relationship (LSER) descriptors and logarithmic hexadecane/air partition coefficients in the UFZ-LSER database.⁸³ For the other 10 substances, the descriptors and log K_HdA(L) values were computed from their SMILES structures using a tool implemented in the UFZ-LSER database and described by Brown.⁸⁴ Sander's⁸⁵ database provides access to a compilation of Henry's law constants for organic and inorganic species in water from literature data. Additional calculated and measured literature data were considered and discussed accordingly. All comparison data extracted from the databases refer to 298 K. Experimental data for other temperatures, when available, are also listed in the SI.

3 Results and discussion

To systematically evaluate the results of the quantum mechanical calculations and the literature values, we based all comparisons on the data generated with COSMO-RS. Accordingly, Table 2 contains all partition coefficients and vapor pressures calculated with COSMO-RS for the 23 target molecules. The corresponding ΔG can be found in the SI, as well as the results of the openCOSMO-RS and SMD calculations and the associated costs for the computational workflow. The results obtained with the other QM methods are also listed, but are not discussed since they do not provide any further insights. With the exception of N₂O, log K_OA, log K_OW, log K_AW (calculated from the Henry solubility) and log P_L were obtained from OPERA. All log K_HdA (L) values were extracted from the UFZ-LSER database. Experimental vapor pressure data were also considered. The literature data are compiled with the respective sources in the SI. As OPERA does not provide data for N₂O, the log K_OW was estimated from the Ostwald coefficients for this substance. Using the data from Makranczy et al.⁸⁶ for 1-octanol and from Gabel and Schultz⁸⁷ for water, log K_OW = 0.55 is obtained at T = 298 K. This is in fair agreement with the values of log K_OW = 0.36 and 0.43 from Hansch et al.^82,88 The comparisons of COSMO-RS results with the other quantum mechanical methods and the literature data for partition coefficients at 298 K are shown in Fig. 3.

Table 2 Partition coefficients and vapor pressures of the subcooled liquid at 298 K for the drug molecules given in Table 1, quantum mechanically calculated using COSMO-RS

Drug	log KOW	log KOA	log KAW	log KHdA	log PL (Pa)
BZP	2.45	7.79	−5.34	6.64	−0.42
HBA	−0.47	7.38	−7.85	3.78	−0.73
AMP	2.15	6.24	−4.08	4.87	1.27
KHAT	4.96	6.93	−1.97	5.3	0.01
COC	1.06	9.37	−8.31	10.3	−3.25
MDMA	2.14	6.26	−4.12	6.22	0.54
FEN	4.47	13.41	−8.93	12	−6.79
HER	3.13	10.51	−7.38	9.36	−4.74
KET	1.63	9.67	−8.04	7.6	−3.56
LSD	3.43	13.91	−10.48	3.44	−7.82
MPD	2.47	7.28	−4.81	6.27	−0.38
MES	1.15	9.85	−8.69	7.83	−2.63
MTD	5.20	10.01	−4.81	9.59	−3.93
METH	2.94	5.73	−2.78	5.25	1.48
MEQ	2.29	10.25	−7.95	9.3	−4.26
MEP	3.94	8.88	−4.94	8.28	−2.16
MOR	2.10	12.02	−9.92	6.97	−5.87
N2O	0.40	0.06	0.34	−0.29	6.63
PET	3.75	8.72	−4.97	8.17	−2.10
PCP	5.19	8.37	−3.18	8.35	−1.96
PSY	−2.97	19.53	−22.5	13.58	−14.80
THC	6.97	14.52	−7.55	12.83	−7.31
TFPP	2.85	7.7	−4.85	6.81	−1.17

---

Fig. 3 Scatter diagrams of the partition coefficients at 298 K calculated using different methods, plotted against COSMO-RS. For the log K_HdA the UFZ-LSER data are used. The log K_AW according to OPERA were calculated from Henry solubilities using eqn (9). Notable points are marked with their respective abbreviation in the figure.

In contrast to our previous work on partition coefficients³² and vapor pressures,³³ the problem here is that there is hardly any reliable experimental data available for the 23 target compounds. Therefore, only computed values can be discussed to estimate their environmental behavior. Since most of the compounds are weak acids or bases, dissociation effects must also be considered.

3.1 Octanol/water partitioning

For the partition coefficients between 1-octanol and water, we generally find good agreement between all methods. openCOSMO-RS predicts slightly larger log K_OW values compared to COSMO-RS, whereas SMD agrees well with the COSMO-RS data. Compared to the QM results, the OPERA data are broader distributed, with a slight tendency to underestimate COSMO-RS. The methods disagree the most for PSY, MPD and KHAT. For PSY, SMD and OPERA predict log K_OW of 0.48 and 1.18, whereas the COSMO-RS variants predict −2.97 and −2.20. The negative values are plausible because the zwitterion is the thermodynamically preferred tautomer and its solubility in water is higher than in organic solvents.⁸⁹ For MPD and KHAT, OPERA predicts a log K_OW around 2–3 log units lower than the quantum chemical methods. Notable is also N₂O, where SMD predicts a log K_OW of around −1.85, whereas all other methods predict values between 0.36 and 0.83.

Assessing the root cause of SMD's misprediction is challenging, since the model is highly empirical and relies on parameters fitted to reproduce experimental data. SMD predicts a higher (i.e. less favourable) solvation free energy for N₂O in 1-octanol of around 3 kcal mol⁻¹, which leads to the wrong sign in the prediction of the log K_OW. Similar shortcomings of SMD have been noted elsewhere in the literature (e.g. in ref. 90–92). However, it should still be noted, that an error of a few kcal mol⁻¹ is still reasonable and often to be expected for implicit solvation models in more difficult cases.^92,93

As mentioned in a previous section, it is difficult to assess the quality of the OPERA data, because no information about the origin of the data is provided. It can only be distinguished between experimental and computed values (see Fig. 3). In this context, “experimental” means that the substance is present in the OPERA dataset. However OPERA provides information about the training set and test set, the confidence level, and the five nearest neighbors, which makes it easier to assess the validity of the calculated value. In the case of morphine (MOR), for example, it can be assumed that the log K_OW = 0.89 given by OPERA is the experimentally determined value published by Avdeef et al.⁹⁴ Moreover, the work of Avdeef et al.⁹⁴ touches the important aspect that most of the 23 target molecules are acids or bases, but the data in Table 2 refer to the undissociated molecule. The pH dependence of log K_OW is discussed in the next section. Furthermore, psilocybin can form a zwitterion through intramolecular proton transfer, and the solvation energy ΔG naturally depends strongly on the tautomeric form.⁸⁹ Here it is obvious that the OPERA algorithm fails, because the five nearest neighbors are molecules without zwitterionic character and OPERA probably cannot take the zwitterionic state into account at all.

3.2 Dependence of log K_OW on the pH value

When considering octanol/water partitioning, it is particularly important to note that at 298 K, a considerable amount of water dissolves in 1-octanol, with a molar fraction of 0.27 – while the solubility of 1-octanol in water is significantly lower.⁹⁵ Therefore, log K_OW usually refers to water-saturated 1-octanol. A further problem arises for ionizable substances.⁹⁶ Hansch and Leo⁹⁷ applied correction terms in their molecular fragment-based CLOGP method for calculating log K_OW to take into account the properties of ionizable compounds and zwitterions. Experimental methods often use buffer solutions, so that the log K_OW of the undissociated acid or base is obtained. To better understand the environmental behavior and bioavailability of ionizable compounds, the distribution coefficient D_OW was introduced, which is defined as the ratio of the neutral and charged species in the lipid and aqueous phase at a given pH. In the pH range where the molecule is predominantly non-ionized, log D_OW = log K_OW. Based on the Henderson–Hasselbalch equation for buffer solutions, Scherrer and Howard⁹⁸ developed relationships between log D_OW and log K_OW for acids (eqn (7)) and bases (eqn (8)) under the assumption that charged species are not present in the octanol phase.

However, due to the solubility of water in 1-octanol, eqn (7) and (8) can only be considered approximate. For an exact treatment, the partition coefficient of the charged species is also required.⁹⁶ Due to the acidic or basic properties of the investigated compounds (see Table 1), reliable pK_a values are needed for the neutral molecule and the corresponding ion to convert log K_OW to log D_OW for a specific pH value. For morphine, Abraham et al.⁹⁹ reported a log D_OW = 0.76 at a pH of 8.9. The value calculated using eqn (8) is 0.89 − 0.06 = 0.83. A commonly used reference value for environmental conditions is pH 7. The agreement between pK_a obtained from OPERA and data from other sources varies considerably. For many substances, the agreement is good^94,100–102 (see also SI), but deviations of approximately two orders of magnitude must be considered in some cases. Examples include tetrahydrocannabinol (7.80 compared to 10.60 (ref. 103)) and methadone (6.76 compared to 8.94 (ref. 104)). Differences are also observed in the zwitterion psilocybin. OPERA yields only one pK_a value of 3.63, but the acidic phosphate group accounts for two dissociation constants, which Richter et al.¹⁰⁵ calculate as 1.87 and 6.21. In addition, the dimethylamine group has basic character with a pK_a of 9.24.¹⁰⁵

3.3 Octanol/air partitioning

This partition coefficient can be used to describe the distribution of an airborne substance between the gas phase and the particle phase,¹⁰⁶ or the distribution between the gas phase and settled house dust.¹⁰⁷ At a particle concentration in the air of 20 μg m⁻³ molecules with log K_OA ≤ 9 are completely in the gas phase, and those with log K_OA ≥ 13 are completely in the particle phase. The curve is sigmoidal and shifts to a higher particle phase/gas phase ratio with higher particle concentrations. Especially in the range of the inflection point of the sigmoidal curve, small differences Δlog K_OA are sufficient to significantly change the particle phase/gas phase ratio.¹⁸ Therefore, it is important to know the log K_OA in this range as precisely as possible.

If log K_OA is taken as a predictor for the partitioning of semi-volatile organic compounds (SVOCs) to aerosols,¹⁰⁶ then with log K_OA ≥ 12 more than 95% of the molecules should be in the particle phase, even at low particle concentrations.¹⁸ However, it is obvious from Fig. 3 that for PSY, MOR, THC, and LSD, the QM and OPERA values vary considerably (note that all OPERA values are computed), where QM tends to significantly higher log K_OA values. The vapor pressure of the subcooled liquid P_L can be used as a support here. According to the theory by Junge and Pankow,¹⁰⁸ molecules with log P_L (Pa) ≤ 6 are almost completely in the particle phase at a particle concentration of 20 μg m⁻³. Thus, it can be assumed that the four substances mentioned above are also essentially attached to particles. Especially for PSY, the difference between the QM methods and OPERA is more than 10 log units. This could be attributed to the zwitterionic nature of psilocybin in 1-octanol. Unless the QSAR tool OPERA fully captures the true zwitterionic nature of PSY, the model will predict an erroneous partitioning between 1-octanol and air. This applies analogously to other binary systems. Also noteworthy is the large difference of almost 4 log units for THC between COSMO-RS (14.52) and OPERA (10.82). In a previous calculation, we obtained 13.76,³² while Askari et al. calculated a value of 12.26.¹⁶ We thus believe the OPERA computed value is too small.

Generally, however, we find good agreement between all QM methods. Note that Parnis and Metcalfe¹⁰⁹ calculated a log K_OA = 9.94 for cocaine using COSMOtherm, which is very close to our value of 9.37. A notable deviation can be seen for N₂O. The log K_OA of N₂O can be estimated from experimental data. For 1-octanol, Makranczy et al.⁸⁶ report an Ostwald coefficient of 2.139 ml ml⁻¹ at 298 K and p = 760 mmHg from which log K_OA = +0.33 is obtained. The calculated values with COSMO-RS and openCOSMO-RS agree well. However, SMD predicts a too small value of −2.28. The reason is the same as stated in Section 3.1 – SMD predicts a too positive solvation free energy in 1-octanol, thus yielding a too small log K_OA.

3.4 Air/water partitioning

The air/water partition coefficient is calculated from the Henry solubility H_S according to eqn (9), with H_S in mol (m⁻³ Pa⁻¹), R = 8.314 J (mol⁻¹ K⁻¹) and T in K.

The conversion unit factor from mol per (atm m³) to mol per (m³ Pa) is 9.87 × 10⁻⁶. In the case of air/water partitioning, similar prerequisites apply to the dissociation of molecules as in octanol/water. The conditions are even simpler, since completely undissociated molecules can be assumed in air. This allows the eqn (7) and (8) with D_AW,acid = αK_AW and D_AW,base = (1 − α)K_AW to be used analogously for the air/water system.

As with the other partition coefficients, hardly any experimental log K_AW values are available for the molecules. For N₂O, Sander⁸⁵ calculated a Henry solubility of H_S = 2.4 × 10⁻⁴ mol (m⁻³ Pa⁻¹) from the data of Weiss and Price.¹¹⁰ Using eqn (9) this results in log K_AW = 0.23, which is in good agreement both with the experimental value log K_AW = 0.22 of Gabel and Schultz⁸⁷ determined from the Ostwald coefficient and with the values calculated with COSMO-RS. Sander's⁸⁵ database includes Henry solubilities for some other of the target substances, but these are not measured values. For cocaine, Parnis and Metcalfe¹⁰⁹ calculated a log K_AW = −6.59 using COSMOtherm and estimated log D_AW values for pH 0 (−14.37), pH 4 (−10.38) and pH 7 (−7.44) with pK_a = 7.78.

We generally observe good agreement between the QM models without any noticeable outliers. All OPERA values are computed and exhibit a significantly larger scatter, particularly for KHAT, FEN, LSD, and MEP, with deviations of about three orders of magnitude compared to COSMO-RS. Furthermore, OPERA predicts a log K_AW for PSY that is ten orders of magnitude higher. Again, it is likely that OPERA does not account for the zwitterionic nature of PSY. This leads to a significantly higher PSY fraction in the gas phase compared to the aqueous phase. Most of the 23 substances in Table 1 show low water solubility. However, it should be noted that this study concerns the actual molecules. Free bases such as cocaine, methamphetamine, etc., are typically further processed into their water soluble hydrochlorides. In the illicit intravenous use of heroin, the free base is dissolved using a weak acid directly before injection.

3.5 Hexadecane/air partitioning

The hexadecane/air partition coefficient log K_HdA ≡ L is useful for characterizing the nonspecific intermolecular interactions of organic chemicals in other partitioning processes and is used in its logarithmic form as a descriptor in LFERs. Due to the hydrophobic properties of n-hexadecane, molecular dissociation processes hardly need to be taken into account. In principle, log K_HdA is a quantity that can be experimentally measured using gas chromatography on nonpolar capillary columns.¹¹¹ However, the data available for the drug molecules listed in Table 1 are limited. An experimental value of log K_HdA = 0.17 for N₂O can be deduced from measurements of the Ostwald coefficient by Makranczy et al.¹¹² In the UFZ-LSER database, values are published for 13 of the 23 compounds, and only for morphine there is a publication by Abraham et al.¹¹³ cited. The authors argue that the descriptors for morphine were generated from experimental values, but the origin of the data is not stated. For the other 12 substances, no literature source is given and it is therefore not clear whether their log K_HdA values are measured or computed. The log K_HdA (L) values of the other 10 substances were calculated from the SMILES structures as described above.

The COSMO-RS methods are in very good agreement, while SMD systematically underestimates the log K_HdA values. Most values from the UFZ-LSER database agree well with COSMO-RS; however, outliers are observed for heroin, morphine and methadone with higher reported values compared to COSMO-RS. Particularly striking, however, is the large deviation of approximately 10 log units for LSD (see Fig. 3). The observed deviation originates from the presence of a relevant tautomer that differs from the structure given by the SMILES string (shown in Fig. 4C). Our QM workflow identified multiple distinct gas-phase tautomers, three of which are shown in Fig. 4, with form A being the most stable. Form B lies about 5.8 kcal mol⁻¹ higher than A, and form C about 12.4 kcal mol⁻¹ higher. This tautomeric effect is effectively the nuclear relaxation contribution to the solvation free energy as defined in eqn (5), and is therefore fully captured by our QM workflow. QSPR methods, in contrast, do not account for tautomerization and the associated energetic changes, leading to substantial deviations in their computed values. This highlights the need for tautomer screening and how a robust QM workflow can yield qualitatively different, and more consistent results. Stenzel et al.¹¹¹ compared experimentally determined log K_HdA values for 387 environmentally relevant compounds with COSMOtherm data and found deviations of up to three orders of magnitude. Assuming that mainly computed values are compared for the drug molecules considered here, most deviations are within the expected range; only for heroin and LSD are there notable discrepancies of five and ten orders of magnitude, respectively.

Fig. 4 Relative gas-phase energies ΔE (kcal mol⁻¹) of three LSD tautomers identified by our QM workflow. Form A is the most stable and set to 0.0 kcal mol⁻¹, while forms B and C are 5.8 and 12.4 kcal mol⁻¹ higher in energy, respectively. Hydrogens involved in tautomerization are highlighted in magenta. The SMILES-based tautomer is C.

3.6 Vapor pressure of the subcooled liquid

Among the vapor pressures obtained from OPERA, the CompTox database indicates that the values for amphetamine (AMP), heroin (HER), and cocaine (COC), are experimental. Verifiable experimental vapor pressure data at 298 K are available from original publications for 12 of the drug molecules considered here. There is agreement between OPERA and the results of Lawrence et al.¹¹⁴ for heroin and cocaine. However, Lawrence et al. extrapolated the vapor pressure of HER at 298 K from Antoine data, which were determined in the range 51 °C to 66 °C. The vapor pressure of cocaine was determined in the range between 21 °C and 41 °C.

The OPERA value for amphetamine roughly matches the value of Thornton et al.¹¹⁵ Since a reasonable experimental data set is available and to avoid confusion, we did not distinguish between experiment and prediction in the OPERA values. Except for PSY, the correlation of COSMO-RS with these experimental values and the OPERA values is very good (see Fig. 5) and within the range of our previous study.³³ As expected, the results for COSMO-RS and openCOSMO-RS are almost identical. The highest deviations can be found for OPERA for MDMA and HER with about two orders of magnitude. In general, the quality of experimental data is difficult to assess. Meng et al.¹¹⁶ used an indirect gas chromatographic method to extrapolate the vapor pressures of drug molecules using di-n-butyl phthalate (DnBP) as a standard. The authors report log(P_L/Pa) = −2.04 for the vapor pressure of DnBP at 298 K, whereas Gobble et al.¹¹⁷ give for the same temperature log(P_L/Pa) = −2.35. The experimentally determined vapor pressures for HER,¹¹⁴ COC,¹¹⁴ AMP,¹¹⁵ FEN,¹¹⁸ METH,¹¹⁹ and PCP¹¹⁷ appear plausible. However, the value for LSD published by Okumus et al.¹²⁰ is not directly verifiable. Uncertainty also exists for THC. COSMO-RS here yields log(P_L/Pa) = −7.31 (see Table 2). A previous calculation without using explicit solution and gas phase optimized geometries (i.e., not including the nuclear relaxation) yielded −4.76. This agrees with the experimental value of Meng et al.¹¹⁶ Additionally, the log K_OA = 12.27 calculated by Askari et al.,¹⁶ also using COSMO-RS, is also noticeably smaller than our calculated log K_OA of 14.52 (see Table 2).

Fig. 5 Left: Scatter diagram of experimental and OPERA calculated vapor pressure data at 298 K plotted against COSMO-RS data. Right: Temperature-dependent vapor pressure curves of the 23 target molecules, calculated with COSMO-RS. Note that the critical temperature of nitrous oxide is 309.5 K.²⁹

In our work, all reported vapor-pressure values include the explicit nuclear relaxation term – i.e., used the optimized geometry of each solute in both the gas and subcooled liquid phases. For comparison, values obtained without nuclear relaxation are provided in the SI, and in the cases of AMP, FEN, METH and PCP these results agree better with experimental vapor pressures. It is important to bear in mind, however, that these calculated quantum-chemical vapor pressure estimate relies on additional approximations – most notably the use of the subcooled liquid state to represent pure-liquid behavior – and that implicit solvation models (like COSMO-RS and openCOSMO-RS) are typically far more extensively parameterized for regular solution phase properties than for vapor pressures.^52,121 Additionally, we approximate the continuum for the geometries of the subcooled liquid to be that one of 1-octanol, for the reasons stated in Section 2.2. Generally, we do not expect a noticeable change in the results when using geometries optimized in their respective subcooled liquid phase, as structural differences across dielectric constants in the typical range for drug-like compounds (ε_r ≈ 8–15) are small and the resulting effect on the calculated vapor pressures is expected to be minor compared to the intrinsic uncertainty of the solvation models themselves.^33,52 The assessment of available publications has shown that some of the experimental vapor pressure determinations of drug molecules are questionable or at least have significant uncertainties. On the other hand, the results of QM calculations also depend on the model used. However, it has been demonstrated here and in previous work that COSMO-RS calculations are in most cases in good agreement with reliable experimental data. We are aware that this data set is limited. However, based on our earlier study,³³ we assume that the uncertainty (standard deviation) of our vapor pressure calculations remains in the range of 0.5 log units.

3.7 Temperature dependence of the partition coefficients and the vapor pressure

In general, the temperature of the environment often differs from the temperature given in the tabulated data. Therefore knowledge of the temperature dependence of the partition coefficients K_αβ and the vapor pressure P_L is essential for many purposes. If both α and β are liquids, the temperature dependence of K_αβ is relatively small. Experimental findings around room-temperature indicate that as a rule of thumb^10,122–127 one can use

In the case of K_OW both an increase as well as a decrease is observed in experiments.¹²² Our findings reveal a similar behavior with the same order of magnitude for dlog K_OW/dT. In the case of METH, direct comparison can be made with the experimental findings of Brodin et al.¹²⁷ (note that the entry of Brodin's result in Table 1 of Sangster's compilation¹²² is wrong). Brodin et al. report log K_OW = 2.13 ± 0.67 at 298.15 K. Our findings of log K_OW = 2.94 ± 0.50 are in fair agreement. The temperature dependence was observed to be¹²⁷ dlog K_OW/dT ≈ 0.01 K⁻¹ whereas our calculations give ≈0.003 K⁻¹. COSMO-RS is parameterized using the Gibbs free energy and does not truly account for enthalpy or entropy directly. Although the model includes some statistical mechanics, it can be limited and can underestimates the quantitative aspect of the temperature dependence.^128–132 The situation is generally much more clear if one phase is a liquid (1-octanol, hexadecane, water) and the other phase is air. As the solubility of substances in a liquid generally decreases with increasing temperature,

should result for α = liquid and β = air. This is also observed in our calculations. From a quantitative perspective, the question arises as to how partition coefficients and vapor pressures determined for a temperature T₁ can be converted to a temperature T₂. This is usually done using the van't Hoff equation, which assumes a constant phase transition enthalpy ΔH_αβ in the temperature interval under consideration. LFER equations based on experimental data sets are often used to predict the required enthalpies.^10,133,134 In contrast, QM methods allow the direct calculation of the solvation energy ΔG_solv at the desired temperature.

The calculations were performed over a temperature interval from 223 K to 333 K for all parameters. The results are fully summarized in the SI. The temperature dependence of log P_L is shown in Fig. 5 for the interval from 283 K to 313 K. Apart from PSY, the largest range occurred for MOR with log P_L(283 K) = −6.26 and log P_L(313 K) = −5.51.

With reference to Fig. 6, we discuss the influence of the solvation energies ΔG_solv for 1-octanol and water on the partition parameters log K_OW, log K_OA and log K_AW using COC as an example in the temperature range between 283 K and 313 K. We chose COC because, among the 23 target compounds, the relatively largest changes were observed here when neglecting PSY. As expected, ΔG_solv increases with temperature, with the effect being more pronounced for water than for 1-octanol. This then leads to the expected decrease in log K_OA and increase in log K_AW. The log K_OW value also increases with increasing temperature. Cocaine is generally a good example to demonstrate the temperature effect on log K_OW. A log coefficient ratio of ≈4 can be observed, with log K_OW(283 K) = 0.38 and K_OW(313 K) = 1.58. The ratio octanol/water is therefore 2.3 at 283 K and 38.0 at 313 K. At higher log K_OW the absolute mass transfer of cocaine into the water phase with decreasing temperature is significantly lower.

Fig. 6 Left: Temperature dependence of the solvation energies ΔG_solv for cocaine in water and 1-octanol. Right: Temperature dependence of the cocaine partition coefficients log K_OW, log K_OA and log K_AW. The calculations were performed using COSMO-RS.

After evaluating the available data, we consider QM calculations to be the best method to obtain reliable temperature-dependent partition coefficients and vapor pressures for the 23 substances. Phase transition enthalpies are only available for some of the experimentally determined vapor pressures. The same applies to the Abraham descriptors for LFER calculations. However, we also note that under ambient conditions, the differences in temperature-dependent partition behavior are small in most cases. Furthermore, a comprehensive analysis must take into account that pK_a values are also temperature-dependent.¹³⁵

3.8 Prediction of the partitioning behavior in the environment

The key question is whether the available data on partition coefficients and vapor pressure allow reasonable predictions about the fate of drug molecules in the environment. Assuming that the partition coefficients depend only on the free solvation energy ΔG_solv, eqn (12) for the relationship between K_OW, K_OA and K_AW is obtained.

log K_OA = log K_OW − log K_AW (12)

Eqn (12) is a simplified approach, which only serves to estimate the preferred accumulation in the hypothetical air/water/octanol system. Since there is a dependence on three parameters, a graphical representation of Mackay et al.¹³⁶ is used. In Fig. 7, log K_OW is plotted against log K_AW, while log K_OA is shown as diagonals.

Fig. 7 Plot of log K_OW versus log K_AW with log K_OA as diagonals at 298 K for the 23 target molecules to estimate which substances preferentially accumulate in which medium (air, water, 1-octanol). The figure is based on space diagrams from a publication by Mackay et al.¹³⁶

For all substances, the partitioning depends on the dimensions of the system, i.e., on the available amounts of air, water, and 1-octanol. Fig. 7 therefore only allows a general statement about which of the three media a substance would prefer. It is clear that the gas N₂O will preferably accumulate in air. Psilocybin (PSY) will be preferentially present in the aqueous phase, which is evident from both the partition parameters and the pK_a value. For the other 21 substances no educated guess can be made. Except for N₂O, all substances have small log K_AW values. The tendency toward accumulation in 1-octanol therefore increases with increasing log K_OA and log K_OW. Wania¹³⁷ uses a similar diagram, but refers to environmental compartments instead of air, water, and 1-octanol. Accordingly, most of these 21 substances would accumulate in soils and sediments or distribute between the compartments as multimedia chemicals.

The Henderson–Hasselbalch eqn (7) and (8) state, that in aqueous media for pH < pK_a, the molecule is predominantly present in the neutral form. Assuming pH 7, this is the case for 14 substances. Five substances have pK_a values between 5 and 7, meaning that both the neutral and ionic species exist at pH 7. For MEQ and PSY with pK_a = 2.23 and pK_a = 3.63, the equilibrium is further shifted towards the ionic form in aqueous media. This means that for a more precise analysis of substances with pK_a < 7 the calculation of log D and the partition coefficients of the ionic species would be necessary, but this is beyond the scope of this work. It must also be noted that more complex models are required to reliably estimate the bioaccumulation and partitioning of drugs in multiple environmental compartments.¹³⁸

4 Conclusion

It is probably pointless to list the drugs available on the illegal market individually. The drug supply is evolving extremely rapidly, and new substances or derivatives of already known substances appear almost daily. The CIA World Factbook,²⁶ for example, defines five categories of illicit drugs: narcotics, stimulants, depressants (sedatives), hallucinogens, and cannabis, in which parent molecules and their derivatives are grouped together. Based on drug consumption statistics, we selected 23 substances that cover the above categories and are structurally different. We included nitrous oxide in the list because it is sold in cans via vending machines and its physical properties are also of interest from a medical perspective.

Since only few experimental data on the partition coefficients of the molecules considered here are available, it is difficult to evaluate the results of the quantum mechanical calculations. However, our previous studies have already shown that quantum mechanics tends to more reliable results than QSAR methods.^32,33 The agreement between log K_OW, log K_AW and log K_OA determined using QSAR (OPERA) and QM (COSMO-RS) varies. With the exception of the zwitterion psilocybin (PSY), which always shows large deviations, various substances are conspicuous (see Fig. 3). There may be several reasons for this; one possibility could be the respective training sets of the QSAR method. In the case of log K_HdA, heroin (HER), morphine (MOR) and LSD show significant discrepancies between LFER and QM. The consideration of tautomers in the calculation of the log K_HdA of LSD demonstrates one of the strengths of a QM-based approach.

For the vapor pressure, a comparison with experimental data is possible, but it is obvious from Fig. 5 that COSMO-RS underestimates the experimental values at low vapor pressures. The temperature dependence of partition coefficients and vapor pressure is important from a practical point of view. As shown for cocaine (COC) in Fig. 6, QM methods allow for direct calculations for a specific temperature. However, Grimme et al.³⁴ state that COSMO-RS can only consider the temperature dependence of the free solvation energy semiquantitatively. Alternatively, the van't Hoff equation and the Clausius–Clapeyron equation can be used if the respective enthalpy of phase transfer is known.

Many of the 23 molecules are acids or bases. This makes the log K_OW values particularly uncertain, since these are determined for the free molecule, but in the aqueous phase, including the aqueous portion of the octanol phase, the pK_a and pH dependent dissociation must be taken into account. Correction formulas for log D are at best approximate. log K_OA, log K_AW and log K_HdA are less influenced by the acid/base properties.

Even in light of all experimental and theoretical uncertainties, the results communicated and discussed in this work contribute to a significantly improved understanding of the physical properties of drug molecules. It is possible to plausibly estimate whether certain substances preferentially accumulate in air, water, or an organic medium (1-octanol). Furthermore, valuable information is obtained on their distribution behavior between the gas and particle phase in aerosols¹³⁹ and on accumulation in house dust.¹⁰⁷

Conflicts of interest

The authors have no conflicts of interest to declare.

Data availability

All data supporting this article have been included as part of the supplementary information (SI). Supplementary information is available. See DOI: https://doi.org/10.1039/d5em00524h.

Acknowledgements

L. W. would like to express his deepest gratitude to Prof. Stefan Grimme for his exceptional support and for granting access to the computational resources at the Mulliken Center for Theoretical Chemistry. L. W. also greatly acknowledges support of the Stiftung Stipendien-Fonds des Verbandes der Chemischen Industrie e.V. through its Kekulé Fellowship program and the granted access to the Marvin cluster hosted by the University of Bonn. T. S. thanks Manuela Lingnau (WK) for her support in preparing the graphic abstract. All authors thank Professor Stefan Grimme for fruitful discussions.

References

1. F. Carod-Artal, Hallucinogenic drugs in pre-Columbian Mesoamerican cultures, Neurologia, 2015, 30, 42–49.
2. W. Sneader, Drug Discovery – A History, John Wiley & Sons Ltd, Chichester, UK, 2005.
3. R. Gerona, Designer Drugs, Elsevier Inc., Amsterdam, 2024.
4. P. Armenian, K. T. Vo, J. Barr-Walker and K. L. Lynch, Fentanyl, fentanyl analogs and novel synthetic opioids: A comprehensive review, Neuropharmacology, 2018, 134, 121–132.
5. E. A. Glasper, Nitrous oxide deaths among children and young people is no laughing matter, Compr. Child Adolesc. Nurs., 2023, 46, 259–261.
6. C. Ort, et al., Spatial differences and temporal changes in illicit drug use in Europe quantified by wastewater analysis, Addiction, 2014, 109, 1338–1352.
7. B. Petrie, R. Barden and B. Kasprzyk-Hordern, A review on emerging contaminants in wastewaters and the environment: Current knowledge, understudied areas and recommendations for future monitoring, Water Res., 2015, 72, 3–27.
8. Y.-C. Lin, T.-C. Hsiao and A. Y.-C. Lin, Urban wastewater treatment plants as a potential source of ketamine and methamphetamine emissions to air, Water Res., 2020, 172, 115495.
9. C. de Dios-Perez, A. Ballesteros-Gomez, J. Ballesteros and S. Rubio, Supramolecular solvent extraction of drugs and pharmaceuticals in indoor dust from public buildings, Sci. Total Environ., 2025, 968, 178859.
10. R. P. Schwarzenbach, P. M. Gschwend and D. M. Imboden, Environmental Organic Chemistry, John Wiley & Sons, Hoboken, 2017.
11. N. Mastroianni, C. Postigo, M. Lopez de Alda, M. Viana, A. Rodriguez, A. Alastuey, X. Querol and D. Barcelo, Comprehensive monitoring of the occurrence of 22 drugs of abuse and transformation products in airborne particulate matter in the city of Barcelona, Sci. Total Environ., 2015, 532, 344–352.
12. M. Viana, X. Querol, A. Alastuey, C. Postigo, M. Lopez de Alda, D. Barcelo and B. Artinano, Drugs of abuse in airborne particulates in urban environments, Environ. Int., 2010, 36, 527–534.
13. G. Morrison, N. V. Shakila and K. Parker, Accumulation of gas-phase methamphetamine on clothing, toy fabrics, and skin oil, Indoor Air, 2015, 25, 405–414.
14. A. Cecinato, C. Balducci and M. Perilli, Illicit psychotropic substances in the air: The state-of-art, Sci. Total Environ., 2016, 539, 1–6.
15. A. Cecinato, P. Romagnoli, M. Perilli and C. Balducci, Psychotropic substances in house dusts: a preliminary assessment, Environ. Sci. Pollut. Res., 2017, 24, 21256–21261.
16. A. Askari, F. Wania and A. W. H. Chan, Modeling the fate and involuntary exposure to tetrahydrocannabinol emitted from indoor cannabis smoking, Environ. Sci.: Atmos., 2023, 3, 760–772.
17. T. Salthammer, The legalization of cannabis may result in increased indoor exposure to Δ⁹-tetrahydrocannabinol, J. Hazard. Mater., 2024, 464, 132949.
18. T. Salthammer, Assessment of methods for predicting physical and chemical properties of organic compounds, Indoor Environ., 2024, 1, 100031.
19. W. Sippl and D. Robaa, QSAR/QSPR, in Applied Chemoinformatics: Achievements and Future Opportunities, ed. T. Engel and J. Gasteiger, WILEY-VCH, 2018, pp. 9–52.
20. I. V. Tetko, A. Yan and J. Gasteiger, Prediction of Physicochemical Properties of Compounds, in Applied Chemoinformatics: Achievements and Future Opportunities, ed. T. Engel and J. Gasteiger, Wiley-VCH, Weinheim, 2018, pp. 53–81.
21. Q. Zang, K. Mansouri, A. J. Williams, R. S. Judson, D. G. Allen, W. M. Casey and N. C. Kleinstreuer, In silico prediction of physicochemical properties of environmental chemicals using molecular fingerprints and machine learning, J. Chem. Inf. Model., 2017, 57, 36–49.
22. K. Mansouri, N. F. Cariello, A. Korotcov, V. Tkachenko, C. M. Grulke, C. S. Sprankle, D. Allen, W. M. Casey, N. C. Kleinstreuer and A. J. Williams, Open-source QSAR models for pKa prediction using multiple machine learning approaches, J. Cheminf., 2019, 11, 60.
23. X. Kang, B. Hu, M. C. Perdana, Y. Zhao and Z. Chen, Extreme learning machine models for predicting the n-octanol/water partition coefficient (K_ow) data of organic compounds, J. Environ. Chem. Eng., 2022, 10, 108552.
24. P. Schossler, T. Schripp, T. Salthammer and M. Bahadir, Beyond phthalates: Gas phase concentrations and modeled gas/particle distribution of modern plasticizers, Sci. Total Environ., 2011, 409, 4031–4038.
25. Bundeskriminalamt, Rauschgiftkriminalität – Bundeslagebild 2023, Bundeskriminalamt, Wiesbaden, 2024.
26. U.S Central Intelligence Agency, The World Factbook – Illicit Drugs, US Central Intelligence Agency – Office of Public Affairs, Washington D.C., 2025.
27. National Institute on Drug Abuse, Commonly Used Drugs Charts, National Institutes of Health, Bethesda, MD, 2023.
28. K. Mansouri, C. M. Grulke, R. S. Judson and A. J. Williams, OPERA models for predicting physicochemical properties and environmental fate endpoints, J. Cheminf., 2018, 10, 10.
29. J. R. Rumble, T. J. Bruno and M. J. Doa, Handbook of Chemistry and Physics – 102st Edition, CRC Press, Boca Raton, 2021.
30. A. Ben-Naim, Solvation Thermodynamics, Springer Science – Business Media, New York, 1987.
31. M. F. Vitha and P. W. Carr, The chemical meaning of the standard free energy of transfer: Use of van der Waals’ equation of state to unravel the interplay between free volume, volume entropy, and the role of standard states, J. Phys. Chem. B, 2000, 104, 5343–5349.
32. T. Salthammer, S. Grimme, M. Stahn, U. Hohm and W.-U. Palm, Quantum chemical calculation and evaluation of partition coefficients for classical and emerging environmentally relevant organic compounds, Environ. Sci. Technol., 2022, 56, 379–391.
33. M. Stahn, S. Grimme, T. Salthammer, U. Hohm and W.-U. Palm, Quantum chemical calculation of the vapor pressure of volatile and semi volatile organic compounds, Environ. Sci.:Processes Impacts, 2022, 24, 2153–2166.
34. S. Grimme, F. Bohle, A. Hansen, P. Pracht, S. Spicher and M. Stahn, Efficient quantum chemical calculation of structure ensembles and free energies for nonrigid molecules, J. Phys. Chem. A, 2021, 125, 4039–4054.
35. A. A. Maryott, Circular of the Bureau of Standards No. 514: Table of Dielectric Constants of Pure Liquids, National Bureau of Standards, Gaithersburg, MD, 1951.
36. Z. Tang, C. Chang, F. Bao, L. Tian, H. Liu, M. Wang, C. Zhu and J. Xu, Feasibility of predicting static dielectric constants of polymer materials: A density functional theory method, Polymers, 2021, 13, 284.
37. M. Jorge and L. Lue, The dielectric constant: Reconciling simulation and experiment, J. Chem. Phys., 2019, 150, 084108.
38. R. M. Redheffer, R. C. Wildman and V. O'Gorman, The computation of dielectric constants, J. Appl. Phys., 1952, 23, 505–508.
39. S. Spicher and S. Grimme, Efficient computation of free energy contributions for association reactions of large molecules, J. Phys. Chem. Lett., 2020, 11, 6606–6611.
40. M. Bursch, J. Mewes, A. Hansen and S. Grimme, Best-practice DFT protocols for basic molecular computational chemistry, Angew. Chem., Int. Ed., 2022, 61, e202205735.
41. J. Gorges, S. Grimme, A. Hansen and P. Pracht, Towards understanding solvation effects on the conformational entropy of non-rigid molecules, Phys. Chem. Chem. Phys., 2022, 24, 12249–12259.
42. A. J. Garza, Solvation entropy made simple, J. Chem. Theory Comput., 2019, 15, 3204–3214.
43. IUPAC, Solvation Energy (ST07102), https://goldbook.iupac.org/terms/view/ST07102.
44. A. V. Marenich, C. J. Cramer and D. G. Truhlar, Universal solvation model based on solute electron density and on a continuum model of the solvent defined by the bulk dielectric constant and atomic surface tensions, J. Phys. Chem. B, 2009, 113, 6378–6396.
45. J. M. Herbert, Dielectric continuum methods for quantum chemistry, WIREs Comput. Mol. Sci., 2021, 11, e1519.
46. A. V. Marenich, C. J. Cramer and D. G. Truhlar, Perspective on foundations of solvation modeling: The electrostatic contribution to the free energy of solvation, J. Chem. Theory Comput., 2008, 4, 877–887.
47. V. Radtke, et al., A unified pH scale for all solvents: part I – intention and reasoning (IUPAC Technical Report), Pure Appl. Chem., 2021, 93, 1049–1060.
48. S. Ehlert, M. Stahn, S. Spicher and S. Grimme, Robust and efficient implicit solvation model for fast semiempirical methods, J. Chem. Theory Comput., 2021, 17, 4250–4261.
49. A. Klamt and G. Schüürmann, COSMO: a new approach to dielectric screening in solvents with explicit expressions for the screening energy and its gradient, J. Chem. Soc., Perkin Trans. 2, 1993, 799–805.
50. V. Barone and M. Cossi, Quantum calculation of molecular energies and energy gradients in solution by a conductor solvent model, J. Phys. Chem. A, 1998, 102, 1995–2001.
51. A. Klamt, Conductor-like screening model for real solvents: A new approach to the quantitative calculation of solvation phenomena, J. Phys. Chem., 1995, 99, 2224–2235.
52. T. Gerlach, S. Müller, A. G. De Castilla and I. Smirnova, An open source COSMO-RS implementation and parameterization supporting the efficient implementation of multiple segment descriptors, Fluid Phase Equilib., 2022, 560, 113472.
53. P. Pracht, C. A. Bauer and S. Grimme, Automated and efficient quantum chemical determination and energetic ranking of molecular protonation sites, J. Comput. Chem., 2017, 38, 2618–2631.
54. P. Pracht, S. Grimme, C. Bannwarth, F. Bohle, S. Ehlert, G. Feldmann, J. Gorges, M. Müller, T. Neudecker, C. Plett, S. Spicher, P. Steinbach, P. A. Wesołowski and F. Zeller, CREST—A program for the exploration of low-energy molecular chemical space, J. Chem. Phys., 2024, 160, 114110.
55. F. Bohle, J. Seibert and S. Grimme, Automated quantum chemistry-based calculation of optical rotation for large flexible molecules, J. Org. Chem., 2021, 86, 15522–15531.
56. A. Katbashev, M. Stahn, T. Rose, V. Alizadeh, M. Friede, C. Plett, P. Steinbach and S. Ehlert, Overview on building blocks and applications of efficient and robust extended tight binding, J. Phys. Chem. A, 2025, 129, 2667–2682.
57. F. Neese, Software Update: The ORCA Program System—Version 6.0, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2025, 15, e70019.
58. N. M. O'Boyle, M. Banck, C. A. James, C. Morley, T. Vandermeersch and G. R. Hutchison, Open Babel: An open chemical toolbox, J. Cheminf., 2011, 3, 33.
59. P. Pracht, F. Bohle and S. Grimme, Automated exploration of the low-energy chemical space with fast quantum chemical methods, Phys. Chem. Chem. Phys., 2020, 22, 7169–7192.
60. N. van Staalduinen and C. Bannwarth, MolBar: a molecular identifier for inorganic and organic molecules with full support of stereoisomerism, Digital Discovery, 2024, 3, 2298–2319.
61. C. Bannwarth, S. Ehlert and S. Grimme, GFN2-xTB - An accurate and broadly parametrized self-consistent tight-binding quantum chemical method with multipole electrostatics and density-dependent dispersion contributions, J. Chem. Theory Comput., 2019, 15, 1652–1671.
62. S. Grimme, A. Hansen, S. Ehlert and J.-M. Mewes, r²SCAN-3c: A “Swiss army knife” composite electronic-structure method, J. Chem. Phys., 2021, 154, 064103.
63. J. W. Furness, A. D. Kaplan, J. Ning, J. P. Perdew and J. Sun, Accurate and numerically efficient r²SCAN meta-generalized gradient approximation, J. Phys. Chem. Lett., 2020, 11, 8208–8215.
64. J. W. Furness, A. D. Kaplan, J. Ning, J. P. Perdew and J. Sun, Correction to: Accurate and numerically efficient r²SCAN meta-generalized gradient approximation, J. Phys. Chem. Lett., 2020, 2020, 9248.
65. L. Wittmann, H. Neugebauer, S. Grimme and M. Bursch, Dispersion-corrected r2SCAN based double-hybrid functionals, J. Chem. Phys., 2023, 159, 224103.
66. E. Caldeweyher, S. Ehlert, A. Hansen, H. Neugebauer, S. Spicher, C. Bannwarth and S. Grimme, A generally applicable atomic-charge dependent London dispersion correction, J. Chem. Phys., 2019, 150, 154122.
67. D. Rappoport and F. Furche, Property-optimized Gaussian basis sets for molecular response calculations, J. Chem. Phys., 2010, 133, 134105.
68. K. A. Peterson, D. Figgen, E. Goll, H. Stoll and M. Dolg, Systematically convergent basis sets with relativistic pseudopotentials. II. Small-core pseudopotentials and correlation consistent basis sets for the post-d group 16–18 elements, J. Chem. Phys., 2003, 119, 11113–11123.
69. D. Andrae, U. Häußermann, M. Dolg, H. Stoll and H. Preuß, Energy-adjusted ab initio pseudopotentials for the second and third row transition elements, Theor. Chim. Acta, 1990, 77, 123–141.
70. G. L. Stoychev, A. A. Auer and F. Neese, Automatic generation of auxiliary basis sets, J. Chem. Theory Comput., 2017, 13, 554–562.
71. F. Neese, F. Wennmohs, A. Hansen and U. Becker, Efficient, approximate and parallel Hartree–Fock and hybrid DFT calculations. A ‘chain-of-spheres’ algorithm for the Hartree–Fock exchange, Chem. Phys., 2009, 356, 98–109.
72. R. Izsák and F. Neese, An overlap fitted chain of spheres exchange method, J. Chem. Phys., 2011, 135, 144105.
73. B. Helmich-Paris, B. De Souza, F. Neese and R. Izsák, An improved chain of spheres for exchange algorithm, J. Chem. Phys., 2021, 155, 104109.
74. A. Klamt, The COSMO and COSMO-RS solvation models, WIREs Comput. Mol. Sci., 2018, 8, e1338.
75. S. G. Balasubramani, et al., TURBOMOLE: Modular program suite for ab initio quantum-chemical and condensed-matter simulations, J. Chem. Phys., 2020, 152, 184107.
76. J. P. Perdew and W. Yue, Accurate and simple density functional for the electronic exchange energy: Generalized gradient approximation, Phys. Rev. B: Condens. Matter Mater. Phys., 1986, 33, 8800–8802.
77. F. Weigend and R. Ahlrichs, Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
78. M. Garcia-Ratés and F. Neese, Effect of the solute cavity on the solvation energy and its derivatives within the framework of the Gaussian charge scheme, J. Comput. Chem., 2020, 41, 922–939.
79. K. Mansouri, C. M. Grulke, A. M. Richard, R. S. Judson and A. J. Williams, An automated curation procedure for addressing chemical errors and inconsistencies in public datasets used in QSAR modelling, SAR QSAR Environ. Res., 2016, 27, 911–937.
80. K. Mansouri, J. T. Moreira-Filho, C. N. Lowe, N. Charest, T. Martin, V. Tkachenko, R. Judson, M. Conway, N. C. Kleinstreuer and A. J. Williams, Free and open-source QSAR-ready workflow for automated standardization of chemical structures in support of QSAR modeling, J. Cheminf., 2024, 16, 19.
81. A. J. Williams, C. M. Grulke, J. Edwards, A. D. McEachran, K. Mansouri, N. C. Baker, G. Patlewicz, I. Shah, J. F. Wambaugh, R. S. Judson and A. M. Richard, The CompTox Chemistry Dashboard: a community data resource for environmental chemistry, J. Cheminf., 2017, 9, 61.
82. C. Hansch, A. Leo and D. Hoeckman, Exploring QSAR – Hydrophobic, Electronic and Steric Constants, American Chemical Society, Washington D.C., 1995.
83. N. Ulrich, S. Endo, T. Brown, N. Watanabe, G. Bronner, M. H. Abraham and K.-U. Goss, UFZ-LSER Database V 3.2.1, Helmholtz Centre for Environmental Research, Leipzig, 2017.
84. T. Brown, Predicting hexadecane–air equilibrium partition coefficients (L) using a group contribution approach constructed from high quality data, SAR QSAR Environ. Res., 2014, 25, 51–71.
85. R. Sander, Compilation of Henry's law constants (version 5.0.0) for water as solvent, Atmos. Chem. Phys., 2023, 23, 10901–12440.
86. J. Makranczy, L. Rusz and K. Balog-Megyery, Solubility of gases in normal alcohols, Hung. J. Ind. Chem., 1979, 7, 41–46.
87. R. A. Gabel and B. Schultz, Solubility of nitrous oxide in water, 20-80 C, Anesthesiology, 1973, 38, 75–81.
88. C. Hansch, A. Vittoria, C. Silipo and P. Y. C. Jow, Partition coefficients and the structure-activity relationship of the anesthetic gases, J. Med. Chem., 1975, 18, 546–548.
89. L. Paul, C. T. Namba-Nzanguim, A. Telesphory, J. Oppong Mensah, D. Mteremko, R. Costa, S. M. Katundu, L. P. Kwiyukwa, N. D. Kambaine, J. Juvenary, S. Mlowe, G. Deogratias, D. M. Shadrack and A. S. Paluch, Investigation of the structure, stability, and relative solubility of psilocybin in water and pure organic solvents: A molecular simulation study, J. Mol. Liq., 2023, 392, 123479.
90. M. Amezcua, L. El Khoury and D. Mobley, SAMPL7 Host-Guest Challenge Overview: Assessing the reliability of polarizable and non-polarizable methods for binding free energy calculations (Version 2), ChemRxiv, 2020, preprint,  DOI:10.26434/chemrxiv.12768353.v2.
91. C. C. Zanith and J. R. Pliego, Performance of the SMD and SM8 models for predicting solvation free energy of neutral solutes in methanol, dimethyl sulfoxide and acetonitrile, J. Comput.-Aided Mol. Des., 2015, 29, 217–224.
92. R. Sure and S. Grimme, Comprehensive benchmark of association (free) energies of realistic host-guest complexes, J. Chem. Theory Comput., 2015, 11, 3785–3801.
93. A. Savoy, E. Paenurk, R. Pollice, P. H. Hünenberger and P. Chen, Solvation free energies of ion dissociations in dichloromethane: En route to accurate computations, J. Phys. Chem. B, 2025, 129, 6276–6288.
94. A. Avdeef, D. A. Barrett, P. N. Shaw, R. D. Knaggs and S. S. Davis, Octanol–, chloroform–, and propylene glycol dipelargonat–water partitioning of morphine–6–glucuronide and other related opiates, J. Med. Chem., 1996, 39, 4377–4381.
95. J. Sangster, Octanol-water Partition Coefficients: Fundamentals and Physical Chemistry, John Wiley & Sons Ltd, Chichester, 1997.
96. M. Kah and C. D. Brown, Log D: Lipophilicity for ionisable compounds, Chemosphere, 2008, 72, 1401–1408.
97. C. Hansch and A. Leo, Exploring QSAR – Fundamentals and Applications in Chemistry and Biology, American Chemical Society, Washington D.C., 1995.
98. R. A. Scherrer and S. M. Howard, Use of distribution coefficients in quantitative structure-activity relations, J. Med. Chem., 1977, 20, 53–58.
99. M. H. Abraham, K. Takacs-Novak and R. C. Mitchell, On the partition of ampholytes: application to blood–brain distribution, J. Pharm. Sci., 1997, 86, 310–315.
100. E. B. Leffler, H. M. Spencer and A. Burger, Dissociation constants of adrenergic amines, J. Am. Chem. Soc., 1951, 73, 2611–2613.
101. M. Rosa Boleda, M. Huerta-Fontela, F. Ventura and M. T. Galceran, Evaluation of the presence of drugs of abuse in tap waters, Chemosphere, 2011, 84, 1601–1607.
102. L. Settimo, K. Bellman and R. M. A. Knegtel, Comparison of the accuracy of experimental and predicted pK_a values of basic and acidic compounds, Pharm. Res., 2014, 31, 1082–1095.
103. E. R. Garrett and C. Hunt, Physicochemical properties, solubility, and protein binding of Δ9 -tetrahydrocannabinol, J. Pharm. Sci., 1974, 63, 1056–1064.
104. A. H. Beckett, Analgesics and their antagonists: Some steric and chemical considerations, J. Pharm. Pharmacol., 1956, 8, 848–859.
105. L. M. Richter, J. Al-Gousous, G. L. B. de Araujo, N. M. Davies and R. Löbenberg, Assessing the utility of in silico tools in early drug development: The case of a pharmaceutically relevant formulation of the prodrug psilocybin, J. Drug Delivery Sci. Technol., 2024, 92, 105305.
106. A. Finizio, D. Mackay, T. Bidleman and T. Harner, Octanol-air partition coefficient as a predictor of partitioning of semi-volatile organic chemicals to aerosols, Atmos. Environ., 1997, 31, 2289–2296.
107. C. J. Weschler and W. W. Nazaroff, SVOC partitioning between the gas phase and settled dust indoors, Atmos. Environ., 2010, 44, 3609–3620.
108. J. F. Pankow, Review and comparative analysis of the theories on partitioning between the gas and aerosol particulate phases in the atmosphere, Atmos. Environ., 1987, 21, 2275–2283.
109. J. M. Parnis and C. D. Metcalfe, The potential for atmospheric transport of contaminants to remote lakes: An assessment using COSMO-RS solvation theory, Chemosphere, 2025, 385, 144555.
110. R. Weiss and B. Price, Nitrous oxide solubility in water and seawater, Mar. Chem., 1980, 8, 347–359.
111. A. Stenzel, S. Endo and K.-U. Goss, Measurements and predictions of hexadecane/air partition coefficients for 387 environmentally relevant compounds, J. Chromatogr. A, 2012, 1220, 132–142.
112. J. Makranczy, K. Megyery-Balog, L. Rusz and L. Patyi, Solubility of gases in normal-alkanes, Hung. J. Ind. Chem., 1976, 4, 269–280.
113. M. H. Abraham, R. Sanchez-Moreno, J. Gil-Lostes, W. E. Acree, J. Enrique Cometto-Muniz and W. S. Cain, The biological and toxicological activity of gases and vapors, Toxicol. in Vitro, 2010, 24, 357–362.
114. A. H. Lawrence, L. Elias and M. Authier-Martin, Determination of amphetamine, cocaine, and heroin vapour pressures using a dynamic gas blending system and gas chromatographic analysis, Can. J. Chem., 1984, 62, 1886–1888.
115. M. Thornton, J. Chickos, I. V. Garist, M. A. Varfolomeev, A. A. Svetlov and S. P. Verevkin, The vaporization enthalpy and vapor pressure of (d)-amphetamine and of several primary amines used as standards at T/K = 298 as evaluated by correlation gas chromatography and transpiration, J. Chem. Eng. Data, 2013, 58, 2018–2027.
116. Y. Meng, A. H. Lichtmann, D. T. Bridgen and B. R. Martin, Inhalation studies with drugs of abuse, NIDA Res. Monogr., 1997, 173, 201–224.
117. C. Gobble, B. Walker and J. S. Chickos, The vaporization enthalpy and vapor pressure of fenpropidin and phencyclidine (PCP) at T/K = 298.15 by correlation gas chromatography, J. Chem. Eng. Data, 2016, 61, 896–902.
118. P. K. Gupta, K. Ganesan, P. K. Gutch, L. Manral and D. K. Dubey, Vapor pressure and enthalpy of vaporization of fentanyl, J. Chem. Eng. Data, 2008, 53, 841–845.
119. M. Thornton, C. Gobble and J. Chickos, The vaporization enthalpy and vapor pressure of S(+)-methamphetamine at T = 298.15 K by correlation gas chromatography, J. Chem. Thermodyn., 2014, 73, 51–56.
120. B. Okumus, A. Metin and I. A. Kariper, Understanding the effects of lysergic acid diethylamide and the importance of its prevention, Health Sci. Rev., 2023, 8, 100107.
121. A. Klamt, V. Jonas, T. Bürger and J. C. W. Lohrenz, Refinement and parametrization of COSMO-RS, J. Phys. Chem. A, 1998, 102, 5074–5085.
122. J. Sangster, Octanol-water partition coefficients of simple organic compounds, J. Phys. Chem. Ref. Data, 1989, 18, 1111–1229.
123. N. P. Bahadur, W.-Y. Shiu, D. G. B. Boocock and D. Mackay, Temperature dependence of octanol water partition coefficient for selected chlorobenzenes, J. Chem. Eng. Data, 1997, 42, 685–688.
124. A. Leo, C. Hansch and D. Elkins, Partition coefficients and their uses, Chem. Rev., 1971, 71, 525–616.
125. Y. D. Lei, S. Baskaran and F. Wania, Measuring the octan-1-ol air partition coefficient of volatile organic chemicals with the variable phase ratio headspace technique, J. Chem. Eng. Data, 2019, 64, 4793–4800.
126. S. Sarraute, H. Delepine, M. F. Costa Gomez and V. Majer, Aqueous solubility, Henry's law constants and air/water partition coefficients of n-octane and two halogenated octanes, Chemosphere, 2004, 57, 1543–1551.
127. A. Brodin, B. Sandin and B. Failerson, Rates of transfer of organic protolytic solutes between an aqueous and an organic phase V. The thermodynamics of mass transfer, Acta Pharm. Suec., 1976, 13, 331–352.
128. J. M. Parnis, D. Mackay and T. Harner, Temperature dependence of Henry’s law constants and K_OA for simple and heteroatom-substituted PAHs by COSMO-RS, Atmos. Environ., 2015, 110, 27–35.
129. Y. S. Sistla and A. Khanna, Validation and prediction of the temperature-dependent Henry’s constant for CO₂–ionic liquid systems using the Conductor-like Screening Model for Realistic Solvation (COSMO-RS), J. Chem. Eng. Data, 2011, 56, 4045–4060.
130. A. Klamt, Prediction of the mutual solubilities of hydrocarbons and water with COSMO-RS, Fluid Phase Equilib., 2003, 206, 223–235.
131. M. Diedenhofen and A. Klamt, COSMO-RS as a tool for property prediction of IL mixtures - A review, Fluid Phase Equilib., 2010, 294, 31–38.
132. D. Zaitsau, M. Trawny, J. Arndt, L. Wittmann, C. Held and R. Siewert, Critical evaluation of amino acids enthalpies of formation and phase transitions with a new approach based on solution calorimetry, quantum cluster growth, COSMO-RS and PC-SAFT, J. Chem. Eng. Data, 2025 , submitted.
133. S. Baskaran and F. Wania, Applications of the octanol–air partitioning ratio: a critical review, Environ. Sci.: Atmos., 2023, 3, 1045–1065.
134. S. Baskaran, A. Podagatlapalli, A. Sangion and F. Wania, Predicting the temperature dependence of the octanol–air partition ratio: a new model for estimating ΔU^o_OA, J. Solution Chem., 2023, 52, 51–69.
135. R. L. Thurlkill, D. A. Cross, J. M. Scholtz and C. N. Pace, pKa of fentanyl varies with temperature: implications for acid-base management during extremes of body temperature, J. Cardiothorac. Vasc. Anesth., 2005, 19, 759–762.
136. D. Mackay, A. K. Celsie and J. M. Parnis, The evolution and future of environmental partition coefficients, Environ. Rev., 2016, 24, 101–113.
137. F. Wania, Assessing the potential of persistent organic chemicals for long-range transport and accumulation in polar regions, Environ. Sci. Technol., 2003, 37, 1344–1351.
138. J. M. Parnis and D. Mackay, Multimedia Environmental Models – the Fugacity Approach, CRC Press, Boca Raton, FL, 2021.
139. B. J. Finlayson-Pitts and J. N. Pitts, Chemistry of the Upper and Lower Atmosphere, Academic Press, San Diego, 2000.

---