API for Mixtures

[ibell]

Hey @jowr, the mixtures are getting to a point that I think it would be logical to think about the API. I have most of the core code working (though there are still plenty of problems to resolve). I got the phase envelopes working (was a stupid bug on my side).

I guess the major question mark for me is how we should handle the DLL wrapping. For the languages that (can) wrap C++ directly (Python, Octave, C#, etc), it's pretty straightforward. All we do is wrap the Mixtures class using SWIG or xdress (for python), which internally handles all the other things in a seamless way.

When we have to pass through the DLL, its more complicated. I think for mixtures we need something like

set_mix_components(const char * fluid_listing) // fluid_listing could be "Methane|Propane"

to set the components, and then something like

T = PropsM('T','P',100,'Q',1,z) 

where z is vector of mole fractions of the components, a std::vector of [0.5, 0.5] or something that can be cast to a std::vector<double>

The obvious limitation of this API is that we can only consider one mixture at a time, and if you had two different mixtures you would have to constantly switch between the components of the two mixtures. But you will always have to provide the mixture mole fractions.

One of the tricky thing is what to do about the phase envelopes, as they only apply to one bulk mole fraction vector. I think that what we really need to do is calculate the phase envelope for all possible mole fractions and then interpolate within them? It's pretty nasty. I'm thinking about an absorption system where the bulk mole fractions are changing throughout the system. It's a different story if you are talking about mixtures in ORC or heat pumps, where the bulk composition is pretty much constant apart from constant pressure phase change.

What do you think?

[jowr]

I have two immediate comments and some more general ideas

• Setting components and their fractions independently might result in inconsistency. We have to make sure to understand the way in which functions get called in the shared library. Just imagine a Modelica simulation with two different mixtures causing calls to the library object with alternating components. What about something like T = PropsM('T','P',100,'Q',1,n,s,z) with 1<=n=s.size()=z.size()+1<=nMax with std::vector<std::string> s and std::vector<double> z. I like the concept of reduced mole fraction information z[n]=1-sum(z[0:n-1]).
• Regarding the name, I would prefer something MProps. This allows us to follow the same naming scheme for unit conversions using PropsU and MPropsU.

Since changes in the bulk mole fractions have almost the same impact as changes in the mixture components (new critical values, new phase envelope, new mole weight, ...), I might make sense to set these two things together in one function. The penalty for parsing some extra strings should not be too bad, or do you disagree?

In general, I would like to restructure some of the code to be more object-oriented and easier to debug. I think it is time for a major clean-up. It is very obvious that CoolProp grew over some years and such an "organic" structure is not easy to maintain as I just noticed when dealing with the CPState, FluidClass and Cache objects. If you don't mind, I would like to think about this topic for a while and come with a suggestion. I attached a first draft of my ideas...

[ibell]

I really like your new idea for the structure of CoolProp into a more object-oriented form. The only major problem I see is that it is going to involve a lot of refactoring of the CoolProp codebase. In principle this ought to not cause major problems and will help for long-term maintainability. You are definitely right that the code has been developed very "organically" - many bad design decisions have been made along the way. I think the worst one was to not use SI units everywhere. Unfortunately I don't see a safe way to make SI units default. Maybe we do that in 5.0, but only if we can safely warn the user.

My proposal is that we finally release 4.0, and then after that we do a major refactor of CoolProp as you propose and also handle the mixtures as part of the refactor.

At the DLL level, I think your solution could work, but we can't pass STL containers, so everything has to be done with double* and char*. The parsing of strings is essentially free compared with doing anything with mixtures. This basically means that we can't do phase envelopes which are critical for being able to approach the critical region. MProps is fine with me.

Essentially the only places that use DLL right now are Excel, MathCAD and Labview (and FORTRAN?). I have been reading about wrapping C++ classes in MATLAB (http://www.mathworks.com/matlabcentral/fileexchange/38964-example-matlab-class-wrapper-for-a-c++-class) and while it is not a very automatic or nice wrapping process, it should be possible. So maybe we just say that anyone that uses the DLL will have to suck it up and take a performance penalty while everyone else who uses advanced programming environments can have multiple mixture instances.

[jowr]

Yes, let's go for a release. I will put my changes in the release notes and we can keep everyone excited about 5.0... I think it is rather difficult to change the default units. However, we can start promoting PropsSI as the new primary function that handles both properties and derivatives.

And yes again, I am aware of the pointer-only concept (double* and char*). That is the reason why I suggested to include additional information regarding the size in parameter n. By the way, I think that Fortran can use the object files directly, but I am not sure. (I'll ask Brian and Leonardo about the wrapper).

[ibell]

Hey @jowr, one late night thought on this topic

What do we do for the case that we have two different mixtures in the system, each with the ability to have different component mole fractions. Suppose in this case you have ammonia/water and R32/R125. As you alternate calls for mixtures, how should we deal with this? Suppose the order of calls is:

1. R32[0.5]/R125[0.5]
2. NH3[0.1]/Water[0.9]
3. R32[0.3]/R125[0.7]
4. NH3[0.9]/Water[0.1]

Does this mean we generate 4 phase envelopes? You are right that we have different reducing temperatures, envelopes etc. if we change the mixture mole fractions.

[jowr]

I am not sure how to handle this. One idea could be to make an list of phase envelopes for a mixture with n components and vary the composition with m predefined step-width.

import numpy as np
# Taken from: http://stackoverflow.com/questions/1208118/using-numpy-to-build-an-array-of-all-combinations-of-two-arrays
def cartesian(arrays, out=None):
    arrays = [np.asarray(x) for x in arrays]
    dtype = arrays[0].dtype

    n = np.prod([x.size for x in arrays])
    if out is None:
        out = np.zeros([n, len(arrays)], dtype=dtype)
    m = n / arrays[0].size
    out[:,0] = np.repeat(arrays[0], m)
    if arrays[1:]:
        cartesian(arrays[1:], out=out[0:m,1:])
        for j in xrange(1, arrays[0].size):
            out[j*m:(j+1)*m,1:] = out[0:m,1:]
    return out

def makeCombinations(numCompo=0,numEntries=0):
    stepWidth  = 1.0 / (entries-1.0)
    theRange   = [x*stepWidth for x in range(numEntries)]
    allRanges  = [theRange for i in range(numCompo)]
    allEntries = cartesian(allRanges)
    result = []
    for i in range(len(allEntries)):
        if np.sum(allEntries[i])==1:
            result.append(allEntries[i])
    return np.array(result)

compone   = 3 # components
entries   = 5 # concentration points

print makeCombinations(compone,entries)

This would give you a list concentrations covering the whole range of possible concentrations:

[[ 0.    0.    1.  ]
 [ 0.    0.25  0.75]
 [ 0.    0.5   0.5 ]
 [ 0.    0.75  0.25]
 [ 0.    1.    0.  ]
 [ 0.25  0.    0.75]
 [ 0.25  0.25  0.5 ]
 [ 0.25  0.5   0.25]
 [ 0.25  0.75  0.  ]
 [ 0.5   0.    0.5 ]
 [ 0.5   0.25  0.25]
 [ 0.5   0.5   0.  ]
 [ 0.75  0.    0.25]
 [ 0.75  0.25  0.  ]
 [ 1.    0.    0.  ]]

Such a matrix could then be used to obtain guess value for a calculation with the actual composition. The advantage is that we speed up the phase envelope calculation without loosing flexibility in terms of composition. The only thing I worry about is the computational efficiency. How costly is it to generate such a matrix and the associated phase envelopes? Maybe we could leave it up to the user to decide whether such a matrix is generated or not.

The above implies that we define mixtures by their components. Your example would consist of 2 different mixtures.

[ibell]

That is an intriguing way of handling the combinations. Check your inbox for an alternative idea for how to handle the phase envelopes for binary mixture at least. There are a lot of tricky issues to work out here. Another idea is to build a really accurate phase envelope for a lot of mole fraction combinations and save it to file, sort of like how we handle the TTSE tables. But the ND interpolation of the phase envelopes is a nasty problem.

It seems that my phase envelopes that I thought were working are not entirely working once I pass the critical point, and I see the problem but I don't understand at all how to remedy it. The basic issue is that once you pass the critical point I don't understand how you switch between considering the vapor and liquid phases if you fix beta. On second thought, maybe you just need a flag to the the x_and_y_from_K function (https://github.com/ibell/coolprop/blob/master/CoolProp/Mixtures.cpp#L1123) that gives the normal or inverted (x=y, y=x) behavior? To be continued.

[jowr]

Done in v5