Pickle serialization of Solution objects

[ischoegl]

Pickle support for Solution objects would facilitate multiprocessing applications, where each thread/process is working with a separate copy of a Solution object. The approach is comparable to constructing the object from an input file plus setting phase and state vector.

Closes Cantera/enhancements#10

The PR illustrates a rudimentary implementation (XML/CTI support only), which should not be considered for adoption (at least not in its present form). … now replaced by a permanent solution using YAML serialization

It allows for the following scenario (obviously, multiprocessing applications are more relevant):

import cantera as ct
import pickle

gas = ct.Solution('h2o2.cti', 'ohmech-multi')
gas.TPX = 500, 500000, 'H2:.75, O2:.25'
pickle.dump(gas, open('gas.pkl', 'wb'))

Recreating the object in a separate session via

import pickle
gas = pickle.load(open('gas.pkl','rb'))
gas()

yields:

 ohmech-multi:

       temperature             500  K
          pressure          500000  Pa
           density         1.14397  kg/m^3
  mean mol. weight          9.5115  amu

                          1 kg            1 kmol
                       -----------      ------------
          enthalpy      6.2484e+05        5.943e+06     J
   internal energy      1.8776e+05        1.786e+06     J
           entropy           16391        1.559e+05     J/K
    Gibbs function     -7.5708e+06       -7.201e+07     J
 heat capacity c_p          3127.2        2.974e+04     J/K
 heat capacity c_v            2253        2.143e+04     J/K

                           X                 Y          Chem. Pot. / RT
                     -------------     ------------     ------------
                H2           0.75         0.158965         -14.8055
                O2           0.25         0.841035           -24.87
     [   +7 minor]              0                0

The current implementation for pickling/unpickling makes use of the static XML tree that is retained after creation of a Solution from CTI/XML files. A more useful (but significantly more involved) implementation would generate YAML strings dynamically.

Changes proposed in this pull request:

• N/A: at least at this point, the PR serves for illustration purposes only

Known issues:

• The current implementation will fail if the Solution object is changed after instantiation (e.g. species are added, etc.)
• Unpickling requires that an empty constructor, i.e. _SolutionBase() is allowed. While this works fine for pickle.load, an attempt for gas = Solution() results in a segmentation fault. ... fixed

[speth]

As I mentioned in the thread on #570, implementation of YAML-based serialization is planned, which will work for objects regardless of how they were created and how they have been manipulated.

With the introduction of the YAML format, deprecation and removal of the CTI and XML formats is also planned, so I wouldn't recommend putting any significant effort into features based around these formats. Beyond that, I'm not sure whether there's a reason to keep this PR open.

[ischoegl]

Agreed. I simply wanted to supply an example for pickling (which was surprisingly simple as the XML tree was available - 40ish lines of code!). As I wrote

The PR illustrates a rudimentary implementation (XML/CTI support only), which should not be considered for adoption (at least not in its present form).

I’d still like to keep this PR open until the core of the code can be switched over to YAML. Ultimately, I hope that this feature can be implemented. See also #696, which i believe is a necessary step in that direction.

[ischoegl]

Linking #984, which will create the basis for a permanent implementation.

[ischoegl]

🎉 ... unsurprisingly, thanks to @speth's recent addition of #984, this became even simpler.

PS: Enabling pickling of SolutionArray looks quite feasible, but I'd rather wait for comments at the moment.

[tsikes]

This would be a wonderful addition. I ran into the problem of not being able to do this when I first wanted to use multiprocessing to speed up calculations for optimization. I know that others have had the same issue based on the Cantera google group. The current method of using Cantera in a multiprocessing application is rather unintuitive and finicky. This PR would greatly simplify this task for other users.

[ischoegl]

Kicking this back to draft - the implementation depends on the ability to create a Solution that is unconfigured, so it can be overwritten by __setstate__. What is implemented thus far works fine for the pickle case, but deactivates errors for calls of Solution() without arguments. Having the _SolutionBase in a safe state (i.e. no Segfaults) will take some more work.

PS: I probably won't work on this unless I hear back from one of the @Cantera/committers, i.e. comments are explicitly welcome.

[speth]

I think it should be possible to instantiate an object of the ThermoPhase base class if no other information is available to instantiate the phase, like we do with the Kinetics and Transport objects. In that case, you should just get a phase with no species where most of the methods raise NotImplementedError, which would probably be fine.

I just tried using this in a modified version of the multiprocessing example, and it seems to work better than I had initially expected. If you use itertools.repeat to pass the Solution object as part of the iterable used by Pool.starmap for example, the default chunking behavior results in the Solution being pickled and unpickled 4 times for each worker process, which is not as bad as I had feared (I thought it might end up doing it once for each iteration). So given the easier user experience, I'll withdraw my earlier objection (Cantera/enhancements#10) to implementing this feature at all.

I'm glad the implementation needed here for standalone Solution objects is pretty straightforward. I'm not sure what's going to happen when you try to serialize a an Interface phase where you have to correctly link the neighboring phases.

[ischoegl]

I think it should be possible to instantiate an object of the ThermoPhase base class if no other information is available to instantiate the phase, like we do with the Kinetics and Transport objects.

@speth Could you clarify? None of these are currently instantiable without arguments. E.g. (using main)

In [3]: kin = ct.Kinetics()
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-3-6f9f3674fb0c> in <module>()
----> 1 kin = ct.Kinetics()

/src/interfaces/cython/cantera/base.pyx in cantera._cantera._SolutionBase.__cinit__()
     64             self._init_parts(thermo, species, kinetics, adjacent, reactions)
     65         else:
---> 66             raise ValueError("Arguments are insufficient to define a phase")
     67 
     68         # Initialization of transport is deferred to Transport.__init__

ValueError: Arguments are insufficient to define a phase

As far as I can tell, all of these are limitations of the implementation of the Python interface, and not necessarily about the underlying C++ structure. The crux here is that half of the time, the required Cxx* references remain uninitialized, which leads to segfaults left and right.

[speth]

Sorry for the confusion. What I meant was setting the thermo attribute of a _SolutionBase object to point to an instance of the base class that has been instantiated, e.g. the equivalent of new ThermoPhase(). This would be similar to how we never leave the _SolutionBase.kinetics or _SolutionBase.transport members as null pointers after initialization (at least for objects that derive from the corresponding Python classes), but have them point to objects of these base classes.

Probably the easiest way to do this would be to add this as a case for ThermoFactory, the same way that using KineticsFactory to create an object with the model none just instantiates the base class, at which point the existing _init_parts method would just about work (although I think we'd still need to require species input except when the thermo model is none).

[ischoegl]

@speth ... interesting indeed. Could you post the code you used? (Edit: I may have overanalyzed the 'serialized' here. If it's just accessing the restored object I can probably take it from here)

[speth]

Yes, that's all I did:

gas0 = ct.Solution('gri30.yaml')
with open('test.pkl', 'wb') as out:
    pickle.dump(gas0, out)

with open('test.pkl', 'rb') as infile:
    gas1 = pickle.load(infile)

print(gas1.thermo_model)
print(gas1.kinetics_model)
print(gas1.transport_model)  # segfaults here

[ischoegl]

@speth ... fixed. The bug was hiding in plain sight ... one of the derived class constructor that Python calls (Transport.__init__) usually takes care of setting the C++ transport object. That this doesn't happen here is obvious in hindsight ... it's probably a good idea to consolidate all those calls in base.pyx once CTI/XML are gone.

I also updated the unit test to ensure that the model is properly set (including "None", which is not the default).

[speth]

Thanks, @ischoegl. I'm glad to see this working based on YAML serialization. I had just a few comments on the implementation.

[ischoegl]

@speth … thank you for the careful review! I believe I took care of everything. I ended up disabling Kinetics() and Transport() as I don’t see a current use case (it was surprisingly hard to catch these due to the __slots__). I also clarified some of the edge cases you pointed out.

[speth]

Thanks for the updates. This looks good to me.