gustavochm / sgtpy

Other
38 stars 14 forks source link

SGTPy.fit n-hexane example: RunTimeWarning #3

Closed JiaYuanChng closed 1 month ago

JiaYuanChng commented 2 years ago

Hello Gustavo.

I just ran the n-hexane parameter estimation example provided in the documentation (page 43). It ran fine, I got the same set of values. However there was a RunTimeWarning alert for the ideal Helmholtz computation:

SAFTpkgs/lib/python3.9/site-packages/sgtpy/vrmie_pure/ideal.py:25: RuntimeWarning: invalid value encountered in log a = np.log(rho * broglie_vol**3) - 1

Is there something that could be done to suppress this alert? I am just wondering if this causes any accuracy problems to the computation.

Also, the solver happens to be rather sensitive to the data input type: there was an error if I had input the attractive exponent without the decimal.

ValueError: Integers to negative integer powers are not allowed.

Thank you.

Best wishes, JY

gustavochm commented 2 years ago

Hi JY,

Thanks for bringing this up. This error: SAFTpkgs/lib/python3.9/site-packages/sgtpy/vrmie_pure/ideal.py:25: RuntimeWarning: invalid value encountered in log a = np.log(rho * broglie_vol**3) - 1.

This error usually means that the density solution procedure diverged to an negative value during the parameter estimation. This could be due to a bad initial guess or because awkward fluid parameters during the minimizing procedure (which are more challenging to converge). I would say that as long as the minimize function converges there should not be any problem. However, I would still verify the final set of parameters. Warnings can be suppressed np.warning but I usually rather see if there is something that it is not working properly.

For the error ValueError: Integers to negative integer powers are not allowed.

Thank you very much! I wasn't aware of that issue, I will update the component object to force the conversion of integers to float numbers.

Regards, Gustavo

JiaYuanChng commented 2 years ago

Hi Gustavo,

Thanks for the info.

Also, thank you for producing such a beautiful package. The ubiquitous numpy scipy aside, this is probably the nicest thing I have used so far in engineering computation!

Two more questions before I leave this page in peace:

  1. As mentioned in the documentation, we need to provide an input geometric parameters for ring molecules. The number of segments are however limited to 3, 4, 5, or 7 segments. Is there an empirical equation that computes the geometric parameter as a function of segments (ms)?
  2. I read that the SAFT gamme Mie object is currently limited to the segments/structures in the database. Would there be a new release that enables one to perform parameter estimations on customised segments?

Best wishes, JY

gustavochm commented 2 years ago

Hey JY,

Many thanks for you kind words, they mean a lot to me!

Regardless your questions.

  1. The extra parameter (ring) was fitted for several geometries using molecular simulation data. Please see this article for further details on the parametrization.

  2. I included the recent published parameters from this article. But the implementation is not limited to those groups, if needed you can edit the database. Please refer to this notebook, essentially the database is a pandas dataframe, so you can set new groups and their interactions parameters.

If you don't have access to the articles, please let me know, I can share them with you.

Regards, Gustavo

JiaYuanChng commented 2 years ago

Hi Gustavo,

Many thanks for the info.

As for my second question, perhaps I should have made my question clearer: In the homonuclear SAFT VR Mie, we could use scipy to optimise the epsilon and sigma. As far as I understood from the tutorials, correct me if I am wrong, we would need to supply a whole set of parameters that describes a group in SAFT gamma Mie, which would then enable fluid properties calculation, but currently it does not support parameter estimation/optimisation?

Best wishes, JY

gustavochm commented 2 years ago

Hi JY,

To be honest, I haven't fitted any particular group for SAFT gamma by myself. But for it could be done using scipy and sgtpy. You would need to define an objective function to fit your group, this can include pure component or mixture data.

from sgtpy import database, component, mixture, saftgammamie

def fobj(inc, database):
    # this is an example of the parameters of the new group, you can also set a given value for any of these (for example Sk=1)
    vk, Sk, sigma, eps, lambda_r, eps_kl, lr_kl = inc
    # adding a new group with overtwrite=True
    database.add_group(name='new_group', vk=vk, Sk=Sk, sigma=sigma, eps=eps, lr=lambda_r, overwrite=True)
    # settting interactions between the new group and existing groups
    database.new_interaction_mie('new_group', 'other_group', eps_kl, lr_kl, overwrite=True)

    fluid = component(GC={'new_group':1, 'other_group':1})
    fluid.saftgammamie(database=database)
    eos = saftgammamie(fluid)

    # compute properties from eos object
    ...
    # compute error against experimental data
    ...
   return error

That function can be optimized by scipy. The body of the calculation will depend on the available experimental data.

Hope that helps,

Gustavo

JiaYuanChng commented 2 years ago

I will try this out and bring this to a close should everything go well.

Thanks again!

JiaYuanChng commented 1 month ago

Problem sloved.