[Cython] Translated some samples to use the new API

interfaces/cython/cantera/examples/multiphase/adiabatic.py

@@ -0,0 +1,91 @@
+"""
+Adiabatic flame temperature and equilibrium composition for a fuel/air mixture
+as a function of equivalence ratio, including formation of solid carbon.
+"""
+
+import cantera as ct
+import numpy as np
+import sys
+import csv
+
+##############################################################################
+# Edit these parameters to change the initial temperature, the pressure, and
+# the phases in the mixture.
+
+T = 300.0
+P = 101325.0
+
+# phases
+gas = ct.Solution('gri30.xml')
+carbon = ct.Solution('graphite.xml')
+
+# the phases that will be included in the calculation, and their initial moles
+mix_phases = [(gas, 1.0), (carbon, 0.0)]
+
+# gaseous fuel species
+fuel_species = 'CH4'
+
+# air composition
+air_N2_O2_ratio = 3.76
+
+# equivalence ratio range
+phi_min = 0.3
+phi_max = 3.5
+npoints = 50
+
+##############################################################################
+
+mix = ct.Mixture(mix_phases)
+
+# create some arrays to hold the data
+phi = np.zeros(npoints)
+tad = np.zeros(npoints)
+xeq = np.zeros((mix.n_species,npoints))
+
+# find fuel, nitrogen, and oxygen indices
+ifuel = gas.species_index(fuel_species)
+io2 = gas.species_index('O2')
+in2 = gas.species_index('N2')
+
+if gas.n_atoms(fuel_species,'O') > 0 or gas.n_atoms(fuel_species,'N') > 0:
+    raise "Error: only hydrocarbon fuels are supported."
+
+stoich_O2 = gas.n_atoms(fuel_species,'C') + 0.25*gas.n_atoms(fuel_species,'H')
+
+for i in range(npoints):
+    phi[i] = phi_min + (phi_max - phi_min)*i/(npoints - 1)
+    X = np.zeros(gas.n_species)
+    X[ifuel] = phi[i]
+    X[io2] = stoich_O2
+    X[in2] = stoich_O2*air_N2_O2_ratio
+
+    # set the gas state
+    gas.TPX = T, P, X
+
+    # create a mixture of 1 mole of gas, and 0 moles of solid carbon.
+    mix = ct.Mixture(mix_phases)
+    mix.T = T
+    mix.P = P
+
+    # equilibrate the mixture adiabatically at constant P
+    mix.equilibrate('HP', solver='gibbs', max_steps=1000)
+
+    tad[i] = mix.T
+    print('At phi = {:12.4g}, Tad = {:12.4g}'.format(phi[i], tad[i]))
+    xeq[:,i] = mix.species_moles
+
+# write output CSV file for importing into Excel
+csv_file = 'adiabatic.csv'
+with open(csv_file, 'w') as outfile:
+    writer = csv.writer(outfile)
+    writer.writerow(['phi','T (K)'] + mix.species_names)
+    for i in range(npoints):
+        writer.writerow([phi[i], tad[i]] + list(xeq[:,i]))
+print('Output written to {}'.format(csv_file))
+
+if '--plot' in sys.argv:
+    import matplotlib.pyplot as plt
+    plt.plot(phi, tad)
+    plt.xlabel('Equivalence ratio')
+    plt.ylabel('Adiabatic flame temperature [K]')
+    plt.show()

interfaces/cython/cantera/examples/onedim/diffusion_flame.py

@@ -6,6 +6,7 @@ import cantera as ct
 import numpy as np
 import csv
 
+# Input parameters
 p = ct.one_atm  # pressure
 tin_f = 300.0  # fuel inlet temperature
 tin_o = 300.0  # oxidizer inlet temperature
@@ -15,6 +16,8 @@ mdot_f = 0.24  # kg/m^2/s
 comp_o = 'O2:0.21, N2:0.78, AR:0.01'  # air composition
 comp_f = 'C2H6:1'  # fuel composition
 
+# Distance between inlets is 2 cm.
+# Start with an evenly-spaced 6-point grid.
 initial_grid = np.linspace(0, 0.02, 6)
 
 tol_ss = [1.0e-5, 1.0e-12]  # [rtol, atol] for steady-state problem
@@ -23,10 +26,17 @@ tol_ts = [5.0e-4, 1.0e-9]  # [rtol, atol] for time stepping
 loglevel = 1  # amount of diagnostic output (0 to 5)
 refine_grid = 1   # 1 to enable refinement, 0 to disable
 
+# Create the gas object used to evaluate all thermodynamic, kinetic, and
+# transport properties.
 gas = ct.Solution('gri30.xml', 'gri30_mix')
 gas.TP = gas.T, p
 
+# Create an object representing the counterflow flame configuration,
+# which consists of a fuel inlet on the left, the flow in the middle,
+# and the oxidizer inlet on the right.
 f = ct.CounterflowDiffusionFlame(gas, initial_grid)
+
+# Set the state of the two inlets
 f.fuel_inlet.mdot = mdot_f
 f.fuel_inlet.X = comp_f
 f.fuel_inlet.T = tin_f
@@ -35,19 +45,27 @@ f.oxidizer_inlet.mdot = mdot_o
 f.oxidizer_inlet.X = comp_o
 f.oxidizer_inlet.T = tin_o
 
+# Set error tolerances
 f.flame.set_steady_tolerances(default=tol_ss)
 f.flame.set_transient_tolerances(default=tol_ts)
 
+# construct the initial solution estimate. To do so, it is necessary
+# to specify the fuel species. If a fuel mixture is being used,
+# specify a representative species here for the purpose of
+# constructing an initial guess.
 f.set_initial_guess(fuel='C2H6')
+
+# First disable the energy equation and solve the problem without
+# refining the grid
 f.energy_enabled = False
 f.solve(loglevel, refine_grid=False)
 
+# Now specify grid refinement criteria, turn on the energy equation,
+# and solve the problem again.
 f.energy_enabled = True
 f.set_refine_criteria(ratio=4, slope=0.2, curve=0.3, prune=0.04)
 f.solve(loglevel, refine_grid=refine_grid)
 f.show_solution()
-
-
 f.save('c2h6_diffusion.xml')
 
 z = f.flame.grid
@@ -55,6 +73,7 @@ T = f.T
 u = f.u
 V = f.V
 
+# write the velocity, temperature, and mole fractions to a CSV file
 with open('c2h6_diffusion.csv', 'w') as csvfile:
     writer = csv.writer(csvfile)
     writer.writerow(['z (m)', 'u (m/s)', 'V (1/s)', 'T (K)', 'rho (kg/m3)'] +
@@ -64,3 +83,5 @@ with open('c2h6_diffusion.csv', 'w') as csvfile:
         writer.writerow([z[n], u[n], V[n], T[n], gas.density] + list(gas.X))
 
 print('solution saved to c2h6_diffusion.csv')
+
+f.show_stats(0)

interfaces/cython/cantera/examples/thermo/critical_properties.py

@@ -0,0 +1,27 @@
+"""
+Print the critical state properties for the fluids for which Cantera has
+built-in liquid/vapor equations of state.
+"""
+
+import cantera as ct
+
+fluids = {'water': ct.Water(),
+          'nitrogen': ct.Nitrogen(),
+          'methane': ct.Methane(),
+          'hydrogen': ct.Hydrogen(),
+          'oxygen': ct.Oxygen(),
+          'carbon dioxide': ct.CarbonDioxide(),
+          'heptane': ct.Heptane(),
+          'hfc134a': ct.Hfc134a()
+          }
+
+print('Critical State Properties')
+print('%20s  %10s  %10s  %10s' % ('Fluid','Tc [K]', 'Pc [Pa]', 'Zc'))
+for name in fluids:
+    f = fluids[name]
+    tc = f.critical_temperature
+    pc = f.critical_pressure
+    rc = f.critical_density
+    mw = f.mean_molecular_weight
+    zc = pc * mw / (rc * ct.gas_constant * tc)
+    print('%20s   %10.4g   %10.4G  %10.4G' % (name, tc, pc, zc))

interfaces/cython/cantera/examples/thermo/isentropic.py

@@ -0,0 +1,71 @@
+import cantera as ct
+import math
+import numpy as np
+
+
+def soundspeed(gas):
+    """The speed of sound. Assumes an ideal gas."""
+
+    gamma = gas.cp / gas.cv
+    return math.sqrt(gamma * ct.gas_constant
+                     * gas.T / gas.mean_molecular_weight)
+
+
+def isentropic(gas=None):
+    """
+    ISENTROPIC  isentropic, adiabatic flow example
+
+    In this example, the area ratio vs. Mach number curve is computed. If a
+    gas object is supplied, it will be used for the calculations, with the
+    stagnation state given by the input gas state. Otherwise, the calculations
+    will be done for a 10:1 hydrogen/nitrogen mixture with stagnation T0 =
+    1200 K, P0 = 10 atm.
+
+    """
+
+    if gas is None:
+        gas = ct.Solution('gri30.xml')
+        gas.TPX = 1200.0, 10.0*ct.one_atm, 'H2:1,N2:0.1'
+
+    # get the stagnation state parameters
+    s0 = gas.s
+    h0 = gas.h
+    p0 = gas.P
+
+    mdot = 1  # arbitrary
+    amin = 1.e14
+
+    data = np.zeros((200,4))
+
+    # compute values for a range of pressure ratios
+    for r in range(200):
+
+        p = p0*(r+1)/201.0
+        # set the state using (p,s0)
+        gas.SP = s0, p
+
+        v2 = 2.0*(h0 - gas.h)      # h + V^2/2 = h0
+        v = math.sqrt(v2)
+        area = mdot/(gas.density*v)    # rho*v*A = constant
+        amin = min(amin, area)
+        data[r,:] = [area, v/soundspeed(gas), gas.T, p/p0]
+
+    data[:,0] /= amin
+
+    return data
+
+
+if __name__ == "__main__":
+    print(isentropic.__doc__)
+    data = isentropic()
+    try:
+        import matplotlib.pyplot as plt
+        plt.plot(data[:,1], data[:,0])
+        plt.ylabel('Area Ratio')
+        plt.xlabel('Mach Number')
+        plt.title('Isentropic Flow: Area Ratio vs. Mach Number')
+        plt.show()
+
+    except ImportError:
+        print('area ratio,   Mach number,   temperature,   pressure ratio')
+        print(data)

interfaces/cython/cantera/examples/thermo/rankine.py

@@ -0,0 +1,75 @@
+"""
+A Rankine vapor power cycle
+"""
+
+import cantera as ct
+
+# parameters
+eta_pump = 0.6     # pump isentropic efficiency
+eta_turbine = 0.8  # turbine isentropic efficiency
+p_max = 8.0e5       # maximum pressure
+
+
+def pump(fluid, p_final, eta):
+    """Adiabatically pump a fluid to pressure p_final, using
+    a pump with isentropic efficiency eta."""
+    h0 = fluid.h
+    s0 = fluid.s
+    fluid.SP = s0, p_final
+    h1s = fluid.h
+    isentropic_work = h1s - h0
+    actual_work = isentropic_work / eta
+    h1 = h0 + actual_work
+    fluid.HP = h1, p_final
+    return actual_work
+
+
+def expand(fluid, p_final, eta):
+    """Adiabatically expand a fluid to pressure p_final, using
+    a turbine with isentropic efficiency eta."""
+    h0 = fluid.h
+    s0 = fluid.s
+    fluid.SP =s0, p_final
+    h1s = fluid.h
+    isentropic_work = h0 - h1s
+    actual_work = isentropic_work * eta
+    h1 = h0 - actual_work
+    fluid.HP = h1, p_final
+    return actual_work
+
+
+def printState(n, fluid):
+    print('\n***************** State {} ******************'.format(n))
+    print(fluid.report())
+
+
+if __name__ == '__main__':
+    # create an object representing water
+    w = ct.Water()
+
+    # start with saturated liquid water at 300 K
+    w.TX = 300.0, 0.0
+    h1 = w.h
+    p1 = w.P
+    printState(1, w)
+
+    # pump it adiabatically to p_max
+    pump_work = pump(w, p_max, eta_pump)
+    h2 = w.h
+    printState(2, w)
+
+    # heat it at constant pressure until it reaches the saturated vapor state
+    # at this pressure
+    w.PX = p_max, 1.0
+    h3 = w.h
+    heat_added = h3 - h2
+    printState(3, w)
+
+    # expand back to p1
+    turbine_work = expand(w, p1, eta_turbine)
+    printState(4, w)
+
+    # efficiency
+    eff = (turbine_work - pump_work)/heat_added
+
+    print('efficiency = ', eff)

interfaces/cython/cantera/examples/thermo/sound_speed.py

@@ -0,0 +1,52 @@
+import cantera as ct
+import math
+
+
+def equilSoundSpeeds(gas, rtol=1.0e-6, maxiter=5000):
+    """
+    Returns a tuple containing the equilibrium and frozen sound speeds for a
+    gas with an equilibrium composition.  The gas is first set to an
+    equilibrium state at the temperature and pressure of the gas, since
+    otherwise the equilibrium sound speed is not defined.
+    """
+
+    # set the gas to equilibrium at its current T and P
+    gas.equilibrate('TP', rtol=rtol, maxiter=maxiter)
+
+    # save properties
+    s0 = gas.s
+    p0 = gas.P
+    r0 = gas.density
+
+    # perturb the pressure
+    p1 = p0*1.0001
+
+    # set the gas to a state with the same entropy and composition but
+    # the perturbed pressure
+    gas.SP = s0, p1
+
+    # frozen sound speed
+    afrozen = math.sqrt((p1 - p0)/(gas.density - r0))
+
+    # now equilibrate the gas holding S and P constant
+    gas.equilibrate('SP', rtol=rtol, maxiter=maxiter)
+
+    # equilibrium sound speed
+    aequil = math.sqrt((p1 - p0)/(gas.density - r0))
+
+    # compute the frozen sound speed using the ideal gas expression as a check
+    gamma = gas.cp/gas.cv
+    afrozen2 = math.sqrt(gamma * ct.gas_constant * gas.T /
+                         gas.mean_molecular_weight)
+
+    return aequil, afrozen, afrozen2
+
+
+# test program
+if __name__ == "__main__":
+    gas = ct.Solution('gri30.xml')
+    gas.X = 'CH4:1.00, O2:2.0, N2:7.52'
+    for n in range(27):
+        T = 300.0 + 100.0 * n
+        gas.TP = T, ct.one_atm
+        print(T, equilSoundSpeeds(gas))