*** empty log message ***
This commit is contained in:
Dave Goodwin
parent
53ceee2aa3
commit
3132e39e1f
6 changed files with 335 additions and 10 deletions
|
|
@ -237,8 +237,13 @@ class Kinetics:
|
|||
def multiplier(self,i):
|
||||
return _cantera.kin_multiplier(self.ckin,i)
|
||||
|
||||
def setMultiplier(self,i,v):
|
||||
return _cantera.kin_setMultiplier(self.ckin,i,v)
|
||||
def setMultiplier(self, value = 0.0, reaction = -1):
|
||||
if reaction < 0:
|
||||
nr = self.nReactions()
|
||||
for i in range(nr):
|
||||
_cantera.kin_setMultiplier(self.ckin,i,value)
|
||||
else:
|
||||
_cantera.kin_setMultiplier(self.ckin,reaction,value)
|
||||
|
||||
def advanceCoverages(self,dt):
|
||||
return _cantera.kin_advanceCoverages(self.ckin,dt)
|
||||
|
|
|
|||
|
|
@ -54,8 +54,10 @@ class BurnerFlame(Stack):
|
|||
Stack.solve(self, loglevel = loglevel, refine_grid = refine_grid)
|
||||
|
||||
|
||||
def setRefineCriteria(self, ratio = 10.0, slope = 0.8, curve = 0.8, prune = 0.0):
|
||||
Stack.setRefineCriteria(self, domain = self.flame, ratio = ratio, slope = slope, curve = curve,
|
||||
def setRefineCriteria(self, ratio = 10.0, slope = 0.8,
|
||||
curve = 0.8, prune = 0.0):
|
||||
Stack.setRefineCriteria(self, domain = self.flame,
|
||||
ratio = ratio, slope = slope, curve = curve,
|
||||
prune = prune)
|
||||
|
||||
def setProfile(self, component, locs, vals):
|
||||
|
|
@ -89,7 +91,7 @@ class BurnerFlame(Stack):
|
|||
for n in range(nsp):
|
||||
nm = self.gas.speciesName(n)
|
||||
y[n] = self.solution(nm, j)
|
||||
self.gas.setState_TPY(self.T(j), self.flame.pressure(), y)
|
||||
self.gas.setState_TPY(self.T(j), self.pressure, y)
|
||||
|
||||
|
||||
|
||||
196
Cantera/python/Cantera/OneD/CounterFlame.py
Normal file
196
Cantera/python/Cantera/OneD/CounterFlame.py
Normal file
|
|
@ -0,0 +1,196 @@
|
|||
"""A counterflow flame."""
|
||||
|
||||
from onedim import *
|
||||
import Numeric
|
||||
import math
|
||||
|
||||
def erfc(x):
|
||||
"""The complementary error function."""
|
||||
exp = math.exp
|
||||
|
||||
p = 0.3275911
|
||||
a1 = 0.254829592
|
||||
a2 = -0.284496736
|
||||
a3 = 1.421413741
|
||||
a4 = -1.453152027
|
||||
a5 = 1.061405429
|
||||
|
||||
t = 1.0 / (1.0 + p*x)
|
||||
erfcx = ( (a1 + (a2 + (a3 +
|
||||
(a4 + a5*t)*t)*t)*t)*t ) * exp(-x*x)
|
||||
return erfcx
|
||||
|
||||
def erf(x):
|
||||
"""The error function."""
|
||||
if x < 0:
|
||||
return -(1.0 - erfc(-x))
|
||||
else:
|
||||
return 1.0 - erfc(x)
|
||||
|
||||
|
||||
class CounterFlame(Stack):
|
||||
"""A non-premixed counterflow flame."""
|
||||
|
||||
def __init__(self, gas = None, grid = None):
|
||||
self.fuel_inlet = Inlet('fuel inlet')
|
||||
self.oxidizer_inlet = Inlet('oxidizer inlet')
|
||||
self.gas = gas
|
||||
self.fuel_inlet.set(temperature = gas.temperature())
|
||||
self.oxidizer_inlet.set(temperature = gas.temperature())
|
||||
self.pressure = gas.pressure()
|
||||
self.flame = AxisymmetricFlow('flame',gas = gas)
|
||||
self.flame.setupGrid(grid)
|
||||
Stack.__init__(self, [self.fuel_inlet, self.flame,
|
||||
self.oxidizer_inlet])
|
||||
self.setRefineCriteria()
|
||||
self._initialized = 0
|
||||
|
||||
|
||||
def init(self, fuel = '', oxidizer = 'O2', stoich = -1.0):
|
||||
"""Set the initial guess for the solution.
|
||||
|
||||
The initial guess is generated by assuming infinitely-fast
|
||||
chemistry."""
|
||||
|
||||
gas = self.gas
|
||||
nsp = gas.nSpecies()
|
||||
wt = gas.molecularWeights()
|
||||
|
||||
# find the fuel and oxidizer species
|
||||
iox = gas.speciesIndex(oxidizer)
|
||||
ifuel = gas.speciesIndex(fuel)
|
||||
|
||||
# if no stoichiometric ratio was input, compute it
|
||||
if stoich < 0.0:
|
||||
if oxidizer == 'O2':
|
||||
nh = gas.nAtoms(fuel, 'H')
|
||||
nc = gas.nAtoms(fuel, 'C')
|
||||
stoich = 1.0*nc + 0.25*nh
|
||||
else:
|
||||
raise CanteraError('oxidizer/fuel stoichiometric ratio must'+
|
||||
' be specified, since the oxidizer is not O2')
|
||||
|
||||
s = stoich*wt[iox]/wt[ifuel]
|
||||
y0f = self.fuel_inlet.massFraction(ifuel)
|
||||
y0ox = self.oxidizer_inlet.massFraction(iox)
|
||||
phi = s*y0f/y0ox
|
||||
zst = 1.0/(1.0 + phi)
|
||||
|
||||
yin_f = Numeric.zeros(nsp, 'd')
|
||||
yin_o = Numeric.zeros(nsp, 'd')
|
||||
yst = Numeric.zeros(nsp, 'd')
|
||||
for k in range(nsp):
|
||||
yin_f[k] = self.fuel_inlet.massFraction(k)
|
||||
yin_o[k] = self.oxidizer_inlet.massFraction(k)
|
||||
yst[k] = zst*yin_f[k] + (1.0 - zst)*yin_o[k]
|
||||
|
||||
gas.setState_TPY(self.fuel_inlet.temperature(), self.pressure, yin_f)
|
||||
mdotf = self.fuel_inlet.mdot()
|
||||
u0f = mdotf/gas.density()
|
||||
t0f = self.fuel_inlet.temperature()
|
||||
|
||||
gas.setState_TPY(self.oxidizer_inlet.temperature(),
|
||||
self.pressure, yin_o)
|
||||
mdoto = self.oxidizer_inlet.mdot()
|
||||
u0o = mdoto/gas.density()
|
||||
t0o = self.oxidizer_inlet.temperature()
|
||||
|
||||
# get adiabatic flame temperature and composition
|
||||
tbar = 0.5*(t0o + t0f)
|
||||
gas.setState_TPY(tbar, self.pressure, yst)
|
||||
gas.equilibrate('HP')
|
||||
teq = gas.temperature()
|
||||
yeq = gas.massFractions()
|
||||
|
||||
# estimate strain rate
|
||||
zz = self.flame.grid()
|
||||
dz = zz[-1] - zz[0]
|
||||
a = (u0o + u0f)/dz
|
||||
diff = gas.mixDiffCoeffs()
|
||||
f = math.sqrt(a/(2.0*diff[iox]))
|
||||
|
||||
x0 = mdotf*dz/(mdotf + mdoto)
|
||||
nz = len(zz)
|
||||
|
||||
y = Numeric.zeros([nz,nsp],'d')
|
||||
t = Numeric.zeros(nz,'d')
|
||||
for j in range(nz):
|
||||
x = zz[j]
|
||||
zeta = f*(x - x0)
|
||||
zmix = 0.5*(1.0 - erf(zeta))
|
||||
if zmix > zst:
|
||||
for k in range(nsp):
|
||||
y[j,k] = yeq[k] + (zmix - zst)*(yin_f[k]
|
||||
- yeq[k])/(1.0 - zst)
|
||||
t[j] = teq + (t0f - teq)*(zmix - zst)/(1.0 - zst)
|
||||
print teq, t[j], t0f, zmix, zst
|
||||
else:
|
||||
for k in range(nsp):
|
||||
y[j,k] = yin_o[k] + zmix*(yeq[k] - yin_o[k])/zst
|
||||
t[j] = t0o + (teq - t0o)*zmix/zst
|
||||
|
||||
zrel = zz/dz
|
||||
self.setProfile('u', [0.0, 1.0], [u0f, -u0o])
|
||||
self.setProfile('V', [0.0, x0/dz, 1.0], [0.0, a, 0.0])
|
||||
self.setProfile('T', zrel, t)
|
||||
for k in range(nsp):
|
||||
self.setProfile(gas.speciesName(k), zrel, y[:,k])
|
||||
|
||||
self._initialized = 1
|
||||
|
||||
|
||||
def solve(self, loglevel = 1, refine_grid = 1):
|
||||
if not self._initialized: self.init()
|
||||
Stack.solve(self, loglevel = loglevel, refine_grid = refine_grid)
|
||||
|
||||
|
||||
def setRefineCriteria(self, ratio = 10.0, slope = 0.8, curve = 0.8,
|
||||
prune = 0.0):
|
||||
Stack.setRefineCriteria(self, domain = self.flame,
|
||||
ratio = ratio, slope = slope, curve = curve,
|
||||
prune = prune)
|
||||
|
||||
def setProfile(self, component, locs, vals):
|
||||
self._initialized = 1
|
||||
Stack.setProfile(self, self.flame, component, locs, vals)
|
||||
|
||||
def set(self, tol = None, energy = '', tol_time = None):
|
||||
if tol:
|
||||
self.flame.setTolerances(default = tol)
|
||||
if tol_time:
|
||||
self.flame.setTolerances(default = tol_time, time = 1)
|
||||
if energy:
|
||||
self.flame.set(energy = energy)
|
||||
|
||||
def T(self, point = -1):
|
||||
"""The temperature [K]"""
|
||||
return self.solution('T', point)
|
||||
|
||||
def u(self, point = -1):
|
||||
"""The axial velocity [m/s]"""
|
||||
return self.solution('u', point)
|
||||
|
||||
def V(self, point = -1):
|
||||
"""The radial velocity divided by radius [s^-1]"""
|
||||
return self.solution('V', point)
|
||||
|
||||
def solution(self, component = '', point = -1):
|
||||
"""The solution for one specified component. If a point number
|
||||
is given, return the value of component 'component' at this
|
||||
point. Otherwise, return the entire profile for this
|
||||
component."""
|
||||
if point >= 0: return self.value(self.flame, component, point)
|
||||
else: return self.profile(self.flame, component)
|
||||
|
||||
def setGasState(self, j):
|
||||
nsp = self.gas.nSpecies()
|
||||
y = Numeric.zeros(nsp, 'd')
|
||||
for n in range(nsp):
|
||||
nm = self.gas.speciesName(n)
|
||||
y[n] = self.solution(nm, j)
|
||||
self.gas.setState_TPY(self.T(j), self.pressure, y)
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
104
Cantera/python/Cantera/OneD/StagnationFlow.py
Normal file
104
Cantera/python/Cantera/OneD/StagnationFlow.py
Normal file
|
|
@ -0,0 +1,104 @@
|
|||
from onedim import *
|
||||
import Numeric
|
||||
|
||||
class StagnationFlow(Stack):
|
||||
"""An axisymmetric flow impinging on a surface at normal incidence."""
|
||||
|
||||
def __init__(self, gas = None, surfchem = None, grid = None):
|
||||
self.inlet = Inlet('inlet')
|
||||
self.gas = gas
|
||||
self.surfchem = surfchem
|
||||
self.inlet.set(temperature = gas.temperature())
|
||||
self.surface = Surface(id = 'surface', surface_mech = surfchem)
|
||||
self.pressure = gas.pressure()
|
||||
self.flow = AxisymmetricFlow('flow',gas = gas)
|
||||
self.flow.setupGrid(grid)
|
||||
Stack.__init__(self, [self.inlet, self.flow, self.surface])
|
||||
self.setRefineCriteria()
|
||||
self._initialized = 0
|
||||
|
||||
|
||||
def init(self):
|
||||
"""Set the initial guess for the solution."""
|
||||
self.getInitialSoln()
|
||||
gas = self.gas
|
||||
nsp = gas.nSpecies()
|
||||
yin = Numeric.zeros(nsp, 'd')
|
||||
for k in range(nsp):
|
||||
yin[k] = self.inlet.massFraction(k)
|
||||
gas.setState_TPY(self.inlet.temperature(), self.pressure, yin)
|
||||
u0 = self.inlet.mdot()/gas.density()
|
||||
t0 = self.inlet.temperature()
|
||||
V0 = 0.0
|
||||
|
||||
tsurf = self.surface.temperature()
|
||||
|
||||
zz = self.flow.grid()
|
||||
dz = zz[-1] - zz[0]
|
||||
|
||||
|
||||
locs = Numeric.array([0.0, 1.0],'d')
|
||||
self.setProfile('u', locs, [u0, 0.0])
|
||||
self.setProfile('V', locs, [V0, V0])
|
||||
self.setProfile('T', locs, [t0, tsurf])
|
||||
for n in range(nsp):
|
||||
self.setProfile(gas.speciesName(n), locs, [yin[n], yin[n]])
|
||||
self._initialized = 1
|
||||
|
||||
|
||||
def solve(self, loglevel = 1, refine_grid = 1):
|
||||
if not self._initialized: self.init()
|
||||
Stack.solve(self, loglevel = loglevel, refine_grid = refine_grid)
|
||||
|
||||
|
||||
def setRefineCriteria(self, ratio = 10.0, slope = 0.8,
|
||||
curve = 0.8, prune = 0.0):
|
||||
Stack.setRefineCriteria(self, domain = self.flow,
|
||||
ratio = ratio, slope = slope, curve = curve,
|
||||
prune = prune)
|
||||
|
||||
def setProfile(self, component, locs, vals):
|
||||
self._initialized = 1
|
||||
Stack.setProfile(self, self.flow, component, locs, vals)
|
||||
|
||||
def set(self, tol = None, energy = '', tol_time = None):
|
||||
if tol:
|
||||
self.flow.setTolerances(default = tol)
|
||||
if tol_time:
|
||||
self.flow.setTolerances(default = tol_time, time = 1)
|
||||
if energy:
|
||||
self.flow.set(energy = energy)
|
||||
|
||||
def T(self, point = -1):
|
||||
return self.solution('T', point)
|
||||
|
||||
def u(self, point = -1):
|
||||
return self.solution('u', point)
|
||||
|
||||
def V(self, point = -1):
|
||||
return self.solution('V', point)
|
||||
|
||||
def solution(self, component = '', point = -1):
|
||||
if point >= 0: return self.value(self.flow, component, point)
|
||||
else: return self.profile(self.flow, component)
|
||||
|
||||
def coverages(self):
|
||||
nsurf = self.surfchem.nSpecies()
|
||||
cov = Numeric.zeros(nsurf,'d')
|
||||
for n in range(nsurf):
|
||||
nm = self.surfchem.speciesName(n)
|
||||
cov[n] = self.value(self.surface, nm, 0)
|
||||
return cov
|
||||
|
||||
def setGasState(self, j):
|
||||
nsp = self.gas.nSpecies()
|
||||
y = Numeric.zeros(nsp, 'd')
|
||||
for n in range(nsp):
|
||||
nm = self.gas.speciesName(n)
|
||||
y[n] = self.solution(nm, j)
|
||||
self.gas.setState_TPY(self.T(j), self.pressure, y)
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
|
@ -2,4 +2,8 @@
|
|||
The classes in this package implement one-dimensional reacting flow problems.
|
||||
"""
|
||||
from onedim import *
|
||||
from BurnerFlame import BurnerFlame
|
||||
from CounterFlame import CounterFlame
|
||||
from StagnationFlow import StagnationFlow
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -22,7 +22,8 @@ class Domain1D:
|
|||
return _cantera.domain_type(self._hndl)
|
||||
|
||||
def index(self):
|
||||
"""Index of this domain in a stack. Returns -1 if this domain is not part of a stack."""
|
||||
"""Index of this domain in a stack. Returns -1 if this domain
|
||||
is not part of a stack."""
|
||||
return _cantera.domain_index(self._hndl)
|
||||
|
||||
def nComponents(self):
|
||||
|
|
@ -36,6 +37,12 @@ class Domain1D:
|
|||
def componentName(self, n):
|
||||
"""Name of the nth component."""
|
||||
return _cantera.domain_componentName(self._hndl, n)
|
||||
|
||||
def componentNames(self):
|
||||
names = []
|
||||
for n in range(self.nComponents()):
|
||||
names.append(self.componentName(n))
|
||||
return names
|
||||
|
||||
def componentIndex(self, name):
|
||||
"""Index of the component with name 'name'"""
|
||||
|
|
@ -248,6 +255,7 @@ class Surface(Bdry1D):
|
|||
Bdry1D.__init__(self)
|
||||
if surface_mech:
|
||||
self._hndl = _cantera.reactingsurf_new()
|
||||
self.setKineticsMgr(surface_mech)
|
||||
else:
|
||||
self._hndl = _cantera.surf_new()
|
||||
if id: self.setID(id)
|
||||
|
|
@ -256,9 +264,12 @@ class Surface(Bdry1D):
|
|||
def setKineticsMgr(self, kin):
|
||||
_cantera.reactingsurf_setkineticsmgr(self._hndl,
|
||||
kin.kinetics_hndl())
|
||||
def enableCoverageEqs(self, onoff=1):
|
||||
_cantera.reactingsurf_enableCoverageEqs(self._hndl, onoff)
|
||||
|
||||
|
||||
def setCoverageEqs(self, onoff='on'):
|
||||
if onoff == 'on':
|
||||
_cantera.reactingsurf_enableCoverageEqs(self._hndl, 1)
|
||||
else:
|
||||
_cantera.reactingsurf_enableCoverageEqs(self._hndl, 0)
|
||||
|
||||
|
||||
class AxisymmetricFlow(Domain1D):
|
||||
|
|
@ -382,6 +393,9 @@ class Stack:
|
|||
_cantera.sim1D_setTimeStep(self._hndl, stepsize,
|
||||
Numeric.asarray(nsteps))
|
||||
|
||||
def getInitialSoln(self):
|
||||
_cantera.sim1D_getInitialSoln(self._hndl)
|
||||
|
||||
def solve(self, loglevel=1, refine_grid=1):
|
||||
return _cantera.sim1D_solve(self._hndl, loglevel, refine_grid)
|
||||
|
||||
|
|
@ -416,7 +430,7 @@ class Stack:
|
|||
for n in range(np):
|
||||
x[n] = self.value(domain, component, n)
|
||||
return x
|
||||
|
||||
|
||||
def workValue(self, dom, icomp, localPoint):
|
||||
idom = dom.index()
|
||||
return _cantera.sim1D_workValue(self._hndl, idom, icomp, localPoint)
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue