andrewlee94
commented
3 years ago Hi @Jackxi1024,
The short answer is that yes you can write a wrapper that will allow Pyomo to call external function, but there are some conditions that need to be met.
The first thing to note is that the ExternalFunction interface you are looking at only supports multiple input-single output type functions - i.e. you need separate function calls for each property of interest. Pyomo does have a new GreyBox model feature that supports multiple input-multiple output type models, but I have no experience working with it yet.
Next, you function needs to return more than just the value of the property at the given state - at a minimum it needs to return the first partial derivatives of each property with respect to each state, and ideally the second partial derivatives as well. Most existing property tools do not support this behaviour, hence why IDAES has created our own tools for this. There are multiple ways to address this, notably analytical differentiation my hand, automatic differentiation tools or numerical methods, but these often involve a bit of work and have different advantages and disadvantages.
Extending from the above however, the solvers we use assume that the functions of the model are twice-continuously differentiable - i.e. they are smooth in both the first and second derivatives. Phase equilibrium phenomena are notably not smooth - phases appear and disappear depending upon the state. When we wrote our external function methods, we put a lot of thought into how to handle this behaviour and smooth these transitions wherever possible. Depending on how the tool you wish to use is structured, this might be a lot of work. I will note that you can sometimes get away with non-smooth functions, but you risk having trouble converging the problem near the phase transition points.
Going back to the example of the Cubic EoS in IDAES - the cubic_roots.so only solves for Z (the roots of the cubic) given the A, B, u and w parameters from the general cubic equation of state (if fact, it is little more than a method for finding the roots of a cubic). The rest of the equations for the cubic EoS are written as Pyomo expressions, including all the departure functions for finding the residual enthalpy, entropy, fugacity, etc.
A better example for what you want can be found in the Helmholtz equation of state code: https://github.com/IDAES/idaes-pse/tree/main/idaes/generic_models/properties/helmholtz. In this case all the calculations for the equation of state are done in the external function, and there is only a thin Python wrapper in place to link this all to Pyomo.
So, to conclude; yes it is possible, but there is a lot of work involved.
Hello, I recently tried to add a new equation of state to idaes. I found a python package to calculate the compress factor, residual enthalpy and fugacity of the equation of state, but I found it difficult for pyomo to call external python functions. Moreover, due to the complexity of the equation of state, I can't write explicit expressions for compressibility factor, residual enthalpy and fugacity (they can be expressed by Constraints of pyomo) Later, I saw the dynamic link library file "cubic_roots.so" in the cubic equation of state file "ceos. py". You use this external file to calculate the compressibility factor of the gas or liquid phase. This makes me very inspired. It seems that pyomo will dynamically call this dynamic link library when solving the model Therefore, I'm wondering whether I can make the external function as a dynamic link library and return it to equation of state of idaes. More in detail, when solving the flash model, I need the fugacity of the gas-liquid components, and I can dynamically get a value by calling the dynamic link library I hope you can provide some suggestions about this DLL. Thank you very much