Reviewer 1 Feedback

[mostafa-razavi]

1_reviewer_attachment_1_1525143789_convrt.pdf

[ramess101]

This manuscript presents the results regarding the applications of the thermodynamics integration (TI) method to calculate the coexistence properties of molecular liquids, i.e. water, ethane, and n-dodecane using canonical NVT simulations. The authors have applied the TI method along isochore and isothermal pathways (the ITIC method). By comparing their results with the NIST data, they show that the ITIC method is valid for reduced temperatures (Tr) of less than 0.85, but it fails for higher Tr. The main conclusion of the manuscript is to perform coexistence calculation with a combination of Monte Carlo (MC) methods such as Gibbs Ensemble MC (GEMC) or Grand Canonical MC (GCMC) and ITIC in order to cover the entire range of industrially relevant temperatures. In addition, the authors conclude that if a single method is preferred, the ITIC can easily be implemented from Tr ranging from 0.85 to 0.45.

The TI method is a well-known method that has been used extensively (see Vega et al. J. Phys.: Condensed Matter 20, 153101, (2008) and Ref. 9 of the current paper). However, most equations in this manuscript are standard and it is hard to find significant development. Regarding Eq. 8, it would be useful to discuss the approximations related to rotational, vibrational, and electronic contributions. Eqs. 14-16, together with the fixed-point iteration method, are used to determine the vapor density and vapor pressure, where the authors discuss the significance of virial coefficients in the fixed-point iteration path.

The authors have applied the ITIC method for only ethane, n-dodecane, and water. However, they have not addressed how the ITIC method works for more systems (see Paluch et al. Ind. Eng. Chem. Res. 48, 4533 (2008), where the authors studied ethane, n-octane, cyclohexane, 2,5- dimethylhexane, 1-propanol, and water). Moreover, a comparison of the ITIC method with the grand-canonical Transition-Matrix Monte Carlo (GC-TMMC) would be useful because the GC-TMMC method has yielded limiting uncertainties in the saturated vapor density and pressure that were significantly smaller, by an order of magnitude in some instances, than those of the GEMC method. In addition, Paluch et al. found that while usually CG-TMMC outperforms GEMC, GEMC can also outperform GC-TMMC. Therefore, the authors have not provided sufficient proof about the advantage of ITIC over other methods such as GC-TMMC.

The authors mention that using flexible bonds results in the systematic discrepancy from exact values (black line in Figure 15) as well as large uncertainties and suggest that one should avoid flexible bonds in simulations at very low densities. However, no sufficient physical insight is presented.

Regarding the application of ITIC for phase coexistence, the authors suggest that at least nine data points are required on the isotherm and three data points on the isochore for each saturation point. Noting that the highest temperature state points serve on both isochores and isotherm, they conclude that a total of 19 state points are required to obtain five coexistence conditions. However, the authors have not discussed this recommendation sufficiently. For instance, if a single ITIC method is preferred, the critical properties may not be obtained via ITIC: the authors have used the law of rectilinear diameter to estimate the critical temperatures and densities. However, the extrapolation to the critical point using the scaling law and the law of rectilinear diameter does not hold well when the Tr is not close to one. Nevertheless, the authors have decided to use this method to estimate the critical properties without appropriate comparison of their results with GCMC or GEMC results (see Binder, Molecular Physics 108, No. 14, 1797–1815 and Dinpajooh et al. J. Chem. Phys. 143, 114113 (2015)). I happened to notice that the critical properties for TIP4P/2005 water predicted by the ITIC method are very far from the NIST Standard Reference Simulation website data. The authors might want to compare the ITIC results with the corresponding GEMC or GCMC results to allow the reader to assess the accuracy of the ITIC method. Therefore, their conclusions about using one preferred single method do not necessarily apply to predicting the critical properties.

Overall, the manuscript is hard to follow due to the large number of figures (about 50, including the supplement) and perhaps tables (11 tables including the supplement) and insufficient discussions. Therefore, I cannot recommend publication in JCP in the current version. An extended revision with detailed discussion might be presented in a more technical journal.

[ramess101]

The TI method is a well-known method that has been used extensively (see Vega et al. J. Phys.: Condensed Matter 20, 153101, (2008) and Ref. 9 of the current paper).

As mentioned in Issue #8, we need to clarify how ITIC is unique.

However, most equations in this manuscript are standard and it is hard to find significant development.

Do we want to move the derivation to an appendix and just focus on the final equations in the text? Since most of these equations are standard, I don't think we want to purport that we are developing the method in this study.

Regarding Eq. 8, it would be useful to discuss the approximations related to rotational, vibrational, and electronic contributions.

I am not sure what he is referring to here. Maybe the ideal gas helmholtz free energy?

Eqs. 14-16, together with the fixed-point iteration method, are used to determine the vapor density and vapor pressure, where the authors discuss the significance of virial coefficients in the fixed-point iteration path.

Again, we need to make it clear that we are solving for Tsat as well.

[ramess101]

The authors have applied the ITIC method for only ethane, n-dodecane, and water. However, they have not addressed how the ITIC method works for more systems (see Paluch et al. Ind. Eng. Chem. Res. 48, 4533 (2008), where the authors studied ethane, n-octane, cyclohexane, 2,5- dimethylhexane, 1-propanol, and water).

As discussed in Issue #1, we are going to present results for other compounds. Again, I have results already for neopentane, and iso-alkanes that were obtained with GROMACS. It could be useful to include those as well.

Moreover, a comparison of the ITIC method with the grand-canonical Transition-Matrix Monte Carlo (GC-TMMC) would be useful because the GC-TMMC method has yielded limiting uncertainties in the saturated vapor density and pressure that were significantly smaller, by an order of magnitude in some instances, than those of the GEMC method. In addition, Paluch et al. found that while usually CG-TMMC outperforms GEMC, GEMC can also outperform GC-TMMC. Therefore, the authors have not provided sufficient proof about the advantage of ITIC over other methods such as GC-TMMC.

This point is addressed in Issue #2. However, I still think we do not want to attempt a rigorous comparison, because then we will need to make sure we are implementing GCMC properly as well. I think we just need to:

1. explain why it is advantageous to avoid insertion moves at low T
2. verify that ITIC yields similar results to GEMC and GCMC at higher T
3. discuss at a superficial level the relative speed of the different methods, i.e. factors that affect how fast one method is compared to the others. Mention that ITIC is easily run in parallel (i.e. multiple jobs), can utilize MD codes (that are highly parallelized), but that GCMC can use smaller systems, etc.
4. demonstrate that the ITIC uncertainties are NOT smaller than GCMC, but they are small enough for most purposes

[ramess101]

The authors mention that using flexible bonds results in the systematic discrepancy from exact values (black line in Figure 15) as well as large uncertainties and suggest that one should avoid flexible bonds in simulations at very low densities. However, no sufficient physical insight is presented.

I can think of three possible points to make:

1. The tested MD packages do not properly account for flexible bonds (hopefully this isn't true)
2. The potential function is different, i.e. the de Pablo harmonic bonds changes the states that are sampled
3. Fluctuations and, therefore, uncertainties are larger because of the additional degree of freedom.

[ramess101]

Regarding the application of ITIC for phase coexistence, the authors suggest that at least nine data points are required on the isotherm and three data points on the isochore for each saturation point. Noting that the highest temperature state points serve on both isochores and isotherm, they conclude that a total of 19 state points are required to obtain five coexistence conditions. However, the authors have not discussed this recommendation sufficiently.

True. Did you ever try fewer state points or more state points? I don't think we need to show this for the simulation examples, just for the REFPROP validation. I think it would suffice to mention some other set-ups that you tested, just briefly.

For instance, if a single ITIC method is preferred, the critical properties may not be obtained via ITIC: the authors have used the law of rectilinear diameter to estimate the critical temperatures and densities. However, the extrapolation to the critical point using the scaling law and the law of rectilinear diameter does not hold well when the Tr is not close to one. Nevertheless, the authors have decided to use this method to estimate the critical properties without appropriate comparison of their results with GCMC or GEMC results (see Binder, Molecular Physics 108, No. 14, 1797–1815 and Dinpajooh et al. J. Chem. Phys. 143, 114113 (2015)). I happened to notice that the critical properties for TIP4P/2005 water predicted by the ITIC method are very far from the NIST Standard Reference Simulation website data. The authors might want to compare the ITIC results with the corresponding GEMC or GCMC results to allow the reader to assess the accuracy of the ITIC method.

As discussed in Issue #7, we might want to remove the critical points.

Therefore, their conclusions about using one preferred single method do not necessarily apply to predicting the critical properties.

However, if we remove the critical point results then we cannot really recommend one preferred single method, unless we clarify that we would not recommend using ITIC for the critical points.

I actually think that we could keep the critical point results and point out that these equations are not valid at low Tr, therefore, the extrapolation to the critical point from ITIC might not be reliable. However, we can show that the critical points are still fairly reasonable in some instances. We could also explain that it is probably best to exclude the lowest Tr ITIC points from the fit, maybe just use the highest three temperatures.

[ramess101]

Overall, the manuscript is hard to follow due to the large number of figures (about 50, including the supplement) and perhaps tables (11 tables including the supplement) and insufficient discussions.

This is discussed in Issue #5.

Therefore, I cannot recommend publication in JCP in the current version. An extended revision with detailed discussion might be presented in a more technical journal.

I think Fluid Phase Equilibria would be a good fit.

[mostafa-razavi]

The TI method is a well-known method that has been used extensively (see Vega et al. J. Phys.: Condensed Matter 20, 153101, (2008)).

This paper is now cited in Introduction