I am arranging a numerical diffusion flame in Matlab, with Incompressible N-S and Cranck-Nicholson. The program is able to make few iterations. After that, it is unstable and in particular the energy residual grows up a lot.

I think it has to do with the very sharp temperature gradient formed in the mixing layer.

have you got any suggestion to stabilize my algorithm?

Thanks.

Do you use relaxation? Solving exothermic reactive flows is known to desabilize the solution. Adding relaxation is a simple and usually quite effective remedy at the expense of more iterations to reach convergence.

Thanks Rami. I've used under-relaxation factors for each equation. I have obtained a weak-reactive solution, but I had to shorten the non-dimensional activation energy (Ea/RTo) to approximately 10. This enhances smaller temperature gradients. But I fix Ea/RTo of the order of 100 (is more accurate), the iteration has a lot of unstabilities and do not reach convergence.

Maybe stronger relaxation is needed. I had similar difficulties some years ago. My relaxation was of the false time step type, and I had to decrease the time steps by an order of magnitude or more together with some additional tricks to achieve convergence. One such useful trick was to start with inert solution to get some mixing, and only then switch on the reaction.

BTW, are you solving simultaneously for all species? I think this also enhances the numerical stability.

Just out of curiosity: Is your flow turbulent? Is the reaction time scale comparable/longer than the mixing time scale (if not, maybe a mixing-controlled reaction is more appropriate)? How many reactions and species are you solving (the less the easier it is to obtain converged solutions)?

I have solving steady equations for incompressible flow (continuity, x-Momentum, Energy, Species-Oxidizer, Species-Fuel). It can be used the mixture fraction, but I want to simulate all this equations for becoming myself stronger doing that. It is a one-step chemical reaction, a very simple mechanism

F+sO2---> products(1+s)+q(Joules)

Due to their parabolic behaviour, I'm using the TDMA algorithm for each abcissae, linearizing the coeficients u,v,and the heat generated by the exponential term in a Crank-Nicholson way. Also, the #Damkhöler is 1, adiabatic temperature 6, Mach<<<1, non-dimensional activation energy of 100. And I have reached under relaxation factors of 0.3.

if your reaction scheme has negative exponent for fuel (as it is typical for 1 step methane equation), you will have big problems with stability.

If this is the case, try another fuel. Otherwise go even lower with relaxation for species and temperature, see where the instability comes from - is it the enthalpy source? Try more robust - implicit then crank-nicholson scheme. Combustion with kinetics is tough.

matej

Thanks matej for your reply. I will try something....