ramess101 / MBAR_GCMC

3 stars 0 forks source link

Enthalpy of vaporization #10

Closed ramess101 closed 5 years ago

ramess101 commented 5 years ago

@msoroush @jpotoff

I am computing the heat of vaporization using the following equation (from Mick et al.):

image

However, there is a clear bias between the MBAR DeltaHv values and those reported by Mick et al. with histogram reweighting (HR):

image

I was hoping you could help me elucidate the cause for this deviation. Recall that the rholiq, rhovap, and Psat values do not demonstrate any bias:

image

So it appears the bias has to be coming from the Uv-Ul term. Just to verify, do these U values include the intramolecular difference (which can be non-negligible for a high density liquid and a low density vapor)? If so, I think the histogram values that Mohammad provided me with only have U for the non-bonded contribution, correct?

If not, is it possible that the small biases in rhov and/or Psat are magnified when computing DeltaHv? In other words, do you think this could just be a propagation of errors scenario?

mrshirts commented 5 years ago

Well, it's possible the histogram values are the ones that are wrong . . . : )

But yes, one would have to look at intramolecular energy differences if there are any as well.

jpotoff commented 5 years ago

@ramess101 that is weird. The histogram data we sent you has just the intermolecular energies. As far as the effect of combining small errors = a bigger error, I don't think you should see that here. We do see small errors combining into large errors in the calculation of the compressibility factor, though. It would be worth calculating Z and seeing what you get.

ramess101 commented 5 years ago

@mrshirts

But yes, one would have to look at intramolecular energy differences if there are any as well.

I agree that the intramolecular differences SHOULD be included, but they are often neglected. The real purpose of this exercise is to compare MBAR and HR, so what we really need to know is whether they WERE included by Mick et al. @msoroush @jpotoff you mentioned that the histogram data only had intermolecular energies, but do you know if the intramolecular energies were also included for computing DeltaHv? If so, are the files with the intramolecular snapshots easily accessible?

It would be worth calculating Z and seeing what you get.

I will look into this.

jpotoff commented 5 years ago

@ramess101 @msoroush @mrshirts In this case, the change in the intramolecular energies between phases was small (< 1%) and they were neglected in the dHv calculation. We know that for larger, more flexible molecules, that approximation isn't going to work, and one must include the intramolecular contributions to get dHv right. Unfortunately, we did not histogram the intramolecular energies, and getting them would require running new simulations.

ramess101 commented 5 years ago

@jpotoff

OK, well I really am just trying to replicate your values. So if you did not include the intramolecular energies, then that is what I will do. No need to run new simulations. But I will need to elucidate what the cause is for this discrepancy.

As for Z, the uncertainties are large but there is no bias, i.e. our results are consistent:

image

jpotoff commented 5 years ago

@ramess101 Just throwing an idea out there, and maybe this is the issue because our equation isn't written exactly right. In our calculation we are taking the difference of the average energy (and molar volumes) in each phase. -+P(-)

ramess101 commented 5 years ago

@jpotoff

That is almost how I am computing it as well. I am using but I am actually computing the volumes as 1/, maybe that is a problem? I kind of doubt it, but I will try it out tomorrow.

Another idea that might impact low density vapor energies, do you include snapshots in the averaging where the number of molecules is 0? I think there are quite a few and this could impact Uvap and/or Vvap. But this wouldnt explain the divergence at high temperatures.

On Wednesday, October 17, 2018, Jeffrey Potoff notifications@github.com wrote:

@ramess101 https://github.com/ramess101 Just throwing an idea out there, and maybe this is the issue because our equation isn't written exactly right. In our calculation we are taking the difference of the average energy (and molar volumes) in each phase. -+P(-)

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/ramess101/MBAR_GCMC/issues/10#issuecomment-430855701, or mute the thread https://github.com/notifications/unsubscribe-auth/AWUvhKX2Xqmb-W_RtsKYkXCLZz7z3o4Tks5ul-ougaJpZM4XmUtf .

ramess101 commented 5 years ago

@jpotoff

I have confirmed that the P DeltaV contribution is not the problem:

image

The DeltaUv term is clearly where the bias in DeltaHv is coming from:

image

ramess101 commented 5 years ago

@mrshirts @jpotoff

OK, I think I found the problem. < Uv > and < Ul > are the ensemble average internal energies per mol. U in the simulation output is internal energy per N molecules (U_perN). Previously, I was dividing U_perN by N to get internal energy per molecule (and then multiple by Avogadro's number) when determining the ensemble average of U. In other words, I had been using < Uv > = < Uv_perN / Nv >. But if I instead take the ensemble average of U_perN and divide that by the ensemble average of N, i.e., < Uv > = < Uv_perN > / < Nv > I get much better agreement:

image

@jpotoff Do you know if Mick et al. used < U_perN / N > or < UperN > / < N >?

mrshirts commented 5 years ago

Thanks for catching that. Still a good question what the right equation is. One quick question though - if you calculate for different numbers of N, is it constant, or does it vary? If it varies, that's a bit of information to use . . .

ramess101 commented 5 years ago

@mrshirts

Still a good question what the right equation is.

True. I don't know which is more correct.

for different numbers of N, is it constant, or does it vary?

I will look into this

ramess101 commented 5 years ago

Yes does vary with respect to N (recall the curvature in the U vs N histogram plots)