Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add missing array properties to SolutionArray and FlameBase #1342

Merged
merged 4 commits into from
Jul 29, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
25 changes: 20 additions & 5 deletions interfaces/cython/cantera/composite.py
Original file line number Diff line number Diff line change
Expand Up @@ -463,7 +463,7 @@ class SolutionArray:

_scalar = [
# From ThermoPhase
'mean_molecular_weight', 'P', 'T', 'density', 'density_mass',
'mean_molecular_weight', 'P', 'T', 'Te', 'density', 'density_mass',
'density_mole', 'v', 'volume_mass', 'volume_mole', 'u',
'int_energy_mole', 'int_energy_mass', 'h', 'enthalpy_mole',
'enthalpy_mass', 's', 'entropy_mole', 'entropy_mass', 'g', 'gibbs_mole',
Expand All @@ -484,18 +484,25 @@ class SolutionArray:
'chemical_potentials', 'electrochemical_potentials', 'partial_molar_cp',
'partial_molar_volumes', 'standard_enthalpies_RT',
'standard_entropies_R', 'standard_int_energies_RT', 'standard_gibbs_RT',
'standard_cp_R',
'standard_cp_R', 'activities', 'activity_coefficients',
# From Transport
'mix_diff_coeffs', 'mix_diff_coeffs_mass', 'mix_diff_coeffs_mole',
'thermal_diff_coeffs'
'thermal_diff_coeffs', 'mobilities', 'species_viscosities',
]

# From Kinetics (differs from Solution.n_species for Interface phases)
_n_total_species = [
'creation_rates', 'destruction_rates', 'net_production_rates',
'creation_rates_ddC', 'creation_rates_ddP', 'creation_rates_ddT',
'destruction_rates_ddC', 'destruction_rates_ddP', 'destruction_rates_ddT',
'net_production_rates_ddC', 'net_production_rates_ddP',
'net_production_rates_ddT'
]

_n_species2 = ['multi_diff_coeffs', 'binary_diff_coeffs']
_n_species2 = [
'multi_diff_coeffs', 'binary_diff_coeffs', 'creation_rates_ddX',
'destruction_rates_ddX', 'net_production_rates_ddX'
]

_n_reactions = [
'forward_rates_of_progress', 'reverse_rates_of_progress',
Expand All @@ -504,6 +511,13 @@ class SolutionArray:
'delta_enthalpy', 'delta_gibbs', 'delta_entropy',
'delta_standard_enthalpy', 'delta_standard_gibbs',
'delta_standard_entropy', 'heat_production_rates',
'forward_rate_constants_ddC', 'forward_rate_constants_ddP',
'forward_rate_constants_ddT', 'forward_rates_of_progress_ddC',
'forward_rates_of_progress_ddP', 'forward_rates_of_progress_ddT',
'net_rates_of_progress_ddC', 'net_rates_of_progress_ddP',
'net_rates_of_progress_ddT', 'reverse_rates_of_progress_ddC',
'reverse_rates_of_progress_ddP', 'reverse_rates_of_progress_ddP',
'reverse_rates_of_progress_ddT', 'third_body_concentrations',
]
_call_scalar = ['elemental_mass_fraction', 'elemental_mole_fraction']

Expand All @@ -519,7 +533,8 @@ class SolutionArray:
'kinetics_species_index', 'reaction', 'reactions', 'modify_reaction',
'multiplier', 'set_multiplier', 'reaction_equations',
'reactant_stoich_coeff', 'product_stoich_coeff',
'reactant_stoich_coeffs', 'product_stoich_coeffs',
'reactant_stoich_coeffs', 'product_stoich_coeffs', 'product_stoich_coeffs3',
'reactant_stoich_coeffs3', 'product_stoich_coeffs_reversible',
# from Transport
'transport_model',
]
Expand Down
33 changes: 29 additions & 4 deletions interfaces/cython/cantera/onedim.py
Original file line number Diff line number Diff line change
Expand Up @@ -711,7 +711,8 @@ def getter(self):
'entropy_mass', 'g', 'gibbs_mole', 'gibbs_mass', 'cv',
'cv_mole', 'cv_mass', 'cp', 'cp_mole', 'cp_mass',
'isothermal_compressibility', 'thermal_expansion_coeff',
'viscosity', 'thermal_conductivity', 'heat_release_rate']:
'viscosity', 'thermal_conductivity', 'heat_release_rate',
'mean_molecular_weight']:
setattr(FlameBase, _attr, _array_property(_attr))
FlameBase.volume = _array_property('v') # avoid confusion with velocity gradient 'V'
FlameBase.int_energy = _array_property('u') # avoid collision with velocity 'u'
Expand All @@ -723,8 +724,13 @@ def getter(self):
'partial_molar_volumes', 'standard_enthalpies_RT',
'standard_entropies_R', 'standard_int_energies_RT',
'standard_gibbs_RT', 'standard_cp_R', 'creation_rates',
'destruction_rates', 'net_production_rates', 'mix_diff_coeffs',
'mix_diff_coeffs_mass', 'mix_diff_coeffs_mole', 'thermal_diff_coeffs']:
'destruction_rates', 'net_production_rates', 'creation_rates_ddC',
'creation_rates_ddP', 'creation_rates_ddT', 'destruction_rates_ddC',
'destruction_rates_ddP', 'destruction_rates_ddT',
'net_production_rates_ddC', 'net_production_rates_ddP',
'net_production_rates_ddT', 'mix_diff_coeffs', 'mix_diff_coeffs_mass',
'mix_diff_coeffs_mole', 'thermal_diff_coeffs', 'activities',
'activity_coefficients', 'mobilities', 'species_viscosities']:
setattr(FlameBase, _attr, _array_property(_attr, 'n_species'))

# Remove misleading examples and references to setters that don't exist
Expand All @@ -738,7 +744,14 @@ def getter(self):
'equilibrium_constants', 'forward_rate_constants', 'reverse_rate_constants',
'delta_enthalpy', 'delta_gibbs', 'delta_entropy',
'delta_standard_enthalpy', 'delta_standard_gibbs',
'delta_standard_entropy', 'heat_production_rates']:
'delta_standard_entropy', 'heat_production_rates',
'third_body_concentrations', 'forward_rate_constants_ddC',
'forward_rate_constants_ddP', 'forward_rate_constants_ddT',
'forward_rates_of_progress_ddC', 'forward_rates_of_progress_ddP',
'forward_rates_of_progress_ddT', 'net_rates_of_progress_ddC',
'net_rates_of_progress_ddP', 'net_rates_of_progress_ddT',
'reverse_rates_of_progress_ddC', 'reverse_rates_of_progress_ddP',
'reverse_rates_of_progress_ddT']:
setattr(FlameBase, _attr, _array_property(_attr, 'n_reactions'))


Expand Down Expand Up @@ -1436,6 +1449,18 @@ def mixture_fraction(self, m):
vals[i] = self.gas.mixture_fraction(Yf, Yo, 'mass', m)
return vals

@property
def equivalence_ratio(self):
Yf = [self.solution(k, 0) for k in self.gas.species_names]
Yo = [self.solution(k, self.flame.n_points-1) for k in self.gas.species_names]

vals = np.empty(self.flame.n_points)
for i in range(self.flame.n_points):
self.set_gas_state(i)
vals[i] = self.gas.equivalence_ratio(Yf, Yo, "mass")
return vals


class ImpingingJet(FlameBase):
"""An axisymmetric flow impinging on a surface at normal incidence."""
__slots__ = ('inlet', 'flame', 'surface')
Expand Down
2 changes: 1 addition & 1 deletion interfaces/cython/cantera/test/test_convert.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,7 @@ def convert(self, inputFile, thermo=None, transport=None,
output.unlink()
ck2yaml.convert_mech(inputFile, thermo_file=thermo,
transport_file=transport, surface_file=surface, out_name=output,
extra_file=extra, **kwargs)
extra_file=extra, quiet=True, **kwargs)
return output

def checkConversion(self, refFile, testFile):
Expand Down
44 changes: 35 additions & 9 deletions interfaces/cython/cantera/test/test_onedim.py
Original file line number Diff line number Diff line change
Expand Up @@ -580,15 +580,35 @@ def test_save_restore_yaml(self):

def test_array_properties(self):
self.create_sim(ct.one_atm, 300, 'H2:1.1, O2:1, AR:5')

for attr in ct.FlameBase.__dict__:
if isinstance(ct.FlameBase.__dict__[attr], property):
if attr in ['u', 'V']:
msg = "Replaced by property"
with self.assertWarnsRegex(DeprecationWarning, msg):
getattr(self.sim, attr)
else:
getattr(self.sim, attr)
grid_shape = self.sim.grid.shape

skip = {
# Skipped because they are constant, irrelevant, or otherwise not desired
"P", "Te", "atomic_weights", "charges", "electric_potential", "max_temp",
"min_temp", "molecular_weights", "product_stoich_coeffs",
"product_stoich_coeffs3", "product_stoich_coeffs_reversible",
"reactant_stoich_coeffs", "reactant_stoich_coeffs3", "reference_pressure",
"state", "u", "v",
# Skipped because they are 2D (conversion not implemented)
"binary_diff_coeffs", "creation_rates_ddX", "destruction_rates_ddX",
"forward_rates_of_progress_ddX", "net_production_rates_ddX",
"net_rates_of_progress_ddX", "reverse_rates_of_progress_ddX"
}

for attr in dir(self.gas):
if attr.startswith("_") or attr in skip:
continue

try:
soln_value = getattr(self.gas, attr)
except (ct.CanteraError, ct.ThermoModelMethodError):
continue

if not isinstance(soln_value, (float, np.ndarray)):
continue

flame_value = getattr(self.sim, attr)
assert flame_value.shape == np.asarray(soln_value).shape + grid_shape

@utilities.slow_test
def test_save_restore_add_species_yaml(self):
Expand Down Expand Up @@ -999,6 +1019,12 @@ def test_mixture_fraction(self):
self.assertTrue(all(Z >= 0))
self.assertTrue(all(Z <= 1.0))

def test_equivalence_ratio(self):
self.create_sim(p=ct.one_atm)
phi = self.sim.equivalence_ratio
assert phi[0] == np.inf
assert np.isclose(phi[-1], 0.0)

class TestCounterflowPremixedFlame(utilities.CanteraTest):
# Note: to re-create the reference file:
# (1) set PYTHONPATH to build/python.
Expand Down
27 changes: 27 additions & 0 deletions interfaces/cython/cantera/test/test_thermo.py
Original file line number Diff line number Diff line change
Expand Up @@ -1891,6 +1891,33 @@ def test_properties_ndim(self):
self.assertArrayNear(self.gas.reverse_rates_of_progress, ropr[i,j,k])
self.assertArrayNear(self.gas.mix_diff_coeffs, Dkm[i,j,k])

def test_array_properties_exist(self):
grid_shape = (7, 3)
states = ct.SolutionArray(self.gas, grid_shape)

skip = {
# Skipped because they are complicated (conversion not implemented)
"forward_rates_of_progress_ddX", "net_rates_of_progress_ddX",
"reverse_rates_of_progress_ddX", "state"
}
skip.update(ct.SolutionArray._passthrough)

for attr in dir(self.gas):
if attr.startswith("_") or attr in skip:
continue

try:
soln_value = getattr(self.gas, attr)
except (ct.CanteraError, ct.ThermoModelMethodError):
continue

if not isinstance(soln_value, (float, np.ndarray)):
continue

assert hasattr(states, attr), attr
array_value = getattr(states, attr)
assert array_value.shape == grid_shape + np.asarray(soln_value).shape, attr

def test_slicing_onedim(self):
states = ct.SolutionArray(self.gas, 5)
states.TPX = np.linspace(500, 1000, 5), 2e5, 'H2:0.5, O2:0.4'
Expand Down