computer simulations 
 exist as a way to connect theories of the 
 microscopic world to experimental results . they help validate our 
 understanding of microscopic phenomena and propose new directions 
 for investigation . such simulations employ molecular models , and the 
 development of molecular models and their component parameters has 
 a long history . one of the primary motivations in their continued development is 
 the increased ability of a molecular model to quantitatively reproduce 
 and/or predict experimental properties of interest . our depiction 
 of molecules in simulations has been grounded by 
 several experimental observables . these include densities , enthalpies , 
 liquid structure , molecular geometry , phase coexistence boundaries , 
 and more . in the case of water , anomalous properties such as the temperature 
 of maximum density provide key motivation . whether the 
 simulation is classical fixed - charge , polarizable , or even quantum 
 mechanical ( qm ) , models used often are unable to reproduce all the 
 observables of interest and , thus , must strike a balance among fitting 
 the selected properties , avoiding too many poorly defined parameters , 
 and keeping the computational cost low . generalized models , such as 
 opls - aa , gaff , and cgenff , typically focus on neat 
 liquid density and heats of vaporization as experimental targets for 
 optimization of lennard  jones ( lj ) parameters coupled with 
 partial charges assigned using fits to hf/6 - 31 g * electrostatic potentials . 
 while models made from these parameter sets can be used to study neat 
 liquids , one of the main intents is for the representation of molecules 
 in other environments , such as a ligand interacting with a protein 
 or a solute transferring between different chemical phases . when using 
 such a model outside of a neat liquid , it is understandable that one 
 may observe a systematic error in a new property of interest . in the 
 case of molecular transfer , researchers have explored a variety of 
 approaches to reduce such systematic errors , and these generally amount 
 to fitting the models to get the answer desired , sometimes coupled to hydrated environment qm calculations , and other times through statistical means by 
 introducing additional water - mixing parameters to the lj potential 
 terms . for researchers interested in 
 simultaneous reproduction of experimental observables in multiple 
 disparate environments , which is sometimes considered impossible for 
 classical fixed - charge models we are interested in both 
 the neat liquid and air - to - water transfer 
 properties of molecules with hydroxyl moieties . recent studies on 
 the hydration of alcohols using common force field parameters have 
 shown systematic errors in ghyd and other thermodynamic quantities , in particular the static dielectric 
 constant , (0 ) . are systematic 
 errors in (0 ) and ghyd coupled 
 to one another in some way ? it has been suggested that functional 
 groups in typical force fields are underpolarized to some degree , though this should not be too surprising given that partial charges 
 are typically fit to reproduce a gas - phase hf/6 - 31 g * electrostatic 
 potential , despite this charge model s fortuitous propensity 
 for overpolarizing the gas - phase electric field . for solution - phase 
 simulations , it would be preferable to have an accurate degree of 
 polarization for the condensed - phase fluid encoded within the atomic 
 partial charges of the molecule while maintaining the typical neat 
 liquid density and enthalpy of vaporization ( hvap ) . recent studies have begun using the experimentally 
 measured (0 ) of a neat liquid as a fitting target for systematically 
 scaling the degree of polarization of classical fixed - charge models . this provides an orthogonal approach to the practice of setting 
 partial charges with a qm electrostatic potential and optimizing lj 
 parameters to fit experimental observables . as it relies on experimental 
 quantities to govern the magnitude of the partial charges , it can 
 potentially provide an independent assessment of the quality of a 
 given qm calculation used to derive partial charges . here , we extend 
 such work with the goal of addressing the following questions : ( 1 ) 
 is it possible to simultaneously capture the density , hvap , and (0 ) for simple alcohol solvents ; 
 ( 2 ) how transferable are the resulting hydroxyl parameters to other 
 more complex alcohols ; and ( 3 ) if an alcohol model reproduces the 
 neat liquid (0 ) , how does that affect the hydration of that 
 alcohol as a solute ? in this work , both extended 
 molecular dynamics simulations and alchemical free energy calculations 
 were used to calculate neat liquid thermodynamic properties and air - to - water 
 transfer free energies for a series of 41 small molecules containing 
 hydroxyl moieties for which both experimental quantities were available . simulations were performed using gromacs 4.5.5 for the neat 
 liquid molecular dynamics and a development version of gromacs 4.6 
 for alchemical free energy calculations on these alcohol systems . hydration free energies were computed using thermodynamic integration 
 ( ti ) , following a combination of previous protocols . initial conformations for starting molecules 
 were taken from structures provided in previous studies , excepting the case of methanol subject to the 
 neat liquid dielectric optimization , which was constructed and minimized 
 using the structure editor in the ucsf chimera visualization system . the gaff small molecule parameters and am1-bcc partial charges were assigned to each structure using the antechamber package ( amber 
 11 version ) . after setup , the resulting 
 structures and topologies were converted to the gromacs format using 
 the acpype python script . for neat liquid 
 calculations , cubic simulation boxes of 256 like small molecules were 
 constructed with gromacs utilities and subject to sets of 150 ps of 
 constant - pressure and constant - temperature equilibration until the 
 target state - point was achieved . again , methanol in the neat liquid 
 optimization was an exception in needing a large cubic simulation 
 box containing 450 molecules in order to accommodate box size fluctuations 
 along the density optimization pathways . for hydration free energy 
 calculations , single copies of the small molecules were solvated in 
 rhombic - dodecahedral tip3p water boxes 
 with at least 1.2 nm of space between any solute atom and a box edge . for all neat liquid calculations , a set of five 50 ns production 
 simulations were performed for each molecule with a unique topology 
 in order to obtain reported (0 ) values and associated standard 
 errors . these calculations were performed at 1 atm and the target 
 temperature from the experimental dielectric constant determination , 
 usually 298.15 k unless otherwise noted . a time step of 2 fs was used with the leapfrog integration algorithm 
 for the equations of motion , and hydrogen atom bond lengths were constrained 
 using the lincs algorithm . rahman 
 barostat and nose  hoover thermostat with time constants of 10 and 1 ps , respectively . smooth particle - mesh 
 ewald with a real - space cutoff of 1.2 
  , spline order of 4 , and energy tolerance of 10 was used for long - ranged electrostatics corrections . lj interactions 
 were switched off smoothly between 8 and 10  , and long - range 
 energy and pressure corrections were applied to these interactions 
 as implemented in gromacs . hydration 
 free energy calculations were performed over a series 
 of separate  windows as described in previous work , with the primary difference being that these 
 simulations were performed at constant pressure as in the neat liquid 
 calculations , but with the langevin thermostat as in previous work . 
 performing the calculations with the parrinello rahman barostat 
 was found to result in numerically equivalent results to previous 
 protocols , while eliminating the need for an extensive constant pressure 
 equilibration and average volume rescaling at the beginning of each 
  window . for our set of 41 alcohols , ti was used to determine 
 both the electrostatic and nonpolar components of solvation rather 
 than the bennett acceptance ratio ( bar ) for technical reasons in the 
 development version of gromacs . for the sampl4 calculations , our standard 
 bar protocol was restored with the use of gromacs 4.6.2 beta 2 . the 
 electrostatic component was calculated from the free energy of turning 
 on the solute partial charges in water , less the free energy of the 
 same charging process in vacuum . the nonpolar component was calculated 
 from the free energy for decoupling the uncharged solute lj interactions 
 from the surrounding water as in previous work . this overall protocol was validated against previous calculations 
 on these 41 small hydroxyl containing solutes and found to be reliable and equivalently accurate . starting with 
 a general model for a solute one would use in a typical alchemical 
 free energy calculation , an iterative optimization process was developed 
 to alter the molecular topology in order to improve agreement with 
 dielectric properties of the neat liquid . this process is a variation 
 of previous work with a focus placed 
 specifically on small molecules containing hydroxyl groups . the simplest 
 of such solutes is methanol , and figure 1 shows 
 an example optimization flow process for methanol . illustration of the iterative 
 optimization process for a methanol 
 model using experimentally measured data . first , the size ( lj  ) 
 is scaled to match experimental density for the neat liquid . this 
 is followed by charge - scaling , to match the experimental dielectric 
 constant , and well - depth ( lj  ) scaling , to match the experimental 
 hvap . the resulting model is tested 
 in a final set of simulations to see how well it matches these three 
 properties , and it is either accepted or returned to the optimization 
 pipeline for further refinement . the goal of this process is to retain or improve the density 
 and 
 hvap agreement that general molecular 
 models are developed upon , while simultaneously improving the agreement 
 with experimental (0 ) values . for simplicity , a three - stage 
 procedure was employed with a linear interpolation in each stage used 
 to target optimal values for these properties . after employing an 
 initial 1 ns simulation to determine the starting model s density , the hydroxyl oxygen  value was scaled up or down by 5% if the 
 density was either greater than or less than the target density window 
 ( 791.4 kg / m to 1.0 kg / m ) respectively . the density of a 1 ns simulation of the resulting model was calculated , 
 and a linear interpolation was performed between these densities to 
 determine a  value that should correspond to the model having 
 the target density . this led to the second stage process 
 of optimizing the (0 ) value ( 33 to 4 ) , calculated using 
 the standard dipole fluctuation formulas , using similar linear interpolations over 15 ns simulations . in order 
 to maintain charge neutrality , all the partial charges were scaled 
 by 5% , effectively scaling the dipole moment of the methanol 
 model in order to raise or lower the (0 ) value . this is illustrated 
 in figure 1 as a saturation of atom and bond 
 color on the methanol molecule . the third stage involved linear 
 interpolation of hvap values to 
 estimate an optimal lj  
 well depth for the hydroxyl oxygen atom . the hvap was calculated from the difference in per - molecule 
 total energies from 10 ns single methanol gas phase and 500 ps liquid 
 phase simulations,1epol is the energetic 
 term corresponding to the work required to polarize the methanol molecule 
 from the gas phase dipole moment to the liquid phase dipole moment,2here , gas and gas were 1.7 d 
 and 3.29  respectively for 
 gaseous methanol . to determine liq from a flexible , yet fixed - charge methanol model , we averaged 
 individual liq values of each molecule from the 
 last frame of the liquid simulation . note that , as the neat liquid 
 optimization system contains 450 alcohol molecules , the eliq for the model is able to converge to 0.1 kcal / mol 
 over this relatively short optimization stage . as in the other stages , 
 this process was iterated until the model reported hvap within the desired tolerance ( 8.95 kcal / mol 
 to 0.59 kcal / mol ) . it should be noted that we only use hvap in this optimization of methanol and do 
 not calculate it for other solutes later in this study . experimental 
 hvap , gas , and 
 gas properties are not available for all subsequent 
 molecules investigated , making such comparisons problematic . after this series of parameter adjustments , a final set of simulations 
 was performed to check if all of these neat liquid properties were 
 within the specified error tolerances for each of their state points . 
 if any of the three properties degraded due to correlations between 
 the parameters and properties not directly associated in the optimization 
 procedure , the series of stages was iterated until the set of properties 
 converged . it should be noted that the experimental values need not 
 be at the same state point . in the case of methanol , the hvap value used was measured at 298.15 k , while 
 the density and (0 ) were measured at 293.15 k. the adjusted topology will simply describe a classical 
 molecular model that captures these experimental quantities at these 
 differing state points . this topology will also not be the only set 
 of parameters to capture the target property values , just the first 
 set found that satisfies the criteria . while the (0 ) 
 optimization procedure could adjust a classical 
 molecular model to reproduce experimental properties for an entire 
 set of molecules , it is not particularly efficient . in our tests developing 
 an optimization approach , we found that this simple linear approach 
 usually only required three to four (0 ) calculations to converge . as converging on (0 ) is the slow step of the process , we actually 
 found the simple line - search strategy to be more efficient than grid - based approaches that required more 
 (0 ) calculations to map out parameter space . one could possibly 
 employ an even more efficient search strategy such as that used in 
 the forcebalance tool , but optimization 
 would still require a computationally intensive search for each molecule 
 of interest or a combined optimization over a set of small molecules . 
 rather than pursue such a specialized route , we chose to test if the 
 knowledge gained from the hydroxyl optimization for methanol was transferrable 
 to other small molecules containing hydroxyl moieties . illustration of the application 
 of optimized methanol polarization 
 to a ( a ) primary alcohol and ( b ) phenolic alcohol . the partial charges 
 on all labeled atoms are scaled by the scaling factor from the methanol 
 optimization with the exception of the nc atoms . in the primary alcohol , 
 the single nc atom acts as a neutralization site . in the phenolic 
 alcohol , the charge difference due to scaling of the hydroxyl group 
 is split over the two nc atoms . to apply the methanol hydroxyl optimization in a general 
 manner , 
 we assume that the hydrogen bonding between hydroxyls in neat methanol 
 is representative of hydroxyl hydrogen bonding in other systems . additionally , 
 we assume that intramolecular perturbations of the electric field 
 about a given hydroxyl group that would alter this behavior are treated 
 by the semiempirical or quantum mechanical approach used in the setting 
 of the partial charges for the model . given these assumptions , we 
 adopt new lj  and  values , taken from the methanol optimization , 
 for each hydroxyl oxygen atom . for the partial charges , we scale all 
 hydroxyl group atoms by the amount observed for methanol , as shown 
 in figure 2 . any non - hydrogen atoms bonded 
 to the  carbon are then treated as neutralization sites ( the 
 nc site in figure 2a ) to compensate for any 
 charge difference introduced by this charge scaling . if there are 
 multiple non - hydrogen atoms bonded to the  carbon ( figure 2b ) , the charge difference between the original and 
 polarized set of hydroxyl atoms is distributed evenly among these 
 bonded atoms . we applied the described optimization protocol to alter 
 the force field parameters for the hydroxyl oxygen in a methanol molecule 
 using the gaff small molecule parameters and am1-bcc partial charges . table 1 shows how the properties for the 
 methanol models successfully balance within the predefined experimental 
 constraints . the density is sitting just within our predefined tolerance 
 ( 1 kg / m ) of the experimental value while the other properties 
 sit well within the kbt window for the hvap and (0 ) 
 = 4 window . the biggest change is seen for the (0 ) value in 
 moving a quantity that has a 40% discrepancy down to a less than 5% 
 difference . to accomplish this , the dipole moment of the methanol 
 model needed an 20% increase in magnitude ( see table 2 ) . these modifications are similar to recently proposed 
 modifications derived from matching kirkwood  buff integrals , 
 though not quite of the same magnitude , and also qualitatively similar to recent changes in the gromos force 
 field for alcohols . our changes indicate 
 that the condensed - phase polarization of methanol is underaccounted 
 for in the am1-bcc charge model and the hf/6 - 31 g * basis - set electrostatic 
 potential upon which it is fit , even 
 given this basis set s 10% to 15% overpolarization relative 
 to the gas phase expectation that is considered fortuitous for condensed - phase 
 simulations . the experimental liquid - phase dipole moment for methanol is currently 
 unknown , so while it is expected to be greater in magnitude than the 
 gas - phase dipole moment , it is unclear if the optimized liquid state 
 electric field is a truer representation of that in real methanol 
 than the field used in other models . we can not expect other 
 properties , such as the density , to remain 
 static in the face of such a large polarization of the liquid - state 
 methanol molecule . thus , the hydroxyl oxygen  expands to counter 
 the increase in density that occurs with such a polarization . this 
 increase in size simultaneously decreases the strength of alcohol 
 model pair interactions and increases the molecular volume . while 
 this change in oxygen atom size acts to counter the increased polarization , 
 this does not mean that the structure of the liquid remains static . figure 3 shows the goo(r ) for methanol before and after the optimization 
 process , alongside a similar curve from neutron - diffraction studies . we see an approximately 0.05  shift of peak maxima to longer 
 distances with the new parameters , in poorer agreement with experimentally 
 supported structural data . we also see that the peaks and troughs 
 are amplified , indicating an increase in the structuring of the liquid , 
 in better agreement with the experimentally supported structural data . these small positive and negative changes , coupled with the fact that 
 the liquid structure can not be unambiguously determined from experiments 
 without using model fits , make it difficult for us to state if the 
 liquid structure is better or worse relative to experiment . we can 
 state that while the liquid density remains the same , the new parameters 
 act to increase the structure of the liquid . if further changes in 
 the structure are desired , one could consider including target data 
 along with condensed - phase geometry optimizations in the fitting process . oxygen  oxygen 
 radial distribution functions for the starting 
 gaff / am1-bcc methanol model ( orange dashes ) , the 
 final dielectric corrected ( dc ) model ( blue ) , and 
 neutron - diffraction studies ( black ) . the topology adjustments in the final dc model make the liquid more 
 structured , but they do little to the radial location of the existing 
 methanol structure of the starting model . in addition to liquid structure , experimental diffusion constant 
 data are available for methanol at standard conditions , 2.28  
 10 cm / s . we computed the neat liquid diffusion constant for this methanol 
 model before and after the optimization process . the original model 
 has a diffusion constant of ( 3.8  0.2 )  10 cm / s , which is considerably faster than experiment . 
 this fast dynamic behavior is quite similar to that seen with the 
 opls - aa force field . after optimization , the diffusion constant drops to ( 1.8  0.2 )  10 cm / s , slower than , but considerably closer to , the experimental 
 value . this change is similar to that seen for the evolution of spc 
 to spc / e water , where the additional polarization in spc / e lowered 
 the diffusion constant closer to the experimental value . as suggested for the model liquid structure , 
 using the diffusion constant , possibly in place of hvap , could be an alternate optimization strategy 
 to further the agreement in dynamical properties of general alcohols . the scaling values listed here were used as general rules in the 
 development of a simple tool for converting such hydroxyl groups in 
 gromacs topologies into similar dielectric corrected ( dc ) forms . the 
 source code for perl and python versions of this utility is provided 
 in the supporting information for interested 
 readers . the optimization improves the dielectric properties of methanol 
 directly but is difficult to apply prospectively to new molecules . if improvement of the (0 ) is desired for models of other molecules 
 with hydroxyl groups , it would be advantageous to not need to undergo 
 such a process for each new molecule . in many cases , such an effort 
 would be impossible due to the lack of relevant experimental data . this is a key concern for researchers wanting to use new molecular 
 models as a step in a computational chemical process for pharmaceutical 
 development , or in any other discovery setting . we therefore decided 
 to test how transferrable these adjustments are to other more complex 
 small molecules . this set , listed in the supporting 
 information , consists of methanol and 40 additional molecules 
 for which the experimental (0 ) and ghyd are known . figure 4a and 4b show the before and after 
 application of dc hydroxyl modifications results for the calculation 
 of the neat liquid (0 ) values at the experimentally measured 
 temperatures . one striking feature of the original general alcohol 
 models is the tight correlation but shallow slope for the calculated 
 (0 ) versus experimental measurements . regression analysis shows 
 the model (0 ) to typically be about two - thirds of the experimental 
 value . applying the dc methanol modifications to all hydroxyl groups scatterplots of the calculated (0 ) values versus 
 experimental 
 quantities for a set of 41 small molecules using ( a ) the gaff small 
 molecule parameters and am1-bcc partial charges and ( b ) those same 
 models after applying the methanol dielectric correction to all hydroxyl 
 groups . for the gaff / am1-bcc models , the data are well - correlated 
 ( r = 0.91 ) to an underpredicting 
 slope of 0.65(3 ) . the dc models show a similar correlation ( r = 0.89 ) but have a slope of 1.07(6 ) , indicating 
 improved agreement with experiment . the general trend for the change in (0 ) values when 
 applying 
 the dc methanol modifications was an approximately 50% increase , closing 
 the gap between the calculations and experimental measurements . there 
 were a few exceptions , the most glaring were cyclohexanol and cycloheptanol . 
 in these neat liquids , the (0 ) changed in the opposite manner , 
 an average 50% decrease . it should be noted that the densities for 
 the starting gaff / am1-bcc model liquids were already 1020 
 kg / m less than experiment , and these values do not get 
  corrected  by the hydroxyl modifications . this indicates 
 a potential issue may reside in the cycloalkane angle and/or torsion 
 parameters , altering the liquid - state dipole orientation distributions 
 in a manner inconsistent with the trends seen for other small molecules . 
 regardless , such occurrences were the exception rather than the norm , 
 as only 4 of the 41 molecules showed a statistically significant decay 
 of (0 ) away from the experimental value , with the noted cycloalcohol 
 pair being the worst offenders . the majority of these alcohols showed 
 dramatic enhancement , indicating that the dc methanol modifications 
 are surprisingly transferable for the purpose of improving model (0 ) 
 values . one of the primary motivations behind improving the dielectric 
 properties of model solutes is the potential for improvement in force 
 field depictions of such solutes in other applications such as hydration 
 and binding free energy calculations . solvation ( and desolvation ) 
 of ligand molecules is a critical aspect of ligand binding processes , 
 and accurate predictions of ghyd are necessary for molecular modeling to be a useful contributor 
 to pharmaceutical design and development . transfer processes are governed 
 by the dielectric permittivity of the mediums in implicit molecular 
 models , and we expect an experimentally relevant depiction of the 
 solute and environment interactions will be necessary for explicit 
 models to be quantitative . the dielectric correction process does 
 not use ghyd information for optimization , 
 but it does target the quality of the model by polarizing hydroxyl 
 groups to a level relevant to a condensed - phase hydrogen - bonding environment . 
 here , we test the dc hydroxyl modifications by comparing alchemical 
 hydration free energy calculations on the original gaff / am1-bcc and 
 modified forms . when performing ghyd calculations on a polarized model , it has been suggested that a 
 work term for the polarization of the molecule ( epol ) be included in final result similar to how it is 
 applied in the hvap calculation . through personal correspondence with the suggesting authors , we 
 were advised not to include such a term when using 
 the am1-bcc partial charges . here , we follow this general advice ; 
 however , because we are systematically polarizing the hydroxyl groups 
 beyond the am1-bcc level , we include an additive per - hydroxyl 
 epol term for polarizing methanol from the am1-bcc 
 partial charge level to that at the end of the optimization , 0.409 
 kcal / mol per - hydroxyl . this term is taken to be additive with the 
 number of hydroxyl groups because we are only modifying the parameters 
 of hydroxyl groups  we leave the remaining model parameters 
 untouched . note that the magnitude of this potential epol term is less than 1 kt , the chosen 
 tolerance for the hvap optimization . this indicates that the new parameters are a valid fit to the computed 
 hvap whether or not the epol term is included . we observe a net benefit whether we 
 apply the epol term or not , though we 
 see the best performance using this correction . this work does not 
 resolve questions on whether or not such a term should be applied 
 in general , but it shows how a simplified consideration 
 of this term can influence hydration free energy calculations . plots 
 and error analysis of all these cases are provided in the supporting information for readers interested 
 in additional comparison clarifications . figure 5a and 5b show the results of ghyd calculations before and after the application 
 of the dc hydroxyl modifications . figure 5a 
 highlights a key concern in hydration calculations of alcohols using 
 gaff and am1-bcc partial charges . when 
 applying the dc hydroxyl modifications , this systematic error disappears , 
 and we are left with accurate calculation results that agree extremely 
 well with experiment . these results indicate that the neat liquid 
 methanol environment from the optimization appears to be a good proxy 
 for an aqueous hydrogen - bonding environment and quantitatively better 
 than the essentially gas - phase environment in an am1-bcc calculation . while we only have two solutes in this set with multiple hydroxyl 
 groups ( 1,2-ethandiol and 1,3-propanediol ) , the dc hydroxyl modifications 
 also appear to avoid potential issues related to nonadditivity of 
 local polarization , taking ghyd values that are nearly 2 and 3 kcal / mol too unfavorable compared 
 to values that are within 0.5 kcal / mol of experiments . as the methanolic environment provided by the dc hydroxyl modifications 
 is beneficial for both neat liquid (0 ) and ghyd values calculated for general alcohols , we are interested 
 in determining if the new lj parameters for the hydroxyl oxygen atom 
 could be used in conjunction with a qm method that better treats the 
 condensed - phase environment s effect on the solute electric 
 field and resulting partial charges . a better accounting of the condensed - phase 
 environment could result in a truer representation of the electrostatics 
 of the solute molecule . ideally , we could set the partial charges 
 for a given solute containing hydroxyl groups by fitting to a qm electrostatic 
 potential calculated using a self - consistent reaction field ( scrf ) 
 or similar environment term , rather than using the modified am1-bcc 
 partial charges . to test this , we applied partial charges from mp2/cc - ptvz 
 scrf calculations to a small subset of 
 phenolic and linear alcohol solutes and compared calculated (0 ) 
 and ghyd values with similar results 
 from gaff / am1-bcc and the dc modified models . scatterplots comparing the calculated to experimental 
 ghyd for 41 alcohols depicted using 
 the ( a ) gaff 
 and am1-bcc parameters and ( b ) the dc hydroxyl modified forms . the 
 gaff / am1-bcc mean signed error ( mse ) of 1.43 kcal / mol is essentially 
 eliminated ( mse of 0.07 kcal / mol ) when using the dc hydroxyl 
 modifications , and the root - mean - square error also drops from nearly 
 3kbt ( 1.52 kcal / mol ) 
 down to under kbt ( 0.45 
 kcal / mol ) . the correlation of the data with experiment improves slightly 
 as well , going from an r of 0.87 to an r of 0.91 when using dc hydroxyls . it is apparent from the improved (0 ) values 
 that the use 
 of an scrf polarizes the alcohols more than the semiempirical am1-bcc 
 approach . still , the more complex qm method does not result in dielectric 
 constants as accurate as our simple application of the polarization 
 changes from a neat liquid optimization of methanol . similarly , the 
 hydration calculation shows improvement when using the more detailed 
 qm method , but agreement with experiment is still not quite to the 
 level of that from our dc hydroxyl modifications . while this is somewhat 
 disappointing from the perspective of the higher - level qm approach , 
 it does indicate that using the scrf improves both dielectric and 
 hydration properties . additionally , the neat liquid optimization rules 
 appear to provide a potential target benchmark for future development 
 of environment terms in qm electric field calculations . comparison 
 of calculated ( a ) (0 ) and ( b ) ghyd with experiment when using gaff / am1-bcc 
 ( red squares ) , the dc hydroxyl form ( blue 
 circles ) , and mp2/cc - ptvz scrf with the dc hydroxyl oxygen 
 lj parameters ( green diamonds ) . by using the new 
 lj parameters , the neat liquid densities for methanol , ethanol , propan-1-ol , 
 phenol , o - cresol , and p - cresol do 
 not degrade when using a condensed - phase quantum calculation . the 
 resulting (0 ) and ghyd are 
 still not quite to the accuracy level of the simple approach outlined 
 in this work . our approach of correcting dielectric 
 properties proved successful for a moderately sized set of small molecules , 
 but our real goal is an approach which is transferable and predictive . thus , it would be helpful to have blind studies to evaluate if this 
 dc hydroxyl modification can be used in predictive situations . to 
 this end , we contributed ghyd results 
 both without and with the dc modifications to the sampl4 hydration 
 prediction event . the results of these contributions are summarized 
 in table 3 . for the actual sampl4 event , the epol contribution to ghyd was not included , but the results with this term are 
 included here for completeness . the reported errors in table 3 were computed as the standard deviation in the 
 mean over 1000 bootstrap trials , with experimental values resampled 
 with added gaussian noise drawn from a distribution with a standard 
 deviation equal to the experimental uncertainty . of the 47 small molecules in the full sampl4 set , 18 contained 
 at least one hydroxyl group . since there are a diverse set of functional 
 groups represented in the sampl molecule sets , we would not expect 
 the application of these dc hydroxyl modifications to result in perfect 
 agreement with experimental values . however , the 0.3 kcal / mol 
 improvement shown in table 3 in hydroxyl containing 
 solutes shown does shift the overall set of predictions closer to 
 experiment by 0.1 to 0.2 kcal / mol . we ran a paired t - test on gaff / am1-bcc versus gaff / am1-bcc dc and found that the latter 
 is significantly better ( p = 0.008 ) . it is interesting 
 to note from the mean signed errors ( mse ) that the dc hydroxyl modifications 
 without the epol term improve the overall 
 results by countering the systematically hydrophobic nature of the 
 general force field , possibly skewing the predictions too hydrophilic . the inclusion of the epol term tempers 
 this effect , yet still results in similar overall prediction quality . while the aim of the sampl experiment is to inform rather than 
 act as a competition between groups , we note that the overall performance 
 of this simple force field + hydroxyl adjustments submission was regularly 
 among the top three predictions across the various performance metrics . these results indicate that further force field 
 improvements using this dc strategy applied to other functional groups 
 could be an effective route toward more accurate predictions of hydration 
 from theoretical methods that rely upon the quality of small molecule 
 force field parameters . we have shown a method 
 to incorporate experimental dielectric response 
 information on neat liquids in classical molecular models of alcohol . 
 by optimizing the lj parameters and polarization of methanol to capture 
 the density , (0 ) , and hvap , we were able to determine simple hydroxyl modulation rules that 
 can be applied in a general manner to other alcohols . we tested how 
 these dc hydroxyl modifications affect the (0 ) for 41 small 
 molecule alcohols , and we found that these changes are transferable 
 and lead to substantially improved model property agreement with experiment . we applied these modifications in the air - to - water transfer of these 
 molecules as well as in a blind prediction event and found that these 
 changes improve agreement with experiment by eliminating an observed 
 systematic hydrophobic bias in the hydration of gaff / am1-bcc alcohols . these simple modifications outperform detailed qm calculations using 
 an scrf environment coupling term and provide a potential benchmark 
 for future developments in qm calculation methods . first , the 
 relative simplicity and general nature of this effort to correct perceived 
 flaws in existing classical molecular mechanics force fields indicates 
 that there is room for improvement in the lj and partial charge representations 
 of standard chemical groups . here , we employed an optimization strategy 
 that utilized accepted experimentally measured thermodynamic properties 
 as target quantities . this work on the simplest of alcohols resulted 
 in a transferable approach that simultaneously improved dielectric 
 and hydration properties of more complex alcohols . additionally , 
 the improvement in ghyd results 
 indicates that the methanolic environment from 
 the hydroxyl optimization is a good analogue for the aqueous hydrogen - bonding 
 condensed - phase environment . using such an environment when optimizing 
 lj and partial charge parameters helps resolve systematic errors in 
 hydroxyl group hydration seen when using accepted classical force 
 fields . we suspect a similar improvement in hydration could occur 
 for other polar functional groups , such as amines and carbonyls , following 
 optimization , as their polar neat liquid environments might also be 
 reasonable analogues to polar aqueous environments . while the 
 optimization process applied here is specific to the 
 originating force field , it appears to be general for a given charge 
 model basis . we have applied a similar optimization to models using 
 the opls - ua force field for alcohols , which are fit to a slightly 
 less detailed basis set , and arrived 
 at similar end - point alterations to those seen here . this indicates 
 that further development in condensed - phase qm electric field calculation 
 methodology , possibly in conjunction with experimentally grounded 
 efforts similar to those presented here , can potentially lead to even 
 more quantitative and insightful future classical molecular simulations .
we present a simple optimization 
 strategy for incorporating experimental 
 dielectric response information on neat liquids in classical molecular 
 models of alcohol . using this strategy , we determine simple and transferable 
 hydroxyl modulation rules that , when applied to an existing molecular 
 parameter set , result in a newly dielectric corrected ( dc ) parameter 
 set . 
 we applied these rules to the general amber force field ( gaff ) 
 to form an initial set of gaff - dc parameters , and we found this to 
 lead to significant improvement in the calculated dielectric constant 
 and hydration free energy values for a wide variety of small molecule 
 alcohol models . 
 tests of the gaff - dc parameters in the sampl4 blind 
 prediction event for hydration show these changes improve agreement 
 with experiment . surprisingly , these simple modifications also outperform 
 detailed quantum mechanical electric field calculations using a self - consistent 
 reaction field environment coupling term . 
 this work provides a potential 
 benchmark for future developments in methods for representing condensed - phase 
 environments in electronic structure calculations .