adiabatic_flame.py

[Jon Tegner]

Hi, I have some questions regarding this code (and pleas excuse me if they ar stupid). We are using adiabatic_flame.py to evaluate different kinetics (hydrogen, propane, methane), and we are often having problem generating a solution. Generally I would be grateful for any hints on how to obtain solutions. Typically, what I do now is varying the initial grid (number of points, length (and previously also tfix)) - but if there is a more "systematic" way of doing this it would certainly be useful.

I also have some more specific questions (on the version of adiabatic_flame.py included in Cantera-2.1).

1. In previous versions (2.0) there was a variable "tfix" which now seems to be missing, is it not needed (or is it specified in a different way)?

2. In the code initial_grid = [0.0, 0.001, 0.01, 0.02, 0.029, 0.03], i.e., grid spacing is refined at the boundaries. If I use a different grid I have noticed (for some cases) that it is very sensitive to HOW the grid is refined there. Is there a "best practice" of how to do this refinement? Would it be advantageous to involve more than one point in this refinement?

3. The flame is solved by calling the "solve" method, is there something written on how this method works? I guess it is using a Newton method (to solve for the steady propagating flame?), but since a time step variable is also specified (f.set_time_step(1e-5, [2, 5, 10, 20, 50, 100]), it seems a time dependent approach is also used, and if so, how do these methods work together? And what does the list above (n_steps [2,5,10,20,50,100]) specify? Doc-string: "Set the sequence of time steps to try when Newton fails." Does it (for this case) take 2 steps, do a Newton, take 5 steps etc.? Would it be possible to get the flame profile after these steps have been taken? It might give a hint whether the computational domain is to large or to small (in the case the procedure fails).

Regards,

/jon

[Jon Tegner]

Continuing on the thread above: for one of the hydrogen-kinetics I'm evaluating I noticed that I managed to obtain solutions for a range of equivalence numbers (really quick) if I omitted the step where the the transport model is sat to 'Multi'. Two questions:

1. What are the implications on the quality of the solutions when using mixture-averaged formulation as opposed to multicomponent?

2. Is it reasonable that it is a lot more difficult to obtain a solution when 'Multi' is used? And if so, why?

Thanks,

/jon

[Ray Speth]

Jon,

The replacement for the "tfix" variable is the "set_fixed_temperature" method of class Sim1D (or FreeFlame). Unlike the old Python interface, there is now a default value for this parameter, halfway between the reactant temperature and the adiabatic flame temperature, so I don't think it should generally be necessary to set a value for this parameter. Please let me know if you find a case where this is not true.

In theory, the initial grid shouldn't matter much, since it's going to be refined algorithmically anyway. The cases where you would want to change it are (1) the solver doesn't converge or is slow to converge on the initial grid, or (2) the solver converges to an undesirable solution (i.e. the nonreacting solution of the counterflow diffusion flame).

There's not much in the way of documentation for how the solver works. After discretizing onto the grid, the steady-state flame is a system of algebraic equations. The nonlinearities in the equations are severe enough that it's non-trivial to find an initial condition for which Newton's method will converge. This is where the timestepping comes in. For a sufficiently small timestep, using Newton's method to solve the time-dependent problem should give a solution that is a better initial guess for solving the steady-state problem.

The list in the "set_time_step" method is how many time steps to take before trying to solve the steady-state problem again. For example, after the initial attempt to solve the steady state problem, The solver will attempt to take 2 timesteps of 1e-5 seconds. It will reduce the size of the timestep until it takes a successful step. Then it will re-attempt the steady state problem. If it fails again, the solver will take 5 timesteps, try to solve thesteady state, take 10 timesteps, etc, until it can successfully solve the steady state problem. After each set of successful timesteps, it will attempt to increase the size of the timestep, as larger timesteps are generally feasible as you approach the steady-state solution.

If you set the 'loglevel' to 7, the solver will output a file 'debug_sim1d.xml' containing the solution after each successful or unsuccessful Newton solve, and after each group of timesteps. If you set 'loglevel' to 8, this file will also contain the residuals for each equation. I can't promise that you'll find this useful -- it's difficult to tell why a particular solution guess isn't good enough. I'm pretty sure there are some bugs in the solver that prevent it from converging in cases where it shouldn't have any particular trouble.

Regards,

Ray

[Bryan W. Weber]

Hi Jon,

Ray gave you some excellent advice about the solvers, and I think I can answer your questions about the diffusion coefficients. You can look at the equations here: http://public.ca.sandia.gov/chemkin/docs/oppdif.pdf on page 8. Equation 6 is for multi-component, Equation 7 is for mixture averaged.

Essentially the difference is that the multi-component diffusion computes the diffusion velocity of every species with respect to every other species - you can see that this problem grows very quickly with the number of species. It is an m choose n, or binomial problem. If we have K species and binary diffusion, then we have a K choose 2 problem. Here is a plot of this problem: http://www.wolframalpha.com/share/clip?f=d41d8cd98f00b204e9800998ecf8427eqa6ko4cp0g You can see that the number of equations to be solved grows very very rapidly as the number of species increases.

The mixture averaged approach, however, computes the diffusion velocity of each species with respect to an average mixture, and thus grows only linearly with the number of species. Here is a plot comparing the number of equations on the scale of a hydrogen mechanism: http://www.wolframalpha.com/share/clip?f=d41d8cd98f00b204e9800998ecf8427e2bbtjp0afe Thus, you can see that it is much more difficult to solve the multi-component diffusion problem.

That said, the multi-component diffusion is generally more accurate. How much more, and how much it matters for your problem, depends on the problem that you are trying to solve. One thing you can do (and this relates to some of your previous questions as well) is use a previously converged solution as the base for a more difficult solution. Thus, you can take your converged mixture-averaged result and switch the transport modeling. The solver will not have as hard a time to solve this problem because it starts with a good guess and it may take less time to solve the multi-component problem than if you started the multi-component problem from scratch. However, the operative word here is *may*. You can see that this is the default behavior in adiabatic_flame.py as well.

I hope that helps,

Bryan

[Jon Tegner]

Ray and Bryan, thank you very much! Clarifications are VERY appreciated!

However, I still do have some questions.

It seems the code adiabatic_flame.py adapts the technique Bryan indicated, i.e., after obtaining a converged mixture-averaged solution switching to the multi-component approach.

However, as it is coded in adiabatic_flame.py supplied with Cantera-2.1 it seems that the "Multi" solution is calculated twice, i.e., *Multi' seems to be the default when the FreeFlame object is created, and even though code states that it gives mixture-averaged flamespeed it in fact prints the multicomponent result (and actually does so twice). Or maybe I'm missing something here?

I guess my question here boils down to how to change this to the intended approach - using "Mix-solution" as initial guess for the "Multi-solution" (or a clarification of my mistakes).

I tried putting in "gas.transport_model = 'Mix'" before the creation of the FreeFlame object, and eveh though this seems to result that in the first run the mixture averaged solution is obtained, I don't get a solution when switching to "Multi".

Again, I really appreciate that you have spent some time helping me with this, thanks!

/jon

--
You received this message because you are subscribed to the Google Groups "Cantera User's Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cantera-user...@googlegroups.com.
To post to this group, send email to canter...@googlegroups.com.
Visit this group at http://groups.google.com/group/cantera-users.
For more options, visit https://groups.google.com/groups/opt_out.

[Bryan W. Weber]

Hi Jon,

You can test which is the default transport model in the following way:

import cantera as ct
gas = ct.Solution('h2o2.xml')
gas.transport_model
f = ct.FreeFlame(gas)
f.transport_model

On my install of 2.1, I got 'Mix' for both of those when I started from scratch. The default transport model is actually set by the phase when it is imported with the Solution class, so it will depend on what the default is in your XML or CTI file. Changing the transport model either by f.transport_model = 'yyyy' or by gas.transport_model = 'yyyy' changes the transport model of the phase associated with the flame or gas object. Thus, if you do this

f.transport_model = 'Multi'
f.transport_model
gas.transport_model
f2 = ct.Freeflame(gas)
f2.transport_model

you should find that they all say 'Multi', even though you've created a whole new flame! This is because the transport model associated with the gas phase object has been changed by f.transport_model = 'Multi' command.

What do you mean by "don't get a solution"? Is there an error that is printed, or does it take a really long time and you can't wait for it to finish? If its the second one, you might want to try either increasing the refine_criteria, turning off grid refinement, or disabling the energy equation (or all three) to see if it helps. Then you can undo whichever changes you've made after a solution is found and try to solve again.

Hope it helps!

Bryan

[Jon Tegner]

Thanks again, didn't realize that it defaults to what's in the CTI-file. This is good to know (and it might be a good idea to set it explicitly in adiabatic_flame.py).

And the code actually crashed, I'll try to set the transport model as you indicated and see if the situation improves.

Thanks!

/jon

[Ray Speth]

Jon,

I agree, it would be more clear to explicitly set the transport model in the example, and more likely to work correctly when users adapt that script to use with their own model.

If you have an example that's actually crashing, could you post your script and mechanism file?

Regards,

Ray

[Jon Tegner]

I had an example which crashed during the multi-component part, however, that seemed to be a result of species present in the definition of the mixture which were neither part of the reaction (one step in this case), nor initialized. That is, in my cti-file I hade the line

species = 'H2 O2 O OH H H2O N2',

in the definition of the mixture, while the reaction part looked like:

reaction( 'H2 + 0.5 O2 => H2O', ...)

So when I removed O, OH and H the multi-component run went through without problem. Don't know if this to be expected, but it does not seem unreasonable to me (but if you want to try them I'd be happy to supply them).

Also tried a mechanism with 7 reaction and 8 species, no problem here either (taking longer time to complete though).

Regards,

/jon

[Ray Speth]

Jon,

Yes, please post your .cti file and modified .py file. I don't think there are any cases where just crashing is a reasonable behavior for the solver.

One thing you can try if you're finding that the multicomponent run converges slowly is to solve it with successively tighter tolerances. I think this idea was suggested by someone else in this group previously, and it seems to work in some cases. For example, do something like this, after solving the mixture-averaged case:

f.transport_model = 'Multi'

f.flame.set_steady_tolerances(default=[1e-3, 1e-12])
f.flame.set_transient_tolerances(default=[1e-2, 1e-12])
f.solve(loglevel, refine_grid)


f.flame.set_steady_tolerances(default=[1e-4, 1e-12])
f.flame.set_transient_tolerances(default=[1e-3, 1e-12])
f.solve(loglevel, refine_grid)


f.flame.set_steady_tolerances(default=[1e-5, 1e-12])
f.flame.set_transient_tolerances(default=[1e-4, 1e-12])
f.solve(loglevel, refine_grid)

Regards,

Ray

[Jon Tegner]

Ray,

the files can be obtained through

wget www.renget.se/asdf/adiabatic_flame.py

and

wget www.renget.se/asdf/mixture.cti

With this mixture transport_model is sat to 'Multi' (I really should have known this...). The code - i principle the same as the one included in 2.1 - calculates a mixture averaged solution and a multi component one. However, since I have defined transport_model='Multi' in mixture.cti it actually prints the multi component solution twice.

Next - for testing - I sat transport_model='Mix' at the last run. This also goes through without problem (even though what the output indicates as "mixture-solution" is actually the "multi-solution" and vice versa).

However, if I set transport_model='Mix' at the first occasion and transport_model='Multi' at the second occasion (according to how the code is supposed to run) I haven't reached a solution in 10 minutes (the cases above finish in a few seconds). The output indicates that the "mixture-solution" is obtained.

Probably doing something stupid here, and would be very grateful for all hints. Thanks!

/jon

[dirk.g...@googlemail.com]

Ray,

the troubles Jon is reporting are very well known to me as well. I am fighting a lot with the 

[dirk.g...@googlemail.com]

Ray,

the troubles Jon is reporting are very well known to me as well. I am struggeling quite a lot to get tolerance settings in a way that the 1D simulations are running, often using tolerances which are by a number of magnitudes smaller than in the exampble codes.

Henry Moffat wrote recently in this forum:
"And, you’ve correctly point out that the 1D capabilities within Cantera aren’t very good. I don’t myself use them in my research, having written better routines within Chemkin 10 years ago. I’ll say more on this later in a separate email. "

It looks like the 1D capabilities are important for a number of Cantera users.  Is there a chance that the 1D Newton solver (this is the critical part as far as I understand ?) within Cantera is getting improved ?  I have no idea how big the effort would be, but I am very interested to get an idea about the pro an cons of the 1D capabilities of Cantera.

Greetings,
Dirk

On Thursday, January 9, 2014 12:57:31 PM UTC+1, Jon wrote:

[ober...@web.de]

Hello guys,

I am also getting Trouble with the tolerances. I have to compute the laminar flame Speed for big amount of different operating points of an internal combustion engine. I am using the adiabatic_flame.py and sometime I get good results but most of the time it ends with this error:  

"CanteraError thrown by OneDim::timeStep:
Time integration failed."

I tried to successively tighter tolerances like you told Ray, but even this doesn't work. Is there anybody who found a good solution for this Problem? Any methods to change the tolerances in a reasonable way?

Please let me know if there is a proper way to do it. Thank you very much in advance, Tobi

[Ray Speth]

Tobi,

What version of Cantera are you using? I believe the underlying issues mentioned in this thread have been resolved in Cantera 2.2.0. If you are using Cantera 2.2.0 and have a specific example that's causing trouble, can you provide the CTI file and a Python script that demonstrates one or more of the troublesome cases?

Regards,

Ray