Skip to content

Commit

Permalink
Correct implementation of FlowReactor component name/index
Browse files Browse the repository at this point in the history
  • Loading branch information
speth authored and ischoegl committed Jun 14, 2023
1 parent 7cab377 commit a605695
Show file tree
Hide file tree
Showing 3 changed files with 60 additions and 9 deletions.
8 changes: 5 additions & 3 deletions include/cantera/zeroD/FlowReactor.h
Expand Up @@ -103,12 +103,14 @@ class FlowReactor : public IdealGasReactor
m_max_ss_error_fails = max_fails;
}

//! Return the index in the solution vector for this reactor of the
//! component named *nm*. Possible values for *nm* are "X" (position),
//! "U", the name of a homogeneous phase species, or the name of a surface
//! Return the index in the solution vector for this reactor of the component named
//! *nm*. Possible values for *nm* are "density", "speed", "pressure",
//! "temperature", the name of a homogeneous phase species, or the name of a surface
//! species.
size_t componentIndex(const string& nm) const override;

string componentName(size_t k) override;

void updateSurfaceState(double* y) override;

protected:
Expand Down
40 changes: 34 additions & 6 deletions src/zeroD/FlowReactor.cpp
Expand Up @@ -342,21 +342,49 @@ void FlowReactor::getConstraints(double* constraints) {

size_t FlowReactor::componentIndex(const string& nm) const
{
// check for a gas species name
size_t k = m_thermo->speciesIndex(nm);
size_t k = speciesIndex(nm);
if (k != npos) {
return k + m_offset_Y;
} else if (nm == "rho" || nm == "density") {
} else if (nm == "density") {
return 0;
} else if (nm == "U" || nm == "velocity") {
} else if (nm == "speed") {
return 1;
} else if (nm == "P" || nm == "pressure") {
} else if (nm == "pressure") {
return 2;
} else if (nm == "T" || nm == "temperature") {
} else if (nm == "temperature") {
return 3;
} else {
return npos;
}
}

string FlowReactor::componentName(size_t k)
{
if (k == 0) {
return "density";
} else if (k == 1) {
return "speed";
} else if (k == 2) {
return "pressure";
} else if (k == 3) {
return "temperature";
} else if (k >= 4 && k < neq()) {
k -= 4;
if (k < m_thermo->nSpecies()) {
return m_thermo->speciesName(k);
} else {
k -= m_thermo->nSpecies();
}
for (auto& S : m_surfaces) {
ThermoPhase* th = S->thermo();
if (k < th->nSpecies()) {
return th->speciesName(k);
} else {
k -= th->nSpecies();
}
}
}
throw CanteraError("FlowReactor::componentName", "Index {} is out of bounds.", k);
}

}
21 changes: 21 additions & 0 deletions test/python/test_reactor.py
Expand Up @@ -1447,6 +1447,27 @@ def test_catalytic_surface(self):
assert X == approx([0.01815402, 0.21603645, 0.045640395])
assert r.thermo.density * r.area * r.speed == approx(mdot)

def test_component_names(self):
surf = ct.Interface('methane_pox_on_pt.yaml', 'Pt_surf')
gas = surf.adjacent['gas']
r = ct.FlowReactor(gas)
rsurf = ct.ReactorSurface(surf, r)
sim = ct.ReactorNet([r])
sim.initialize()

assert r.n_vars == 4 + gas.n_species + surf.n_species
assert sim.n_vars == r.n_vars

for i in range(r.n_vars):
name = r.component_name(i)
assert r.component_index(name) == i

with pytest.raises(IndexError, match="No such component: 'spam'"):
r.component_index('spam')

with pytest.raises(ct.CanteraError, match='out of bounds'):
r.component_name(200)


class TestSurfaceKinetics(utilities.CanteraTest):
def make_reactors(self):
Expand Down

0 comments on commit a605695

Please sign in to comment.