Suppose you have developed some useful extensions to Cantera for a particular application. Your code probably consists of a library of C++ classes and functions, and one or more C++ application programs that use the library. Since you are the developer of the library, you are most likely comfortable working in C++, but even so, for solving problems quickly you may wish you could knock out a short Python script or MATLAB m-file, instead of having to write, compile, and link a C++ driver program every time. Or perhaps you are happy to work in C++, but you have shown your work to others, and now they want to use your extensions too, but they are not C++ programmers.
Since Cantera provides bindings for Python, MATLAB, and Fortran, you can create interfaces for your library by following the same procedures used for the rest of Cantera. Here we'll describe step-by-step how to do this.
Before we start, we should say upfront that this process is a bit involved, and is not automated. The main reason for the complexity is that Cantera supports multiple user environments. If we only wanted to generate Python bindings for C++ classes, we could use SWIG (http://www.swig.org), which reads in a C++ class definition, and writes out a corresponding Python "wrapper" class and C++ interfacing code so that the methods of the Python class invoke those of the C++ class.
There are two problems with using SWIG, however. The first is that while SWIG supports multiple output languages, MATLAB and Fortran 90 are not among them. The second problem is that the SWIG-generated Python classes function (by design) very much like the C++ classes. They do not take advantage of things you can do in Python that are not possible in C++.
For example, Cantera classes that represent phases of matter have a method called getMoleFractions that can be used like this:
ThermoPhase* gas = importPhase("gas.cti");
...
int nsp = gas.nSpecies();
double x[nsp];
gas.getMoleFractions(x); Since C++ cannot return arrays as return values, the array x must first be created with an appropriate size, and then the method gtMoleFractions is called to write the mole fraction values into x. Mole fraction values are written to the aray for all species, in the order in which they were declared in the input file. If only selected values are desired, several lines of additional code would be required to find and select the desired species:
int i_H2 = has.speciesIndex("H2");
int i_CO = has.speciesIndex("CO");
int i_CH3 = has.speciesIndex("CH3")
double xh2 = x[i_H2];
double xco = x[i_CO];
double xch3 = x[i_CH3]; But in Python, this can all be done much more compactly:
gas = importPhase("gas.cti")
...
x = gas.moleFractions(["H2", "CO", "CH3"])
In the Python class, there is a method moleFractions that returns an array of mole fractions, so there is no need to manually allocate the array as in C++ - Python takes care of the memory allocation, and deletes it when it is no longer used. The moleFractions method takes an optional argument that is a list of species names, in which case the array returned contains mole fractions for just those species, in the order listed.
This example illustrates why Cantera Python classes are not copies of the underlying C++ class. Instead, they are designed to operate in a way that uses features of Python to enhance usability. This means that the task of writing the interfacing code is more complex, however, since it requires though and judgement, which automated tools are not very good at.
To see how this works, consider the class Droplet, shown below, that models a single evaporating droplet. This class is intentionally simplified for use as an example, but all of the essential features are the same as for more complex classes.
#ifndef DROPLET_H
#define DROPLET_H
// parameters taken from Example 3.2, p. 103, in "An Introduction to
// Combustion" by S.R. Turns
#include <cantera/Cantera.h>
namespace DropletNamespace {
const double Dab_0 = 8.1e-6; // m2/s
const double Tboil = 489.5; //
const double molWt = 170.337;
const double hfg = 2.56e5; // J/kg
const double rho_vapor = 0.4267;
const double rho_liquid = 749.0;
class Droplet {
public:
Droplet(double D=-1.0, double T=300.0) : m_diam(D),
m_temp(T),
m_yinf(0.0) {}
virtual ~Droplet() {}
void setTemperature(double T) { m_temp = T; }
void setYinf(double yinf = 0.0) { m_yinf = yinf; }
double diam() const { return m_diam; }
double mass() const { return rho_liquidPi*m_diam*m_diam*m_diam/6.0}
void evaporate(double time) {
if (m_diam > 0.0) {
double d2 = m_diam*m_diam - evapConstant*time;
if (d2 < 0.0) m_diam = 0.0;
else m_diam = sqrt(d2);
}
}
double evapRate() const {
return 4.0*Pi*rho_vapor*Dab()*log(1.0+B());
}
double B() const {
double ys = surfaceMassFraction();
return (ys - m_yinf)/(1.0 - ys);
}
double surfaceMassFraction() const {
double ps = Psat();
return ps*molWt/(ps*molWt + (OneAtm - ps)*28.014);
}
double Psat() const { return OneAtm*exp(-hfg*molWt*
(1.0/m_temp - 1.0/Tboil)/GasConstant); }
double evapConstant() const { return (8.0*rho_vapor*Dab()/rho_liquid)
*log(1.0 + B());}
double lifetime() const { return m_diam*m_diam/evapConstant();}
protected:
double Dab() { return Dab_0*pow(800.0/399.0, 1.5); }
double m_diam, m_temp, m_yinf;
};
}
#endif
We would like to write bindings for Python, MATLAB, and Fortran,so that we can do thigs like this:
class Spray:
def __init__(self):
self._drops = []
def addDroplet(self, d):
self._drops.append(d)
def evaporate(self, time):
for d in self._drops:
d.evaporate(time)
def mass(self):
m = 0.0
for d in self._drops:
m += d.mass()
def d32(self):
sumd3 = 0.0
sumd2 = 0.0
for drop in self._drops:
d = drop.diam()
sumd3 += d*d*d
sumd2 += d*d
return sumd3/sumd2
dbar = 1.0e-4
sigma = 1.0e-5
cloud = []
for i in range(20):
d = random.gauss(dbar, sigma)
d = Droplet(d, T)
cloud.append(d)
print T, d.lifetime()
1.6.3