ramess101 / MBAR_ITIC

Merging the MBAR and ITIC methodologies to optimize a TraPPE force switch model
0 stars 2 forks source link

Virial Coefficients #3

Open ramess101 opened 7 years ago

ramess101 commented 7 years ago

Once we start implementing the true ITIC optimization we will need to be more careful with respect to the Virial coefficients. We will either want to run those simulations in Gromacs/Cassandra/Towhee, use DIPPR/REFPROP, or the Kofke approach.

ramess101 commented 7 years ago

@jrelliottoh discussed virial coefficients in Issue #1

For future discussion on how we want to get virials, let's continue here

mostafa-razavi commented 7 years ago

Here is a B2 and B3 comparison which could help us determine virial coefficients more carefully:

Note: For each temperature, I simulated 400 TraPPE-UA ethane molecules in GOMC at these 4 low densities: 1) rho=HighRho/28 2) rho=HighRho/21 3) rho=HighRho/14 4) rho=HighRho/7 "Lin" means linear regression and "Quad" means quadratic regression. For example, Lin1,3 means I used points 1 and 3 to do a linear regression.

b23

Conclusion1 Kofke and Siepmann results are in perfect agreement. One way to solve B2 and B3 issue completely is to just use Kofke's code, make a correlation for B2 and B3 and use it inside ITIC code. The benefit of this approach is that it's fast and it gives us true simulated B2 and B3 from Tr=1.2 all the way to Tr=0.45.

An alternative approach is to get B2 and B3 at Tr=1.2 and Tr=0.95 from (Z-1)/rho regression and fit to some sort of correlation and see if that will predict B2 and B3 at TrSat's of interest. (Virial coefficients are not require at all temperatures. I will post a detailed comment about "At what TrSat's do we need B2 and/or B3?")

@jrelliottoh , @ramess101 What do you guys think about these two approaches?

Conclusion2 If we want to use regression to get B2 and B3, the best way seems to be a linear regression using 1,4 or 1,3 or 2,4 points.

ramess101 commented 7 years ago

Conclusion 1: Coming at this from a force field optimization stand point, I would like for this approach to be amenable to MBAR. For example, with the linear regression method we can use MBAR to predict Z-1/rho vs rho for a given force field. Is the Kofke approach as simple and fast for reevaluating B2 and B3 with a different force field? If so, then I am in favor of Kofke. Otherwise, at least for rapid force field parameterization, I am in favor of the regression approach or even just using the NIST values. Once the force field is optimized, then I would be in favor of whatever is most accurate, likely the Kofke approach.

Conclusion 2: I would suspect that the best approach would be to do a linear regression of all four points. However, a plot of Z-1/rho vs rho would help visualize this. Also, this plot could help understand why we see such sporadic results in B2 and B3 for different linear/quadrature 1,2,3,4 combinations. So could you actually plot the Z-1/rho vs rho results? If you can give us an idea of error bars for the simulations that would be helpful as well. Then, could you also include the Kofke B2 points as the intercept with the Kofke B3 points as the slope? This way we could see if the problem is our assumption of linearity or if it is the simulation output being noisy or if the Z values are just not accurate or some other reason. So I am envisioning a plot where for each Tr we have the 4 simulation points with respect to rho and then we have the Kofke point and line for each Tr as well. Does that make sense?

jrelliottoh commented 7 years ago

Hi Rich, FYI, we have an update about B2, B3. The B2 we get from low density extrapolation agree quite well with Kofke's values (less than 1% deviation). Unfortunately, the B3's are about 25% too high. This happens because the simulations at hiRho/7 are significantly influenced by B4. We can get accurate estimates of Helmholtz energy using our B3eff estimates up to hiRho/7, but B3eff is not equal to B3. The good news is that we do not need B3 for saturation values at Tr<0.86. So we are focusing our ITIC paper on Tr<0.86 for now. Hopefully it will be wrapped up soon.

ramess101 commented 7 years ago

That is encouraging that B2 agrees well with Kofke's values. Just to make sure I understand you correctly, although B3eff is not the same as B3 it is still useful for the ITIC approach, right? In other words, because it is an effective parameter that accounts for the B4 term as well we should still get accurate saturation values. So then wouldn't we be able to use the ITIC approach above 0.86? I understand that below ~0.86 you don't need B3, but if B3eff gives accurate estimates of rhov, Pv, rhol above 0.86 that could be important information as well.