[Python] Make ReactorSurface available from Python

This commit is contained in:
Ray Speth 2016-05-16 22:23:53 -04:00
parent 3f766370b9
commit c3b9a0c821
4 changed files with 142 additions and 66 deletions

View file

@ -467,6 +467,7 @@ cdef extern from "cantera/equil/MultiPhase.h" namespace "Cantera":
cdef extern from "cantera/zeroD/ReactorBase.h" namespace "Cantera":
cdef cppclass CxxWall "Cantera::Wall"
cdef cppclass CxxReactorSurface "Cantera::ReactorSurface"
cdef cppclass CxxFlowDevice "Cantera::FlowDevice"
cdef cppclass CxxReactorBase "Cantera::ReactorBase":
@ -490,6 +491,7 @@ cdef extern from "cantera/zeroD/Reactor.h":
string componentName(size_t) except +
size_t neq()
void getState(double*)
void addSurface(CxxReactorSurface*)
void addSensitivityReaction(size_t) except +
void addSensitivitySpeciesEnthalpy(size_t) except +
@ -530,6 +532,19 @@ cdef extern from "cantera/zeroD/Wall.h":
size_t nSensParams(int)
cdef extern from "cantera/zeroD/ReactorSurface.h":
cdef cppclass CxxReactorSurface "Cantera::ReactorSurface":
CxxReactorSurface()
double area()
void setArea(double)
void setKinetics(CxxKinetics*)
void setCoverages(double*)
void setCoverages(Composition&) except +
void syncCoverages()
void addSensitivityReaction(size_t) except +
size_t nSensParams()
cdef extern from "cantera/zeroD/flowControllers.h":
cdef cppclass CxxFlowDevice "Cantera::FlowDevice":
CxxFlowDevice()
@ -933,6 +948,10 @@ cdef class IdealGasConstPressureReactor(Reactor):
cdef class FlowReactor(Reactor):
pass
cdef class ReactorSurface:
cdef CxxReactorSurface* surface
cdef Kinetics _kinetics
cdef class WallSurface:
cdef CxxWall* cxxwall
cdef object wall

View file

@ -83,12 +83,9 @@ upstream = ct.Reservoir(gas, name='upstream')
# this reservoir is irrelevant.
downstream = ct.Reservoir(gas, name='downstream')
# use a 'Wall' object to implement the reacting surface in the reactor.
# Since walls have to be installed between two reactors/reservoirs, we'll
# install it between the upstream reservoir and the reactor. The area is
# set to the desired catalyst area in the reactor, and surface reactions
# are included only on the side facing the reactor.
w = ct.Wall(upstream, r, A=cat_area, kinetics=[None, surf])
# Add the reacting surface to the reactor. The area is set to the desired
# catalyst area in the reactor.
rsurf = ct.ReactorSurface(surf, r, A=cat_area)
# The mass flow rate into the reactor will be fixed by using a
# MassFlowController object.

View file

@ -398,6 +398,83 @@ cdef class WallSurface:
self.cxxwall.addSensitivityReaction(self.side, m)
cdef class ReactorSurface:
"""
Represents a surface in contact with the contents of a reactor.
:param kin:
The `Kinetics` or `Interface` object representing reactions on this
surface.
:param A:
The area of the reacting surface [m^2]
"""
def __cinit__(self):
self.surface = new CxxReactorSurface()
def __dealloc__(self):
del self.surface
def __init__(self, kin=None, Reactor r=None, *, A=None):
if kin is not None:
self.kinetics = kin
if r is not None:
self.install(r)
if A is not None:
self.area = A
def install(self, Reactor r):
r.reactor.addSurface(self.surface)
property area:
""" Area on which reactions can occur [m^2] """
def __get__(self):
return self.surface.area()
def __set__(self, A):
self.surface.setArea(A)
property kinetics:
"""
The `InterfaceKinetics` object used for calculating reaction rates on
this surface.
"""
def __get__(self):
return self._kinetics
def __set__(self, Kinetics k):
self._kinetics = k
self.surface.setKinetics(self._kinetics.kinetics)
property coverages:
"""
The fraction of sites covered by each surface species.
"""
def __get__(self):
if self._kinetics is None:
raise Exception('No kinetics manager present')
self.surface.syncCoverages()
return self._kinetics.coverages
def __set__(self, coverages):
if self._kinetics is None:
raise Exception("Can't set coverages before assigning kinetics manager.")
if isinstance(coverages, (dict, str, unicode, bytes)):
self.surface.setCoverages(comp_map(coverages))
return
if len(coverages) != self._kinetics.n_species:
raise ValueError('Incorrect number of site coverages specified')
cdef np.ndarray[np.double_t, ndim=1] data = \
np.ascontiguousarray(coverages, dtype=np.double)
self.surface.setCoverages(&data[0])
def add_sensitivity_reaction(self, int m):
"""
Specifies that the sensitivity of the state variables with respect to
reaction *m* should be computed. *m* is the 0-based reaction index.
The Surface must be installed on a reactor and part of a network first.
"""
self.surface.addSensitivityReaction(m)
cdef class Wall:
r"""
A Wall separates two reactors, or a reactor and a reservoir. A wall has a

View file

@ -743,10 +743,12 @@ class TestConstPressureReactor(utilities.CanteraTest):
C = np.zeros(self.interface1.n_species)
C[0] = 0.3
C[4] = 0.7
self.w1.left.kinetics = self.interface1
self.w2.left.kinetics = self.interface2
self.w1.left.coverages = C
self.w2.left.coverages = C
self.surf1 = ct.ReactorSurface(self.interface1, A=0.2)
self.surf2 = ct.ReactorSurface(self.interface2, A=0.2)
self.surf1.coverages = C
self.surf2.coverages = C
self.surf1.install(self.r1)
self.surf2.install(self.r2)
self.net1 = ct.ReactorNet([self.r1])
self.net2 = ct.ReactorNet([self.r2])
@ -783,8 +785,7 @@ class TestConstPressureReactor(utilities.CanteraTest):
self.assertNear(self.r1.T, self.r2.T, rtol=1e-5)
self.assertNear(self.r1.thermo.P, self.r2.thermo.P, rtol=1e-6)
if surf:
self.assertArrayNear(self.w1.left.coverages,
self.w2.left.coverages,
self.assertArrayNear(self.surf1.coverages, self.surf2.coverages,
rtol=1e-4, atol=1e-8)
def test_closed(self):
@ -854,7 +855,7 @@ class TestFlowReactor(utilities.CanteraTest):
self.assertNear(r.speed, v, 1e-3)
class TestWallKinetics(utilities.CanteraTest):
class TestSurfaceKinetics(utilities.CanteraTest):
def make_reactors(self):
self.net = ct.ReactorNet()
@ -870,34 +871,24 @@ class TestWallKinetics(utilities.CanteraTest):
self.r2.volume = 0.01
self.net.add_reactor(self.r2)
self.w = ct.Wall(self.r1, self.r2)
self.w.area = 1.0
def test_coverages(self):
self.make_reactors()
self.w.left.kinetics = self.interface
surf1 = ct.ReactorSurface(self.interface, self.r1)
self.w.left.coverages = {'c6HH':0.3, 'c6HM':0.7}
self.assertNear(self.w.left.coverages[0], 0.3)
self.assertNear(self.w.left.coverages[1], 0.0)
self.assertNear(self.w.left.coverages[4], 0.7)
surf1.coverages = {'c6HH':0.3, 'c6HM':0.7}
self.assertNear(surf1.coverages[0], 0.3)
self.assertNear(surf1.coverages[1], 0.0)
self.assertNear(surf1.coverages[4], 0.7)
self.net.advance(1e-5)
C_left = self.w.left.coverages
self.assertEqual(self.w.right.kinetics, None)
with self.assertRaises(Exception):
self.w.right.coverages
C_left = surf1.coverages
self.make_reactors()
self.w.right.kinetics = self.interface
self.w.right.coverages = 'c6HH:0.3, c6HM:0.7'
self.assertNear(self.w.right.coverages[0], 0.3)
self.assertNear(self.w.right.coverages[4], 0.7)
self.assertEqual(self.w.left.kinetics, None)
with self.assertRaises(Exception):
self.w.left.coverages
surf2 = ct.ReactorSurface(self.interface, self.r2)
surf2.coverages = 'c6HH:0.3, c6HM:0.7'
self.assertNear(surf2.coverages[0], 0.3)
self.assertNear(surf2.coverages[4], 0.7)
self.net.advance(1e-5)
C_right = self.w.right.coverages
C_right = surf2.coverages
self.assertNear(sum(C_left), 1.0)
self.assertArrayNear(C_left, C_right)
@ -907,14 +898,14 @@ class TestWallKinetics(utilities.CanteraTest):
self.make_reactors()
self.r1.energy_enabled = False
self.r2.energy_enabled = False
self.w.left.kinetics = self.interface
surf1 = ct.ReactorSurface(self.interface, self.r1)
C = np.zeros(self.interface.n_species)
C[0] = 0.3
C[4] = 0.7
self.w.left.coverages = C
self.assertArrayNear(self.w.left.coverages, C)
surf1.coverages = C
self.assertArrayNear(surf1.coverages, C)
data = []
test_file = 'test_coverages_regression1.csv'
reference_file = '../data/WallKinetics-coverages-regression1.csv'
@ -922,7 +913,7 @@ class TestWallKinetics(utilities.CanteraTest):
for t in np.linspace(1e-6, 1e-3):
self.net.advance(t)
data.append([t, self.r1.T, self.r1.thermo.P, self.r1.mass] +
list(self.r1.thermo.X) + list(self.w.left.coverages))
list(self.r1.thermo.X) + list(surf1.coverages))
np.savetxt(test_file, data, delimiter=',')
bad = utilities.compareProfiles(reference_file, test_file,
@ -932,14 +923,14 @@ class TestWallKinetics(utilities.CanteraTest):
def test_coverages_regression2(self):
# Test with energy equation enabled
self.make_reactors()
self.w.left.kinetics = self.interface
surf = ct.ReactorSurface(self.interface, self.r1)
C = np.zeros(self.interface.n_species)
C[0] = 0.3
C[4] = 0.7
self.w.left.coverages = C
self.assertArrayNear(self.w.left.coverages, C)
surf.coverages = C
self.assertArrayNear(surf.coverages, C)
data = []
test_file = 'test_coverages_regression2.csv'
reference_file = '../data/WallKinetics-coverages-regression2.csv'
@ -947,7 +938,7 @@ class TestWallKinetics(utilities.CanteraTest):
for t in np.linspace(1e-6, 1e-3):
self.net.advance(t)
data.append([t, self.r1.T, self.r1.thermo.P, self.r1.mass] +
list(self.r1.thermo.X) + list(self.w.left.coverages))
list(self.r1.thermo.X) + list(surf.coverages))
np.savetxt(test_file, data, delimiter=',')
bad = utilities.compareProfiles(reference_file, test_file,
@ -992,16 +983,14 @@ class TestReactorSensitivities(utilities.CanteraTest):
r2 = ct.IdealGasReactor(gas2)
net.add_reactor(r2)
w = ct.Wall(r1, r2)
w.area = 1.5
w.left.kinetics = interface
surf = ct.ReactorSurface(interface, r1, A=1.5)
C = np.zeros(interface.n_species)
C[0] = 0.3
C[4] = 0.7
w.left.coverages = C
w.left.add_sensitivity_reaction(2)
surf.coverages = C
surf.add_sensitivity_reaction(2)
r2.add_sensitivity_reaction(18)
for T in (901, 905, 910, 950, 1500):
@ -1155,17 +1144,11 @@ class TestReactorSensitivities(utilities.CanteraTest):
rB = ct.IdealGasReactor(gas2)
if order % 2 == 0:
wA = ct.Wall(rA, rB)
wB = ct.Wall(rB, rA)
surfX = ct.ReactorSurface(interface, rA, A=0.1)
surfY = ct.ReactorSurface(interface, rA, A=10)
else:
wB = ct.Wall(rB, rA)
wA = ct.Wall(rA, rB)
wA.left.kinetics = interface
wB.right.kinetics = interface
wA.area = 0.1
wB.area = 10
surfY = ct.ReactorSurface(interface, rA, A=10)
surfX = ct.ReactorSurface(interface, rA, A=0.1)
C1 = np.zeros(interface.n_species)
C2 = np.zeros(interface.n_species)
@ -1174,8 +1157,8 @@ class TestReactorSensitivities(utilities.CanteraTest):
C2[0] = 0.9
C2[4] = 0.1
wA.left.coverages = C1
wB.right.coverages = C2
surfX.coverages = C1
surfY.coverages = C2
if order // 2 == 0:
net.add_reactor(rA)
@ -1184,7 +1167,7 @@ class TestReactorSensitivities(utilities.CanteraTest):
net.add_reactor(rB)
net.add_reactor(rA)
return rA,rB,wA,wB,net
return rA, rB, surfX, surfY, net
def integrate(r, net):
net.advance(1e-4)
@ -1193,16 +1176,16 @@ class TestReactorSensitivities(utilities.CanteraTest):
S = []
for order in range(4):
rA,rB,wA,wB,net = setup(order)
for (obj,k) in [(rB,2), (rB,18), (wA.left,2),
(wA.left,0), (wB.right,2)]:
rA, rB, surfX, surfY, net = setup(order)
for (obj,k) in [(rB,2), (rB,18), (surfX,2),
(surfY,0), (surfY,2)]:
obj.add_sensitivity_reaction(k)
integrate(rB, net)
S.append(net.sensitivities())
rA,rB,wA,wB,net = setup(order)
for (obj,k) in [(wB.right,2), (wA.left,2), (rB,18),
(wA.left,0), (rB,2)]:
rA, rB, surfX, surfY, net = setup(order)
for (obj,k) in [(surfY,2), (surfX,2), (rB,18),
(surfX,0), (rB,2)]:
obj.add_sensitivity_reaction(k)
integrate(rB, net)