Note

Go to the end to download the full example code.

Simulate a combustor using a coupled equivalent reactor network with stream recycling#

This example shows how to set up and solve a “coupled” ERN (equivalent reactor network) in PyChemkin.

An ERN is employed as a reduced-order model to simulate the steady-state combustion process inside a gas turbine combustor chamber. This reduced-order reactor network model retains the complexity of the combustion chemistry by sacrificing details of the combustor geometries, the spatial resolution, and the mass and energy transfer processes. An ERN usually comprises PSRs (perfectly-stirred reactors) and PFRs (plug-flow reactors). The network configuration and connectivity, the reactor parameters, and the mass flow rates can be determined from “hot” steady-state CFD simulation results and/or from observations and measured data of the actual/similar devices. Once a reactor network is calibrated against the experimental data of a gas combustor, it becomes a handy tool for quickly estimating the emissions from the combustor when it is subjected to certain variations in the fuel compositions.

This figure shows a proposed ERN model of a fictional gas turbine combustor.

System Message: WARNING/2 (/home/runner/work/pychemkin/pychemkin/doc/source/examples/reactor_network/PSRnetwork_coupled.rst, line 45)

Cannot scale image! Could not get size from "examples/reactor_network/combustor_ERN.png": [Errno 2] No such file or directory: 'examples/reactor_network/combustor_ERN.png'

The primary fuel is mixed with the incoming air to form a fuel-lean mixture before entering the chamber through the primary inlet. Additional air, the primary air, is introduced to the combustion chamber separately through openings surrounding the primary inlet. The secondary air is entrained into the combustion chamber through well-placed holes on the liners at a location slightly downstream from the primary inlet.

The first PSR (reactor #1) represents the mixing zone around the main injector where the cool premixed fuel-air stream and the primary air stream are preheated by mixing with the hot combustion products from the recirculation zone.

The configurations of a “coupled” ERN (the same as the reactor network when you use the Chemkin GUI) has one major difference in requirements. The “hybrid” ERN in PyChemkin allows all PSRs to have an external outlet. However, for the “coupled” reactor network in PyChemkin and in the Chemkin GUI, ONLY the “last” PSR of the network can have an external outlet. Since the “recirculation zone* does not have any downstream outflow (through flow), it cannot be the “last” reactor of a “coupled” ERN. Consequently, even both the “PSRnetwork” and the current “PSRnetwork_coupled” represent exactly the same ERN, their configurations are different. In the current “coupled” ERN configuration, the recirculation zone becomes PSR #2 and the flame zone becomes the last reactor of the ENR, PSR #3.

Downstream from PSR #1, PSR #3, the flame zone is where the combustion of the heated fuel-air mixture takes place. The secondary air is injected here to cool down the combustion exhaust before it exits the combustion chamber. A portion of the exhaust gas coming out of PSR #3 does not leave the combustion chamber directly and is diverted to PSR #2, the recirculation zone. The external outlet of the ERN must be the direct downstream from this reactor.

The majority of the outlet flow from PSR #2 is recirculated back to the flame zone to sustain the fuel-lean premixed flame there. The rest of the hot gas from PSR #2 travels further back to PSR #1, the mixing zone, to preheat the fuel-lean mixture just entering the combustion chamber. Finally, the cooled flue gas leaves the chamber in a stream-like manner. Typically, a PFR is applied to simulation of the outflow.

The reactors in the “coupled” network are solved simultaneously. Therefore, there is no need to define any tear stream.

Import PyChemkin packages and start the logger =======================================+======

from pathlib import Path
import time

import ansys.chemkin.core as ck  # Chemkin
from ansys.chemkin.core import Color
from ansys.chemkin.core.inlet import Mixture
from ansys.chemkin.core.inlet import Stream  # external gaseous inlet
from ansys.chemkin.core.logger import logger

# Chemkin PSR model (steady-state)
from ansys.chemkin.core.stirreactors.PSR import PSRSetResTimeEnergyConservation as Psr
from ansys.chemkin.core.stirreactors.PSRcluster import PSRCluster as Ern
import numpy as np  # number crunching

# check working directory
current_dir = str(Path.cwd())
logger.debug("working directory: " + current_dir)
# set verbose mode
ck.set_verbose(True)
# set interactive mode for plotting the results
# interactive = True: display plot
# interactive = False: save plot as a PNG file
global interactive
interactive = True

Create a chemistry set#

The mechanism to load is the GRI 3.0 mechanism for methane combustion. This mechanism and its associated data files come with the standard Ansys Chemkin installation in the /reaction/data directory.

# set mechanism directory (the default Chemkin mechanism data directory)
data_dir = Path(ck.ansys_dir) / "reaction" / "data"
mechanism_dir = data_dir
# create a chemistry set based on the GRI mechanism
MyGasMech = ck.Chemistry(label="GRI 3.0")
# set mechanism input files
# including the full file path is recommended
MyGasMech.chemfile = str(mechanism_dir / "grimech30_chem.inp")
MyGasMech.thermfile = str(mechanism_dir / "grimech30_thermo.dat")

Preprocess the gasoline chemistry set#

# preprocess the mechanism files
ierror = MyGasMech.preprocess()

Set up gas mixtures based on the species in this chemistry set#

Create the “fuel” and the “air” mixtures to initialize the external inlet streams. The fuel for this case is pure methane.

# fuel is pure methane
fuel = Mixture(MyGasMech)
fuel.temperature = 650.0  # [K]
fuel.pressure = 10.0 * ck.P_ATM  # [atm] => [dyne/cm2]
fuel.x = [("CH4", 1.0)]

# air is modeled as a mixture of oxygen and nitrogen
air = Mixture(MyGasMech)
air.temperature = 650.0  # [K]
air.pressure = 10.0 * ck.P_ATM
air.x = ck.Air.x()  # mole fractions

Create external inlet streams from the mixtures#

Create the fuel and the air streams/mixtures before setting up the external inlet streams. The fuel in this case is pure methane. The X_by_Equivalence_Ratio() method is used to form the main premixed inlet stream to the mixing zone reactor. The external inlets to the first reactor, the mixing zone, and the second reactor, the flame zone, are simply air mixtures with different mass flow rates and temperatures.

Note

PyChemkin has air redefined as a convenient way to set up the air stream/mixture in the simulations. Use the ansys.chemkin.Air.x() or ansys.chemkin.Air.y() method when the mechanism uses O2 and N2 for oxygen and nitrogen. Use the ansys.chemkin.Air.x('L') or ansys.chemkin.Air.y('L') method when oxygen and nitrogen are represented by o2 and n2.

# primary fuel-air mixture
# products from the complete combustion of the fuel mixture and air
products = ["CO2", "H2O", "N2"]
# species mole fractions of added/inert mixture.
# (You can also create an additives mixture here.)
add_frac = np.zeros(MyGasMech.kk, dtype=np.double)  # no additives: all zeros

# create the unburned fuel-air mixture
premixed = Stream(MyGasMech)
# mean equivalence ratio
equiv = 0.6
ierror = premixed.x_by_equivalence_ratio(
    MyGasMech, fuel.x, air.x, add_frac, products, equivalenceratio=equiv
)
# check fuel-oxidizer mixture creation status
if ierror != 0:
    print("Error: Failed to create the fuel-oxidizer mixture.")
    exit()

# list the composition of the unburned fuel-air mixture
premixed.list_composition(mode="mole")

# set premixed fuel-air inlet temperature and mass flow rate
premixed.temperature = fuel.temperature
premixed.pressure = fuel.pressure
premixed.mass_flowrate = 500.0  # [g/sec]

# primary air stream to mix with the primary fuel-air stream
primary_air = Stream(MyGasMech, label="Primary_Air")
primary_air.x = air.x
primary_air.pressure = air.pressure
primary_air.temperature = air.temperature
primary_air.mass_flowrate = 50.0  # [g/sec]

# secondary bypass air stream
secondary_air = Stream(MyGasMech, label="Secondary_Air")
secondary_air.x = air.x
secondary_air.pressure = air.pressure
secondary_air.temperature = 670.0  # [K]
secondary_air.mass_flowrate = 100.0  # [g/sec]

# find the species index
ch4_index = MyGasMech.get_specindex("CH4")
o2_index = MyGasMech.get_specindex("O2")
no_index = MyGasMech.get_specindex("NO")
co_index = MyGasMech.get_specindex("CO")

Define reactors in the reactor network#

Set up the PSR for each zone one by one with external inlets only. For PSRs, use the set_inlet() method to add the external inlets to the reactor. A PFR always requires one external inlet when it is instantiated. From upstream to downstream, they are premix zone, recirculation zone, and flame zone. The recirculation zone does not have an external inlet.

Note

PyChemkin requires that the first reactor/zone must have at least one external inlet. Because the rest of the reactors have at least the through-flow from the immediate upstream reactor, they do not require an external inlet.

Note

The stream parameter used to instantiate a PSR is used to establish the guessed reactor solution and is modified when the network is solved by the ERN.

# PSR #1: mixing zone
mix = Psr(premixed, label="mixing zone")
# use different guess temperature
mix.set_estimate_conditions(option="TP", guess_temp=800.0)
# set PSR residence time (sec): required for PSRSetResTimeEnergyConservation model
mix.residence_time = 0.5 * 1.0e-3
# add external inlets
mix.set_inlet(premixed)
mix.set_inlet(primary_air)

# PSR #2: recirculation zone
recirculation = Psr(premixed, label="recirculation zone")
# use the equilibrium state of the inlet gas mixture as the guessed solution
recirculation.set_estimate_conditions(option="TP", guess_temp=1600.0)
# set PSR residence time (sec): required for PSRSetResTimeEnergyConservation model
recirculation.residence_time = 1.5 * 1.0e-3

# PSR #3: flame zone
flame = Psr(premixed, label="flame zone")
# use the equilibrium state of the inlet gas mixture as the guessed solution
flame.set_estimate_conditions(option="TP", guess_temp=1600.0)
# set PSR residence time (sec): required for PSRSetResTimeEnergyConservation model
flame.residence_time = 1.5 * 1.0e-3
# add external inlet
flame.set_inlet(secondary_air)

Create the reactor network#

Create a coupled reactor network named PSRcluster with the list of the PSRs in the correct order. For a simple chain network, you do not need to define the connectivity among the reactors. The reactor network model automatically figures out the through-flow connections.

Note

The order of the reactor in the list is important as it dictates the solution sequence and thus the convergence rate.

# instantiate the PSR network as a coupled reactor network
PSR_list = [mix, recirculation, flame]
PSRcluster = Ern(PSR_list, label="combustor_cluster")

Define the PSR connectivity#

Because the current reactor network contains recycling streams, for example, from PSR #3 to PSR #1 and PSR #2, you must explicitly specify the network connectivity. To do this, you use a split dictionary to define the outflow connections among the PSRs in the network. The key of the split dictionary is the label of the originate PSR. The value is a list of tuples consisting of the label of the target PSR and its mass flow rate fraction.

{ "originate PSR" : [("target PSR1", fraction), ("target PSR2", fraction)],
  "another originate PSR" : [("target PSR1", fraction), ... ], ... }

For example, the outflow from reactor PSR1``is split into two streams. 90% of the mass flow rate goes to ``PSR2 and the rest is diverted to PSR5. The split dictionary entry for the outflow split of PSR1 is as follows:

"PSR1" : [("PSR2", 0.9), ("PSR5", 0.1)]

Note

The outlet flow from a reactor that is leaving the reactor network must be labeled as EXIT>> when you define the outflow splitting.

# PSR #1 outlet flow splitting
split_table = [(flame.label, 1.0)]
PSRcluster.set_recycling_stream(mix.label, split_table)
# PSR #2 outlet flow splitting
# PSR #2 does not have an external outlet.
split_table = [(mix.label, 0.15), (flame.label, 0.85)]
PSRcluster.set_recycling_stream(recirculation.label, split_table)
# PSR #3 outlet flow splitting
# part of the outlet flow from PSR #3 exits the reactor network
split_table = [(recirculation.label, 0.2)]
PSRcluster.set_recycling_stream(flame.label, split_table)

Solve the reactor network#

Use the run() method to solve the entire reactor network. The PSRCluster solves the PSRs simultaneously.

# set the start wall time
start_time = time.time()
# solve the reactor network iteratively
status = PSRcluster.run()
if status != 0:
    print(Color.RED + "Failed to solve the reactor network." + Color.END)
    exit()

# compute the total runtime
runtime = time.time() - start_time
print()
print(f"Total simulation duration: {runtime} [sec].")
print()

Postprocess the reactor network solutions#

Since the reactors in the same cluster are solved together, you must get use the process_cluster_solution() method first, and extract the solution stream of a specific PSR outlet by using the get_reactor_stream() method with the (1-based) reactor index. Once you get the solution as a stream, you can use any stream or mixture method to further manipulate the solutions.

# postprocess the solutions of all PSRs in the cluster
iErr = PSRcluster.process_cluster_solution()

# verify the mass flow rate in and out of the PSR cluster
print(
    f"net external inlet mass flow rate = "
    f"{PSRcluster.total_inlet_mass_flow_rate} [g/sec]."
)
print(
    f"net outlet mass flow rate = {PSRcluster.get_cluster_outlet_flowrate()} [g/sec]."
)

# display the reactor solutions
print("=" * 10)
print("reactor/zone")
print("=" * 10)
for index in range(PSRcluster.numb_psrs):
    id = index + 1
    name = PSRcluster.get_reactor_label(id)
    # get the solution stream of the reactor
    sstream = PSRcluster.get_reactor_stream(name)
    print(f"Reactor: {name}.")
    print(f"Temperature = {sstream.temperature} [K].")
    print(f"Mass flow rate = {sstream.mass_flowrate} [g/sec].")
    print(f"CH4 = {sstream.x[ch4_index]}.")
    print(f"O2 = {sstream.x[o2_index]}.")
    print(f"CO = {sstream.x[co_index]}.")
    print(f"NO = {sstream.x[no_index]}.")
    print("-" * 10)

# clean up
ck.done()

Gallery generated by Sphinx-Gallery