Comparison of selected polarizable and nonpolarizable water models in molecular dynamics simulations of ice I(h)

We present a molecular dynamics simulation study in which we determined the melting point of ice I(h) for the polarizable SWM4-NDP water model (Lamoureux et al., Chem. Phys. Lett., 2006, 418, 245-249) and compared the performance of several popular water force fields, both polarizable and nonpolarizable, in terms of melting temperature, stability and orientational structuring of ice. The simulations yield the melting temperature of SWM4-NDP ice as low as T(m) = 185 ± 10 K, despite the quadrupole moment of a SWM4-NDP water molecule being close to the experimental gas phase value. The results thus show that the dependence of T(m) on the molecular quadrupole, observed for the three- and four-site water models, is generally lost if polarization is explicitly included. The study also shows that adding polarizability to a planar three-charge water model increases orientational disorder in hexagonal ice. In addition, analysis of the tetrahedral order in bulk ice reveals a correlation between the pre-existing degree of orientational disorder in ice simulated using different polarizable and nonpolarizable models and the melting temperature of the models. Our findings thus suggest some new considerations regarding the role of polarization forces in a crystalline solid that may guide future development of reliable polarizable water models for ice.     

Efficient sampling of ice structures by electrostatic switching

An electrostatic switching procedure is introduced that enables the systematic generation of high-quality ice configurations at various temperatures. Proton disordered ice Ih configurations were generated for the TIP4P water model at temperatures from 50 to 240 K, for the SPC/Fw water model from 100 to 240 K, and for the DC97 water model at 240 K. The resulting configurations were found to properly sample the canonical ensemble. The dielectric constant of ice Ih was determined from the net dipole fluctuation of the ice configurations. The calculated dielectric constant compares favorably with the study by Rick and Haymet [J. Chem. Phys. 2003, 118, 9291]. However, our method gives smaller error bars, especially at lower temperatures. At temperatures above 200 K, a type of hydrogen-bond defect is identified that cannot be categorized as a D or L type defect.     

Homogeneous nucleation of water between 200 and 240 K: new wave tube data and estimation of the Tolman length

We have measured homogeneous nucleation rates of water at 200-240 K in the carrier gas helium, in the range of 10(13) - 10(17) m(-3) s(-1) using an expansion wave tube. The rates agree well with the results of Wolk and Strey [J. Phys. Chem. B 105, 11683 (2001)] in the range of overlap (220-240 K), and are summarized by the empirical fit J = S exp[4.6 + 0.244T-(906.8 - 2.914T)(ln S)(2)], with J the nucleation rate in m(-3) s(-1), S the supersaturation, and T the temperature in K. We find that the supersaturation dependence of both our rates and those of Wolk and Strey is lower than classical theory predicts, and that the critical cluster is smaller than the classical critical size. These deviations are explained in the framework of the Tolman theory for surface tension, and the "Tolman length" is estimated from our experimental results. We find a positive Tolman length that increases with decreasing temperature, from about 0.1 Angstrom at 260 K to (0.6 +/- 0.4) Angstroms at 200 K. We present a nucleation rate expression that takes the Tolman length into account and show that both the supersaturation and temperature dependence are improved, compared to the classical theory.     

Morphology and crystallization kinetics of compact (HGW) and porous (ASW) amorphous water ice

An investigation of porosity and isothermal crystallization kinetics of amorphous ice produced either by background water vapour deposition (ASW) or by hyperquenching of liquid droplets (HGW) is presented. These two types of ice are relevant for astronomical ice research (Gálvez et al., Astrophys. J., 2010, 724, 539) and are studied here for the first time under comparable experimental conditions. From CH(4) isothermal adsorption experiments at 40 K, surface areas of 280 ± 30 m(2) g(-1) for the ASW deposits and of 40 ± 12 m(2) g(-1) for comparable HGW samples were obtained. The crystallization kinetics was studied at 150 K by following the evolution of the band shape of the OD stretching vibration in HDO doped ASW and HGW samples generated at 14 K, 40 K and 90 K. Comparable rate constants of ～7 × 10(-4) s(-1) were obtained in all cases. However a significant difference was found between the n Avrami parameter of the samples generated at 14 K (n～ 1) and that of the rest (n > 2). This result hints at the possible existence of a structurally different form of amorphous ice for very low generation temperatures, already suggested in previous literature works.     

Quantum effects in liquid water and ice: model dependence

This paper explores the influence of choice of potential model on the quantum effects observed in liquid water and ice. This study utilizes standard rigid models and a more formal context for the rigid-body centroid molecular dynamics methodology used to perform the quantum simulations is provided. Quantum and classical molecular dynamics simulations are carried out for liquid water and ice Ih at 298 and 220 K, respectively, with the simple point charge/extended and TIP4P-Ew water models. The results obtained for equilibrium and dynamical properties are compared with those recently reported on TIP4P [L. Hernandez de la Pena and P. G. Kusalik, J. Chem. Phys. 121, 5992 (2004); L. Hernandez de la Pena et al., J. Chem. Phys 123, 144506 (2005)]. For the liquid, an energy shift of about 8% and an average molecular uncertainty of about 11 degrees were found independently of the water model. The self-diffusion coefficient consistently increases by more than 50% when going from the classical to the quantum system and quantum dynamics are found to reproduce the experimental isotopic shifts with the models examined. The ice results compare remarkably well with those previously reported for the TIP4P water model; they confirm that quantum effects are considerable and that the quantum mechanical uncertainty and the energy shifts due to quantization are smaller in ice than in liquid water. The relevance of these findings in the context of the construction of water models is briefly discussed.     

Homogeneous ice nucleation from aqueous inorganic/organic particles representative of biomass burning: water activity, freezing temperatures, nucleation rates

Homogeneous ice nucleation plays an important role in the formation of cirrus clouds with subsequent effects on the global radiative budget. Here we report on homogeneous ice nucleation temperatures and corresponding nucleation rate coefficients of aqueous droplets serving as surrogates of biomass burning aerosol. Micrometer-sized (NH(4))(2)SO(4)/levoglucosan droplets with mass ratios of 10:1, 1:1, 1:5, and 1:10 and aqueous multicomponent organic droplets with and without (NH(4))(2)SO(4) under typical tropospheric temperatures and relative humidities are investigated experimentally using a droplet conditioning and ice nucleation apparatus coupled to an optical microscope with image analysis. Homogeneous freezing was determined as a function of temperature and water activity, a(w), which was set at droplet preparation conditions. The ice nucleation data indicate that minor addition of (NH(4))(2)SO(4) to the aqueous organic droplets renders the temperature dependency of water activity negligible in contrast to the case of aqueous organic solution droplets. The mean homogeneous ice nucleation rate coefficient derived from 8 different aqueous droplet compositions with average diameters of ～60 μm for temperatures as low as 195 K and a(w) of 0.82-1 is 2.18 × 10(6) cm(-3) s(-1). The experimentally derived freezing temperatures and homogeneous ice nucleation rate coefficients are in agreement with predictions of the water activity-based homogeneous ice nucleation theory when taking predictive uncertainties into account. However, the presented ice nucleation data indicate that the water activity-based homogeneous ice nucleation theory overpredicts the freezing temperatures by up to 3 K and corresponding ice nucleation rate coefficients by up to ～2 orders of magnitude. A shift of 0.01 in a(w), which is well within the uncertainty of typical field and laboratory relative humidity measurements, brings experimental and predicted freezing temperatures and homogeneous ice nucleation rate coefficients into agreement. The experimentally derived ice nucleation data are applied to constrain the water activity-based homogeneous ice nucleation theory to smaller than ±1 order of magnitude compared to the predictive uncertainty of larger than ±6 orders of magnitude. The atmospheric implications of these findings are discussed.     

Molecular dynamics simulations of ice nucleation by electric fields

Molecular dynamics simulations are used to investigate heterogeneous ice nucleation in model systems where an electric field acts on water molecules within 10-20 ? of a surface. Two different water models (the six-site and TIP4P/Ice models) are considered, and in both cases, it is shown that a surface field can serve as a very effective ice nucleation catalyst in supercooled water. Ice with a ferroelectric cubic structure nucleates near the surface, and dipole disordered cubic ice grows outward from the surface layer. We examine the influences of temperature and two important field parameters, the field strength and distance from the surface over which it acts, on the ice nucleation process. For the six-site model, the highest temperature where we observe field-induced ice nucleation is 280 K, and for TIP4P/Ice 270 K (note that the estimated normal freezing points of the six-site and TIP4P/Ice models are ～289 and ～270 K, respectively). The minimum electric field strength required to nucleate ice depends a little on how far the field extends from the surface. If it extends 20 ?, then a field strength of 1.5 × 10(9) V/m is effective for both models. If the field extent is 10 ?, then stronger fields are required (2.5 × 10(9) V/m for TIP4P/Ice and 3.5 × 10(9) V/m for the six-site model). Our results demonstrate that fields of realistic strength, that act only over a narrow surface region, can effectively nucleate ice at temperatures not far below the freezing point. This further supports the possibility that local electric fields can be a significant factor influencing heterogeneous ice nucleation in physical situations. We would expect this to be especially relevant for ice nuclei with very rough surfaces where one would expect local fields of varying strength and direction.     

Comment on "HCl adsorption on ice at low temperature: a combined X-ray absorption, photoemission and infrared study" by P. Parent, J. Lasne, G. Marcotte and C. Laffon, Phys. Chem. Chem. Phys., 2011, 13, 7142

On the basis of NEXAFS, photoemission and FTIR spectra of ice films with low doses of adsorbed HCl, the authors of the PCCP paper "HCl adsorption on ice at low temperature: a combined X-ray absorption, photoemission and infrared study", Phys. Chem. Chem. Phys., 2011, 13, 7142, have come to conclusions regarding the behavior of submonolayer amounts of HCl at 50 K that contradict published results of the authors of this Comment. Our purpose is to argue that the conclusion, attributed going forward to PLML (authors' initials), that nearly 100% of HCl ionizes for dosage levels near to 0.16 monolayer (ML) or 0.3 Langmuir (L) at 50 K is questionable. Rather, we reaffirm our conclusions of much lower levels of ionization for similar temperatures and HCl dosages based on reactive ion scattering (RIS) and low energy sputtering (LES) data for ice films and FTIR spectra of ice nanocrystals. A second current paper by Ayotte et al., J. Phys. Chem. A, 2011, 115, 6002, that largely parallels in method and results the RAIR spectroscopy of PLML, is also given special notice.     

Crystallization, melting, and structure of water nanoparticles at atmospherically relevant temperatures

Water nanoparticles play an important role in atmospheric processes, yet their equilibrium and nonequilibrium liquid-ice phase transitions and the structures they form on freezing are not yet fully elucidated. Here we use molecular dynamics simulations with the mW water model to investigate the nonequilibrium freezing and equilibrium melting of water nanoparticles with radii R between 1 and 4.7 nm and the structure of the ice formed by crystallization at temperatures between 150 and 200 K. The ice crystallized in the particles is a hybrid form of ice I with stacked layers of the cubic and hexagonal ice polymorphs in a ratio approximately 2:1. The ratio of cubic ice to hexagonal ice is insensitive to the radius of the water particle and is comparable to that found in simulations of bulk water around the same temperature. Heating frozen particles that contain multiple crystallites leads to Ostwald ripening and annealing of the ice structures, accompanied by an increase in the amount of ice at the expense of the liquid water, before the particles finally melt from the hybrid ice I to liquid, without a transition to hexagonal ice. The melting temperatures T(m) of the nanoparticles are not affected by the ratio of cubic to hexagonal layers in the crystal. T(m) of the ice particles decreases from 255 to 170 K with the particle size and is well described by the Gibbs-Thomson equation, T(m)(R) = T(m)(bulk) - K(GT)/(R - d), with constant K(GT) = 82 ± 5 K·nm and a premelted liquid of width d = 0.26 ± 0.05 nm, about one monolayer. The freezing temperatures also decrease with the particles' radii. These results are important for understanding the composition, freezing, and melting properties of ice and liquid water particles under atmospheric conditions.     

Melting temperature of ice Ih calculated from coexisting solid-liquid phases

We carried out molecular-dynamics simulations by using the two-phase coexistence method with the constant pressure, particle number, and enthalpy ensemble to compute the melting temperature of proton-disordered hexagonal ice I(h) at 1-bar pressure. Four models of water were considered, including the widely used TIP4P [W. L. Jorgensen, J. Chandrasekha, J. D. Madura, R. W. Impey, and M. L. Klein, J. Chem. Phys.79, 926 (1983)] and TIP5P [M. W. Mahoney and W. L. Jorgensen J. Chem. Phys.112, 8910 (2000)] models, as well as recently improved TIP4P and TIP5P models for use with Ewald techniques-the TIP4P-Ew [W. Horn, W. C. Swope, J. W. Pitera, J. C. Madura, T. J. Dick, G. L. Hura, and T. Head-Gordon, J. Chem. Phys.120, 9665 (2004)] and TIP5P-Ew [S. W. Rick, J. Chem. Phys.120, 6085 (2004)] models. The calculated melting temperature at 1 bar is T(m) = 229 +/- 1 K for the TIP4P and T(m) = 272.0 +/- 0.6 K for the TIP5P ice I(h), both are consistent with previous simulations based on free-energy methods. For the TIP4P-Ew and TIP5P-Ew models, the calculated melting temperature is T(m) = 257.0 +/- 1.1 K and T(m) = 253.9 +/- 1.1 K, respectively.     

Binary nucleation rates for ethanol/water mixtures in supersonic Laval nozzles

Although the conditions corresponding to the onset of condensation of aqueous-alcohol mixtures have been measured in supersonic nozzles [B. E. Wyslouzil et al., J. Chem. Phys. 113, 7317 (2000)], the true nucleation rates have not. Here, we propose a new analytical method to estimate the temperature, the concentrations of condensable species in both the vapor and the liquid phases, and the amount of the condensate using only the measured static pressure profiles in the nozzle. We applied the method to ethanol/water (CH(3)CH(2)OH/D(2)O or CH(3)CH(2)OD/D(2)O) mixtures and confirmed that the aerosol volume fractions derived from pressure measurements and small angle neutron scattering measurements are in very good agreement when this method is used. Combining the results from the pressure measurements with the number densities of the condensed droplets, measured either by small angle neutron or small angle x-ray scattering, we determined the first quantitative ethanol/water binary nucleation rates in the supersonic nozzle at a temperature of 229±1?K.     

Collisional relaxation of O2(X3Σg(-), υ = 1) and O2(a1Δg, υ = 1) by atmospherically relevant species

Laboratory measurements are reported of the rate coefficient for collisional removal of O(2)(X(3)Σ(g)(-), υ = 1) by O((3)P), and the rate coefficients for removal of O(2)(a(1)Δ(g), υ = 1) by O(2), CO(2), and O((3)P). A two-laser method is employed, in which the pulsed output of the first laser at 285 nm photolyzes ozone to produce oxygen atoms and O(2)(a(1)Δ(g), υ = 1), and the output of the second laser detects O(2)(a(1)Δ(g), υ = 1) via resonance-enhanced multiphoton ionization. The kinetics of O(2)(X(3)Σ(g)(-), υ = 1) + O((3)P) relaxation is inferred from the temporal evolution of O(2)(a(1)Δ(g), υ = 1), an approach enabled by the rapid collision-induced equilibration of the O(2)(X(3)Σ(g)(-), υ = 1) and O(2)(a(1)Δ(g), υ = 1) populations in the system. The measured O(2)(X(3)Σ(g)(-), υ = 1) + O((3)P) rate coefficient is (2.9 ± 0.6) × 10(-12) cm(3) s(-1) at 295 K and (3.4 ± 0.6) × 10(-12) cm(3) s(-1) at 240 K. These values are consistent with the previously reported result of (3.2 ± 1.0) × 10(-12) cm(3) s(-1), which was obtained at 315 K using a different experimental approach [K. S. Kalogerakis, R. A. Copeland, and T. G. Slanger, J. Chem. Phys. 123, 194303 (2005)]. For removal of O(2)(a(1)Δ(g), υ = 1) by O((3)P), the upper limits for the rate coefficient are 4 × 10(-13) cm(3) s(-1) at 295 K and 6 × 10(-13) cm(3) s(-1) at 240 K. The rate coefficient for removal of O(2)(a(1)Δ(g), υ = 1) by O(2) is (5.6 ± 0.6) × 10(-11) cm(3) s(-1) at 295 K and (5.9 ± 0.5) × 10(-11) cm(3) s(-1) at 240 K. The O(2)(a(1)Δ(g), υ = 1) + CO(2) rate coefficient is (1.5 ± 0.2) × 10(-14) cm(3) s(-1) at 295 K and (1.2 ± 0.1) × 10(-14) cm(3) s(-1) at 240 K. The implications of the measured rate coefficients for modeling of atmospheric emissions are discussed.     

Computer simulation study of metastable ice VII and amorphous phases obtained by its melting

Molecular dynamics simulations of metastable ice VII and cubic ice Ic are carried out in order to examine (1) the ability of commonly used water interaction potentials to reproduce the properties of ices, and (2) the possibility of generating low-density amorphous (LDA) structures by heating ice VII, which is known to transform to LDA at approximately 135 K at normal pressure [S. Klotz, J. M. Besson, G. Hamel, R. J. Nelmes, J. S. Loveday, and W. G. Marshall, Nature (London) 398, 681 (1999)]. We test four simple empirical interaction potentials of water: TIP4P [W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. Impey, and M. L. Klein, J. Chem. Phys. 79, 926 (1983)], SPC/E [H. J. C. Berendsen, J. R. Grigera, and T. P. Straatsma, J. Phys. Chem. B 91, 6269 (1987)], TIP5P [M. W. Mahoney and W. L. Jorgensen, J. Chem. Phys. 112, 8910 (2000)], and ST2 [F. H. Stillinger and A. Rahman, J. Chem. Phys. 60, 1545 (1974)]. We have found that TIP5P ice VII melts at 210 K, TIP4P at 90 K, and SPC/E at 70 K. Only TIP5P water after transition has a structure similar to that of LDA. TIP4P and SPC/E have almost identical structures, dissimilar to any known water or amorphous phases, but upon heating both slowly evolve towards LDA-like structure. ST2 ice VII is remarkably stable up to 430 K. TIP4P and SPC/E predict correctly the cubic ice collapse into a high-density amorphous ice (HDA) at approximately 1 GPa whereas TIP5P remains stable up to approximately 5 GPa. The densities of the simulated ice phases differ significantly, depending on the potential used, and are generally higher than experimental values. The importance of proper treatment of long-range electrostatic interactions is also discussed.     

A cold and slow molecular beam

Employing a two-stage cryogenic buffer gas cell, we produce a cold, hydrodynamically extracted beam of calcium monohydride molecules with a near effusive velocity distribution. Beam dynamics, thermalization and slowing are studied using laser spectroscopy. The key to this hybrid, effusive-like beam source is a "slowing cell" placed immediately after a hydrodynamic, cryogenic source [Patterson et al., J. Chem. Phys., 2007, 126, 154307]. The resulting CaH beams are created in two regimes. In one regime, a modestly boosted beam has a forward velocity of v(f) = 65 m s(-1), a narrow velocity spread, and a flux of 10(9) molecules per pulse. In the other regime, our slowest beam has a forward velocity of v(f) = 40 m s(-1), a longitudinal temperature of 3.6 K, and a flux of 5 × 10(8) molecules per pulse.     

High resolution overtone spectroscopy of the acetylene van der Waals dimer, ((12)C2H2)2

CW-cavity ring down spectroscopy was used to record in a free jet expansion the spectrum of the absorption band in ((12)C(2)H(2))(2) with origin at 6547.6 cm(-1). It is a perpendicular band and corresponds to 2CH excitation in the hat unit of the T-shaped dimer. Calibration (better than ±1 × 10(-3) cm(-1) accuracy) and ring-down time (130 μs) were improved compared to a previous contribution (Didriche et al. Mol. Phys., 2010, 108, 2158-2164). A line-by-line analysis was achieved. Three series of lines were identified involving levels with A(1)(+), E(+) and B(1)(+) ground state tunneling symmetries, confirming the spectral and symmetry analyses reported in the literature for the 1CH excitation band (Fraser et al. J. Chem. Phys., 1988, 89, 6028-6045). 164 vibration-rotation-tunneling lines were assigned in the K(a)' - K(a)' = 2 - 3, 0 - 1, 2 - 1 and 4 - 3 sub-bands and effective rigid rotor vibration-rotation constants were obtained by simultaneously fitting 1CH and 2CH lines with the same symmetry series. Perturbations affecting the K(a) stacks, in particular, are reported. The tunneling frequency in 2CH is estimated to be ν(tunn)(2CH) = 270 MHz for the K(a) = 0 stack. The rotational temperature is determined to be 23 K from relative line intensities and the lifetime of the dimer in the 2CH hat state is estimated to be 1 ns from individual line widths.     

Melting points and thermal expansivities of proton-disordered hexagonal ice with several model potentials

A method of free energy calculation is proposed, which enables to cover a wide range of pressure and temperature. The free energies of proton-disordered hexagonal ice (ice Ih) and liquid water are calculated for the TIP4P [J. Chem. Phys. 79, 926 (1983)] model and the TIP5P [J. Chem. Phys. 112, 8910 (2000)] model. From the calculated free energy curves, we determine the melting point of the proton-disordered hexagonal ice at 0.1 MPa (atmospheric pressure), 50 MPa, 100 MPa, and 200 MPa. The melting temperatures at atmospheric pressure for the TIP4P ice and the TIP5P ice are found to be about T(m)=229 K and T(m)=268 K, respectively. The melting temperatures decrease as the pressure is increased, a feature consistent with the pressure dependence of the melting point for realistic proton-disordered hexagonal ice. We also calculate the thermal expansivity of the model ices. Negative thermal expansivity is observed at the low temperature region for the TIP4P ice, but not for the TIP5P ice at the ambient pressure.     

Photocatalytic reactions at the graphite/ice interface

We have studied the photodissociation of water molecules at the ice/graphite interface in the presence of submonolayer amounts of potassium. The ice films were grown at cryogenic temperatures and ultrahigh vacuum conditions. They are transparent for 240-900 nm photons, but the strong light absorption in the uppermost layers of the graphite substrate generates energetic charge carriers that may drive photoreactions at the interface. Similar schemes have been demonstrated and investigated before in the monolayer regime, without an ice layer (D. Chakarov, M. Gleeson and B. Kasemo, J. Chem. Phys., 2001, 115, 9477). Here, using ice films of tens of monolayer thickness and with different morphologies, we have investigated the confinement effects due to the ice layer, and the ice permeability for reaction products.     

Fluorescence studies on phenylene moieties embedded in a framework of periodic mesoporous organosilica

The excited state characteristics of phenylene (Ph)-bridged periodic mesoporous organosilica (PMO) powders with crystal-like and amorphous wall structures are investigated. Crystal-like Ph-PMO has a molecular ordering of the bridging organic moieties with intervals of 0.76 and 0.44 nm parallel and perpendicular to the mesochannel direction, respectively, whereas amorphous Ph-PMO has no molecular-level periodicity in the wall. Fluorescence from the exciton and excimer of the Ph moieties and the defect center in the silicate network were detected at room temperature, but fluorescence from the excimer and the defect center were not detected at 77 K for crystal-like Ph-PMO dispersed in a methanol/ethanol mixed solvent. The decay curve of the exciton fluorescence of crystal-like Ph-PMO at room temperature was analyzed successfully using a one-dimensional diffusion model quenched by the defect center and the excimer site. The results were discussed in comparison with those for the crystal-like biphenylene-bridged PMO reported in the preceding paper (Yamanaka et al., Phys. Chem. Chem. Phys., 2010, 12, 11688-11696). The existence of excited states with various conformations including ground state dimers or aggregates of the Ph moieties was suggested for amorphous Ph-PMO. It was clearly apparent that the differences in the excited state dynamics reflected the differences in the molecular-level structure in the wall.     

Nitrogen nucleation in a cryogenic supersonic nozzle

We follow the vapor-liquid phase transition of N(2) in a cryogenic supersonic nozzle apparatus using static pressure measurements. Under our operating conditions, condensation always occurs well below the triple point. Mean field kinetic nucleation theory (MKNT) does a better job of predicting the conditions corresponding to the estimated maximum nucleation rates, J(max) = 10(17±1) cm(-3) s(-1), than two variants of classical nucleation theory. Combining the current results with the nucleation pulse chamber measurements of Iland et al. [J. Chem. Phys. 130, 114508-1 (2009)], we use nucleation theorems to estimate the critical cluster properties. Both the theories overestimate the size of the critical cluster, but MKNT does a good job of estimating the excess internal energy of the clusters.     

Molecular dynamics simulations of surface-initiated melting of nitromethane

The melting of nitromethane initiated at solid-vacuum interfaces has been investigated using molecular dynamics nvt simulations with a realistic force field [D. C. Sorescu et al., J. Phys. Chem. B 104, 8406 (2000)]. The calculated melting point (251+/-5 K) is in good agreement with experiment (244.73 K) and values obtained previously (approximately 255.5 and 266.5+/-8 K) using other simulation methods [P. M. Agrawal et al., J. Chem. Phys. 119, 9617 (2003)]. Analyses of the molecular orientations and diffusion during the simulations as functions of the distance from the exposed surfaces show that the melting is a direct crystal-to-liquid transition, in which the molecules first gain rotational freedom, then mobility. There is a slight dependence of the melting temperature on the exposed crystallographic face.