Integral constant for internal energy of an ionic liquid

[hosein]

Integral constant for internal energy of ionic liquid

I have a question, and I will be really grateful if someone helps me. I have a polynomial equation for internal energy which I calculated by integration an equation of state formula, which is based on density. But, because I calculated this using integration one integration constant which is temperature dependent( based on other articles) that I don't know how can I formulate it to have its magnitude to calculate internal energy at other range of data. My simulation box has contained 200 molecules of ionic liquid with one negative ion( PF6) and a positive one( butyl methyl imidazolium). Because according to internal energy equation at zero density internal energy is equal to the integration constant, we considered it as ionic liquid internal energy at ideal gas state. With all those in mind, how can I use a degree of freedom of rotational, vibrational, and translational to formulate this integration constant dependent of temperature to use it in other range of data? Or, is there any other method to formulate it?

the main equation is this:

([Zth + Zin] - 1)V^2 = e +f/rho+ g*rho^2
in=internal
th=thermal
Z=compressibility factor
(Zth - 1)V^2 = eth +fth/rho+ gth*rho^2
(Zin)V^2 = ein +fin/rho+ gin*rho^2

Ein =∫Pin/rho^2 drho+ F(T) = RT[(ein(T)/2)*rho^2+ fin(T)*rho + (gin(T)/4)*rho^4]+ F(T)
F(T)?

the main equation can be fitted to experimental data, but the Ein cannot.
Actually, I considered the assumption of zero density(ideal gas state for my ionic liquid), but it will show that Ein=F(T), somehow it will help that maybe F(T) can be Ein at the ideal gas state, and if I wanted to formulate it with degree of freedom it will be 3RT( polyatomic ideal gas).

But I calculated Cv and Cp from Ein and when I fitted them( Cv and Cp) to experimental data, their constant(F'(T) or their ideal part contribution which is the derivation of Ein F(T)) shows temperature dependency which means F(T) should have a power 2 or more temperature in its formula to produce a derivation which has temperature dependency.

Actually, the only thing that I can guess is that F(T) could be formulated properly if I calculate the degree of freedom classically for my ionic liquid ([BMIM][PF6]), which I don't know how.

thanks in advance

 

[hosein]

I did not find any formulation of DoF, but I decided to calculate the non-ideal part instead. So, there is no need of the ideal part.

 

[Roark]

It seems like you are doing molecular dynamics simulations? Or is it Monte Carlo?

In any case, let's figure out what you are doing.

First, what is the point of dividing the energy into the internal and thermal components. And an even more serious questions is how do you even subdivide it that way? Internal energy, by definition, contains thermal energy. Do you mean to say the potential energy and the thermal energy instead?

Larger question. How come you are so concerned with the integration constant? Is it really an integration constant after all if you believe that it depends on the temperature? Anyway usually we are concerned in the change of thermodynamic properties between states. Integration constants are usually not important because they always depend on the reference state that we use.

 

[hosein]

Thank you, Roark, for your good questions. I am doing molecular dynamics simulation. Internal energy in my system means total energy because it is thermodynamics system. I need integration constant because I want to report the thermodynamic properties, not their change. I want to explain it more in detail, so you can help me.
Ein=∫Pin/rho^2 d rho +F(T)= RT{ (e*rho^2)+(f*rho)+(g*rho^4)}+F(T)
Cp=R{ (a*rho^2)+(b*rho)+(c*rho^4)}+F'(T)
in these formulas: rho= density. e,f,g,a,b,c = temperature dependant coefficients. Pin= internal pressure. Ein= internal energy=total energy. F(T)= integration constant. F'(T)= temperature derivation of integration constant.

I calculate enthalpy from Ein and then calculate Cp. My main problem is that according to the article that I used, these integration constants, which are y-intercept,( F(T) and F'(T)) should be an ideal contribution of corresponding thermodynamical properties, but in MD and its numerical solutions, these parameters are just y-intercepts. So if I have enough data in MD (from high rho to zero rho) these y-intercepts would be an ideal contribution( because they are in rho=0) and it would be easy for me to have them. In my case, there is no way to calculate all the rho until zero( I work between 4.8-4.9 mol/lit), so I just extrapolate the thermodynamical properties against rho and then find y-intercept which is not ideal contribution ( and not temperature dependent and fluctuational to temperature) because this charts cannot properly extrapolate and find the exact ideal contribution( the thermodynamical properties at zero rho) because every extrapolation in different temperature would find the easiest way to reach y-intercept. So I cannot say that this is an ideal contribution at all, I need these y-intercepts( in my case not ideal contribution) to calculate the properties correctly, but the only way could be some mathematical methods to calculate F(T) and F'(T). Do you know any methods? or is there anything which you can advice me to do?
Thanks in advance

 

[Roark]

I am still unsure of what you are trying to achieve, but my advice is to just treat F(T) and F'(T) as zero. At the zero density limit (empty vacuum) I imagine the thermodynamic internal energy to be equal to zero, and the same with the heat capacity.

Do you have a published article that you are using as a reference? If I can see an example it might help me help you.

 

[hosein]

Dear Roark,
I uploaded the article. As you can see in this article, the author used "3/2 RT" for F(T) because it was monoatomic gas, so in my case, I should consider the F(T) equal to 3RT( it is an ionic liquid fluid which is polyatomic)( right?). But the problem is in my case I don't have the data to zero density like the author of this article had. I mean, his F(T) which is y-intercept accurately became ideal contribution (3/2 RT) because he traced the data until zero density and y-intercept matched the exact amount of F(T). In my case, I used Molecular Dynamics Simulation and because of force field limitation, I cannot have the data to zero density. So consequently, the F(T) in my formula would become a y-intercept which is not equal to really ideal contribution or 3RT because the data( density= 4-5 mol/lit) would be extrapolated to find their corresponding y-intercept, which ideally should be ideal contribution, but here, they are ideal contribution plus an error for extrapolation ( or not having the exact data to zero density). I know, I am a little confused because of so much overanalyzing.
If you were in my place, and in the above equation for Cp you had everything but F'(T) how would you calculate that( considering the article)??
thank you very much for your time and consideration.

 

[Roark]

OK, thanks for the article! I just realized that Z is the compressibility factor. I hope I can give you something useful out of this conversation, since you've done a lot of explaining already.

Basically, you calculated the total energy of the system as a function of temperature. And then you did a polynomial fit with respect to density (rho) and temperature-dependent variables (e,f,g,F). You want to know the zero-density limit value of F, without actually doing the simulation. So you turned to using fundamental relationships (like the degree of freedom relationship to heat capacity) to approximate F.

I am sorry, but I don't think you are going to be able to calculate F (or F') without doing the simulations. Easy fundamental relationships like the DoF to heat capacity is not going to give you the temperature dependence observed in experiments. Even if you could find a way to approximating this value, as a reader I would be skeptical of the consistency of different models (since you are mixing different theories).

So I want to discuss with you exactly why you can't do ideal gas simulations. Is it because the simulation software you are using demands charge neutrality? Also, if you could simulate ideal gas ionic liquid, would it exist as pairs or separate ions? I am thinking perhaps there is a way to do this because I think it will yield better results.

 

[hosein]

Thank you, dear Roark, for your comment,
I can simulate at the ideal gas state( using just a pair of the ion at no periodic boundary condition=zero pressure). But my problem is not that, my problem is that the F'(T) in my calculations is not exactly the ideal contribution. I fitted the experimental data for heat capacity at constant pressure with the Cp formula in the article, But instead of something logically dependent on temperature, I got some fluctuated number for F'(T) as follows:
T-----------F'(T)
298.15 -0.1304
303.15 0.5617
308.15 -2.756
313.15 0.7692
318.15 0.1357
323.15 -1.044
I think there is an ideal contribution in these F'(T) plus some random error created by extrapolation error( better to say caused by lack of enough data at a wide range near to zero density,which are in gas state, which is a limitation for my calculation because of force field accuracy at one phase, my main state was fluid,). But I think this could be possible that I calculate the values at the ideal gas state to prevent this extrapolation error. I mean by adding one point near to zero density, which I calculated at the ideal gas state, to a far range ( 4.8-4.9 mol/lit) and then fitting the data to calculate coefficients. With this approach, I would prevent the extrapolation error which is caused by lack of data in a wide range near zero density point. Am I right? do you think one point is enough to prevent that?
Really thanks for your time

 

[Chestermiller]

I don't see what the problem is. The basic equation you are working with is $$dU=C_vdT+RT\frac{Z_I}{\rho}d\rho$$where $Z_I$ is the internal pressure compressibility factor, given by: $$Z_I=f\rho + e\rho^2+g\rho^4$$If we combine these two equations, we obtain:$$dU=C_vdT+RT(f+e\rho+g\rho^3)d\rho$$If we take as a reference state the gas at ideal gas conditions of temperature $T_{ref}$ and density approaching zero, we can use Hess' Law to express the internal energy at temperature T and density $\rho$ as:$$U(T,\rho)=U(T_{ref},0)+\int_{T_{ref}}^T{C_v^{IG}(T')dT'}+RT\rho\left(f+\frac{e}{2}\rho+\frac{g}{4}\rho^3\right)$$where $C_v^{IG}(T)$ is the ideal gas heat capacity at temperature T. In deriving this relationship, the density integration has been performed as a definite integral from density zero (where ideal gas behavior applies) to density $\rho$. So where does this F(T) come in?

 

[Roark]

So where does this F(T) come in?

I am still wrapping my head around it, but I think F(T) is the ideal gas heat capacity and temperature T. We are talking about finding out what this value really is.

T-----------F'(T)
298.15 -0.1304
303.15 0.5617
308.15 -2.756
313.15 0.7692
318.15 0.1357
323.15 -1.044

Are these fitted values from your non-ideal gas state simulations? I am going to assume that they are. Please let me know if I have misunderstood.

Thank you, dear Roark, for your comment,
I think there is an ideal contribution in these F'(T) plus some random error created by extrapolation error( better to say caused by lack of enough data at a wide range near to zero density,which are in gas state, which is a limitation for my calculation because of force field accuracy at one phase, my main state was fluid,). But I think this could be possible that I calculate the values at the ideal gas state to prevent this extrapolation error. I mean by adding one point near to zero density, which I calculated at the ideal gas state, to a far range ( 4.8-4.9 mol/lit) and then fitting the data to calculate coefficients. With this approach, I would prevent the extrapolation error which is caused by lack of data in a wide range near zero density point. Am I right? do you think one point is enough to prevent that?

In any sort of polynomial fit, extrapolation outside the fitted range have artifacts. Is one point enough? Probably not. If you have 100 points at other densities and 1 point near the zero density limit, it won't do anything to help mitigate the error. If you want to do it this way, I think you would need to do many more simulations at different densities than just one.

I think I see why you are hesitant to use "ideal gas" simulations to get F(T). The zero-density limit is exactly what its names suggests - it is a limit and does not exist in real life. But my opinion is that it will probably be OK to use ideal gas simulations to approximate the real F(T).

I'll illustrate my experiences which will explain why I think it is OK to use ideal gas simulations for F(T). There is a free energy value called the solvation free energy. It is the free energy difference between a solute molecule in the ideal gas state and its free energy in a solvent (like water). Molecular simulations people like you and I will use one molecule in vacuum and one molecule in solution to estimate the solvation free energy all the time. All I am trying to say is that simulated infinite dilution and simulated ideal gases are actually very good approximates for the actual infinite dilution limit and the actual zero density limit thermodynamics.

If I were in your shoes, I would calculate F'(T) for one molecule as a function of temperature. Hopefully this will have less error than your current results.

If you are worried about the force-field, look at how it was parameterized in the first place. Was it thermodynamics based or quantum-chemistry based? If it is quantum-chemistry based, it was most likely parameterized assuming non-interacting molecules, although with ionic liquids it might be different. So your ideal gas simulations should be consistent with the bulk simulations. If it was optimized using thermodynamic values, that is an inconsistency and could be discussed further.

 

[hosein]

I don't see what the problem is. The basic equation you are working with is $$dU=C_vdT+RT\frac{Z_I}{\rho}d\rho$$where $Z_I$ is the internal pressure compressibility factor, given by: $$Z_I=f\rho + e\rho^2+g\rho^4$$If we combine these two equations, we obtain:$$dU=C_vdT+RT(f+e\rho+g\rho^3)d\rho$$If we take as a reference state the gas at ideal gas conditions of temperature $T_{ref}$ and density approaching zero, we can use Hess' Law to express the internal energy at temperature T and density $\rho$ as:$$U(T,\rho)=U(T_{ref},0)+\int_{T_{ref}}^T{C_v^{IG}(T')dT'}+RT\rho\left(f+\frac{e}{2}\rho+\frac{g}{4}\rho^3\right)$$
Dear Chestermiller,

The first two sentence of right side of above equation is F(T) which we considered as ideal contribution of internal ( total) energy and F'(T) is its derivation to temperature.

where $C_v^{IG}(T)$ is the ideal gas heat capacity at temperature T. In deriving this relationship, the density integration has been performed as a definite integral from density zero (where ideal gas behavior applies) to density $\rho$. So where does this F(T) come in?

 

[Chestermiller]

What do e, f, and g look like as functions of T?

 

[hosein]

I am still wrapping my head around it, but I think F(T) is the ideal gas heat capacity and temperature T. We are talking about finding out what this value really is.

Dear Raork,

thank you for your time and consideration. F(T) is ideal gas state internal ( total) energy and F'(T) is ideal gas state heat capacity contribution which is temperature derivation of F(T).

Are these fitted values from your non-ideal gas state simulations? I am going to assume that they are. Please let me know if I have misunderstood.

No, they are fitted values of experimental data, and I chose F'(T) to show you that they don't have any reasonable relation with temperature which I consider this as extrapolation artifact. But using this F'(T) and fitted coefficients I can retrieve experimental Cp with a small error of less than 1. So the equations work even when the F'(T) has no reasonable relation with temperature( because of extrapolation).

In any sort of polynomial fit, extrapolation outside the fitted range have artifacts. Is one point enough? Probably not. If you have 100 points at other densities and 1 point near the zero density limit, it won't do anything to help mitigate the error. If you want to do it this way, I think you would need to do many more simulations at different densities than just one.

I think I see why you are hesitant to use "ideal gas" simulations to get F(T). The zero-density limit is exactly what its names suggests - it is a limit and does not exist in real life. But my opinion is that it will probably be OK to use ideal gas simulations to approximate the real F(T).

I'll illustrate my experiences which will explain why I think it is OK to use ideal gas simulations for F(T). There is a free energy value called the solvation free energy. It is the free energy difference between a solute molecule in the ideal gas state and its free energy in a solvent (like water). Molecular simulations people like you and I will use one molecule in vacuum and one molecule in solution to estimate the solvation free energy all the time. All I am trying to say is that simulated infinite dilution and simulated ideal gases are actually very good approximates for the actual infinite dilution limit and the actual zero density limit thermodynamics.

If I were in your shoes, I would calculate F'(T) for one molecule as a function of temperature. Hopefully this will have less error than your current results.
It was really useful, and I would do it.
If you are worried about the force-field, look at how it was parameterized in the first place. Was it thermodynamics based or quantum-chemistry based? If it is quantum-chemistry based, it was most likely parameterized assuming non-interacting molecules, although with ionic liquids it might be different. So your ideal gas simulations should be consistent with the bulk simulations. If it was optimized using thermodynamic values, that is an inconsistency and could be discussed further.

Thank you very much.

 

[hosein]

What do e, f, and g look like as functions of T?

They are coefficients of the equation of state and temperature dependent, but I should do some calculation to find their relation to T which is the very center of my work. My main problem was finding the ideal contribution values {F(T) for internal energy and F'(T) for heat capacities}. Roark proposed to do simulation at ideal gas state( zero pressure and system temperature), but it is a little risky for ion pair. I want to do it and see what happens.

 

[Chestermiller]

I would have done the fit to the model results very differently. I would have expressed e, f, and g as low order polynomials in T, such as $$e(T)=e_0+e_1(T-310.65)+e_2(T-310.65)^2$$
Similarly for f and g. I would also have expressed F(T) as a low order polynomial in T:
$$F(T)=F_0+F_1(T-310.65)+F_2(T-310.65)^2$$Then I would have had 12 constant parameters to fit to the global set of data. This would have caused much less noise in the fit to F(T).

 

[hosein]

I would have done the fit to the model results very differently. I would have expressed e, f, and g as low order polynomials in T, such as $$e(T)=e_0+e_1(T-310.65)+e_2(T-310.65)^2$$
Similarly for f and g. I would also have expressed F(T) as a low order polynomial in T:
$$F(T)=F_0+F_1(T-310.65)+F_2(T-310.65)^2$$Then I would have had 12 constant parameters to fit to the global set of data. This would have caused much less noise in the fit to F(T).

Dear Chestmiller,
I would normally fit them to ordinary polynomial, and you suggest that when I compute coefficients and simulated ideal gas state for F(T) instead of fitting them to ordinary polynomial I should fit them to these low order polynomial, right? But how did you choose 310.65? please explain in more detail why they are better because I am not familiar with this equation, or at least, give me some website to learn about this concept. I think I can comprehend that they are better, but I have no logical mathematic reason for it. Also, one of the riskiest steps for me is a derivation of all these coefficients to temperature, and if this new approach would cause less noise, it is a great improvement for me.
Thank you very muchThank you very much

 

[Chestermiller]

I would have done the fit to the model results very differently. I would have expressed e, f, and g as low order polynomials in T, such as $$e(T)=e_0+e_1(T-310.65)+e_2(T-310.65)^2$$
Similarly for f and g. I would also have expressed F(T) as a low order polynomial in T:
$$F(T)=F_0+F_1(T-310.65)+F_2(T-310.65)^2$$Then I would have had 12 constant parameters to fit to the global set of data. This would have caused much less noise in the fit to F(T).

Dear Chestmiller,
I would normally fit them to ordinary polynomial, and you suggest that when I compute coefficients and simulated ideal gas state for F(T) instead of fitting them to ordinary polynomial I should fit them to these low order polynomial, right? But how did you choose 310.65? please explain in more detail why they are better because I am not familiar with this equation, or at least, give me some website to learn about this concept. I think I can comprehend that they are better, but I have no logical mathematic reason for it. Also, one of the riskiest steps for me is a derivation of all these coefficients to temperature, and if this new approach would cause less noise, it is a great improvement for me.
Thank you very muchThank you very much

I chose 310.65 because it is right at the middle of the range of temperatures you indicated in post #9. Here's why using a lower order polynomial has the effect of smoothing out the noise. Imagine having 9 data points of y vs x, and you fit a 8th degree polynomial to the data for a function you expect to be relatively smooth in the region of the data. The high order polynomial will pass through all 9 data points, but will be very noisy because of the variability and uncertainty in the data. But, if you find the best least squares fit for a 2nd degree polynomial representation of the data, it will of course be very smooth, and will smooth out the noise. In the fitting that I am recommending for your problem, you will be determining the 12 constant coefficients that give the best least squares fit to your generated internal energy data, for all the data at all the temperatures and densities.

 

[hosein]

I chose 310.65 because it is right at the middle of the range of temperatures you indicated in post #9. Here's why using a lower order polynomial has the effect of smoothing out the noise. Imagine having 9 data points of y vs x, and you fit a 8th degree polynomial to the data for a function you expect to be relatively smooth in the region of the data. The high order polynomial will pass through all 9 data points, but will be very noisy because of the variability and uncertainty in the data. But, if you find the best least squares fit for a 2nd degree polynomial representation of the data, it will of course be very smooth, and will smooth out the noise. In the fitting that I am recommending for your problem, you will be determining the 12 constant coefficients that give the best least squares fit to your generated internal energy data, for all the data at all the temperatures and densities.

Thank you very much, the 310.65 part was really interesting, and my main question is about that. Is it going to affect the real values? I mean, for coefficients and the derivation to temperature, am I going to get real values or they would be different because of this part? I need them to calculate Cp, Cv, E and etc.

Best regards

 

[Chestermiller]

Thank you very much, the 310.65 part was really interesting, and my main question is about that. Is it going to affect the real values? I mean, for coefficients and the derivation to temperature, am I going to get real values or they would be different because of this part? I need them to calculate Cp, Cv, E and etc.

Best regards

No, this choice of reference temperature couldn't affect anything.

 

[hosein]

No, this choice of reference temperature couldn't affect anything.

thank you very much.
I will implement that. That is a great help.
Best regards

 

[hosein]

Dear Roark,
I hope you have a good time. Excuse me for taking your time again. I wanted to do ideal state simulation, but there is something bothering me. I took on ion pair in nonperiodic boundary condition and used fix NVE and a thermostat to calculate F(T) (energy in the ideal gas state). I have a thermodynamic question: as you know the absolute total energy in Molecular Dynamics simulation means nothing(because of the different reference state), but its differences have to mean. Do you think, in my simulation when I take the average of production run total energy, is this can be the answer of F(T)(total energy in the ideal gas state)?